## Abstract

We obtained likelihoods in the lower mantle for long-wavelength models of bulk sound and shear wave speed, density, and boundary topography, compatible with gravity constraints, from normal mode splitting functions and surface wave data. Taking into account the large uncertainties in Earth's thermodynamic reference state and the published range of mineral physics data, we converted the tomographic likelihoods into probability density functions for temperature, perovskite, and iron variations. Temperature and composition can be separated, showing that chemical variations contribute to the overall buoyancy and are dominant in the lower 1000 kilometers of the mantle.

To understand the nature of mantle convection, it is essential to quantify thermal and compositional contributions to the density variations that drive the solid-state flow. Although seismic tomography is probably the best probe for Earth's three-dimensional structure, its main constraint is on wave speeds rather than density. It has therefore been common practice in tomography to prescribe a scaling between density and velocity variations (*1*) and invert for velocity only. Such a scaling is justified if a single cause is responsible for the observed variations. Temperature-induced ratios of relative density to relative shear wave speed variations between 0.2 to 0.4 have been measured (*2*, *3*) and are compatible with geodynamic data, combined with specific viscosity profiles (*4*, *5*). This, together with evidence from seismology that slabs penetrate deep into the mantle (*6*, *7*), led to the view that mantle dynamics is dominated by thermally driven whole-mantle convection (*8*). Chemical buoyancy (*9*) was introduced mainly to explore the possible thermochemical nature of D″ in terms of a primordial layer (*10*–*14*), subducted oceanic crust (*14*–*16*), or chemical reactions with the core (*17*). These Boussinesq calculations, however, are not realistic, because the simulated high-density contrasts are not compatible with the observed seismic velocities and a plausible mineralogic model (*18*). When an extended Boussinesq or compressible calculation is used, the required density contrasts are reduced (*19*–*21*). More interestingly, in models where thermal expansivity decreases with depth, thermochemical superplumes are seen to develop (*20*–*22*), not unlike those found under Africa and the Pacific in tomography (*23*, *24*). With improving resolution of seismic velocities and, especially, of density, indirect evidence has emerged suggesting that compositional heterogeneity is present in the lower mantle (*5*, *25*–*32*). In an effort to reconcile evidence from various research fields, dynamical models with a strong compositional component (*33*–*35*) have challenged the classic view of thermally driven mantle convection.

Owing to trade-offs between temperature and composition, wave speeds alone are not sufficient to infer their variations, and density constraints should be included (*5*, *31*, *36*). Normal modes require weak and/or negative correlations between density and shear wave speed variations throughout most of the lower mantle (*28*), but amplitudes of density are difficult to infer (*37*). We represent the seismic constraints with more complete likelihoods, rather than individual models, and have extended the work of Resovsky and Trampert (*32*) to spherical harmonic degree 6 for relative variations of bulk sound (*d*ln*V*_{Φ}) and shear wave speed (*d*ln*V*_{s}), density (*d*lnρ), and topography at the 670-km discontinuity and at the core-mantle boundary (CMB). In addition to providing a full uncertainty analysis (errors and trade-offs), representing the data as likelihoods of seismic parameters allows a subsequent incorporation of additional data constraints. Most often, geodynamic data are jointly inverted with the seismologic data (*5*, *28*), but in our approach, it is more efficient to filter a posteriori the purely seismic models by retaining only those that fit the gravity field within its error bars (*37*). The final likelihoods for wave speeds and density variations are close to Gaussian and can thus be represented by a mean and standard deviation (Table 1). They are a complete and compact representation of all long-period seismic data, compatible with the observed gravity field, and robust (defined here as >1 SD) mean variations of wave speeds agree with previous work (*28*, *29*) (Fig. 1).

Layer (km depth) | dT (K) | dPv (%) | dFe (%) | dlnV_{s} (%) | dlnV_{Π} (%) | dlnρ (%) |
---|---|---|---|---|---|---|

670-1200 | 180 | 5.5 | 0.75 | 0.16 | 0.34 | 0.26 |

1200-2000 | 112 | 3.0 | 0.55 | 0.12 | 0.22 | 0.28 |

2000-2891 | 198 | 3.6 | 0.86 | 0.12 | 0.26 | 0.48 |

Although correlations and average amplitudes can be a good indicator of chemical heterogeneities (*37*), we directly inverted the seismic likelihoods for variations of temperature (*dT*) and composition. Describing the chemical variations by the relative variations of total perovskite (*dPv*) and total iron (*dFe*) content in the lower mantle (*37*), the seismic likelihoods are related to probability density functions (pdfs) for temperature and composition by: (1) (2) (3) where the partial derivatives are the sensitivities of velocities and density to temperature and composition. We calculated sensitivities (fig. S1) using available mineral physics data and a reasonable range for the thermal and chemical reference state of the mantle (*37*). This leads to uncertainties in the sensitivities that are also close to Gaussian.

Solving the algebraic system (Eqs. 1 to 3) would be trivial if some quantities were not pdfs. No routine mathematical tools are available to solve such a system. We solved the system for fixed sensitivities, which implies that the lateral variations in *dT, dPv*, and *dFe* are Gaussian distributed. In this case, the system can be written in vector form as *d* = *Gm*, where *d* represents the mean of the seismic likelihoods and *G* the partial derivatives in the system (Eqs. 1 to 3). The mean thermochemical model is found by *m̄* = *G*^{–1}*d* and the variance is given by *C*_{m̄} = *G*^{–1}*C _{d}*(

*G*

^{–1})* (

*38*), where

*C*

_{d}represents the variance of the seismic likelihoods and * denotes the matrix transpose. Of course,

*G*is not a constant matrix, but each sensitivity is a pdf itself. We therefore solved the system a million times, drawing randomly in the Gaussian distributions of the partial derivatives. We then averaged the mean models

*m*and their variances

*C*

_{m}and determined the corresponding spread. This allows us to distinguish between contributions from the widths of the seismic likelihoods and from uncertainties in the sensitivities in the final model. The average model variance is taken as the sum of all variances; however, the variance in all cases is dominated by

*C*

_{d}. Uncertainties in the sensitivities contribute less than 10% to the total model uncertainties.

For each step (seismic models, sensitivities, thermo-chemical models), we determine complete uncertainties, providing the tools to quantify how meaningful the results are. Robust mean variations of temperature, perovskite, and iron are shown in Fig. 2. Uncertainties in these maps are uniform within a given layer and are listed in Table 1. The large uncertainties in temperature are mainly due to the large uncertainties in density. Nevertheless, density is indispensable to infer the compositional variations. Previous estimates of thermal and chemical parameters without density (*5*, *31*) were not robust (*31*). We obtained robust temperature and compositional anomalies in many places (Fig. 2). Iron and temperature variations are largest in the bottom 1000 km of the mantle. Perovskite variations are smallest in our mid-mantle layer. These characteristic changes with depth are consistent with previous suggestions based on seismologic (*27*) and rheologic (*5*) properties of the mantle. Our estimates are valid for long wavelengths (spherical harmonic degrees 2, 4, 6, and vertical layers of ≈1000 km). Because current tomographic models do not show much power beyond degree 6 (*39*), we do not expect significant lateral changes, but a finer vertical resolution could concentrate the thermochemical signal over smaller vertical length scales. The recent discovery of a postperovskite phase at the bottom of the mantle (*40*) is not likely to change our inferences, because the expected change in elastic properties (*41*) falls inside our range of input parameters used for the calculation of the sensitivities. In most interpretations of tomography, it is assumed that wave speeds, and in particular shear wave speeds, can be scaled to temperature. Our results show that this is not the case (*37*) (fig. S2) and explain why thermochemical inferences obtained without density (*5*, *31*) are so different from those presented here.

Density anomalies generate buoyancy forces that drive mantle flow. Because we determined independent likelihoods for thermochemical variations in the mantle, we are in a position to separate the driving force into thermal and chemical contributions. We resampled the pdfs for *dT, dPv*, and *dFe* and the corresponding density sensitivities to obtain separate likelihoods for thermal and compositional parts of the density (Fig. 3). In the lower 1000 km of the mantle, thermal buoyancy is weak compared to chemical buoyancy, even though temperature variations are highest in this layer. This is because the thermal expansivity decreases with increasing depth (fig. S1). The Pacific and African superplumes, identified in global tomography (*23*, *24*), are dense and have a chemical origin as previously suggested (*20*–*22*, *28*, *29*, *42*). Unless complex, as-yet-unmodeled processes are at play, our findings rule out that superplumes are thermally buoyant, as has often been proposed (*5*, *43*). In the mid- and upper-lower mantle layers, thermal and chemical buoyancies are equally important. Other important concepts of mantle flow are slabs and hotspots. Very few slabs, identified from seismic tomography (*44*), continuously plot on heavier-than-average material down to the CMB, and hardly any hotspots, recently classified as coming from the CMB (*45*, *46*), continuously plot on buoyant material throughout the lower mantle.

Our results demonstrate that long-wavelength chemical heterogeneities exist and play an important role throughout the mantle. Our ability to separate the total buoyancy into thermal and chemical components further shows that compositional variations play a first-order role in large-scale mantle dynamics and cannot be ignored. Recent experimental (*47*) and numerical (*35*) models of thermochemical convection show a variety of structures, depending on input parameters. Whether the structures identified in our models are chemically stable (*33*), organized in piles (*20*), or in a doming regime (*47*) will crucially depend on the input parameters of the simulations: the initial composition of Earth; depth-dependent parameters such as thermal expansivity, thermal conductivity, viscosity; and many other parameters.

**Supporting Online Material**

www.sciencemag.org/cgi/content/full/306/5697/853/DC1

Materials and Methods

Figs. S1 and S2

Table S1

References