Insights into the nature of the transition zone from physically constrained inversion of long-period seismic data

  1. Fabio Cammarano* and
  2. Barbara Romanowicz
  1. Berkeley Seismological Laboratory, University of California, 215 McCone Hall, Berkeley, CA 94720
  1. Edited by Ho-kwang Mao, Carnegie Institution of Washington, Washington, DC, and approved February 7, 2007 (received for review January 4, 2007)

Abstract

Imposing a thermal and compositional significance to the outcome of the inversion of seismic data facilitates their interpretation. Using long-period seismic waveforms and an inversion approach that includes constraints from mineral physics, we find that lateral variations of temperature can explain a large part of the data in the upper mantle. The additional compositional signature of cratons emerges in the global model as well. Above 300 km, we obtain seismic geotherms that span the range of expected temperatures in various tectonic regions. Absolute velocities and gradients with depth are well constrained by the seismic data throughout the upper mantle, except near discontinuities. The seismic data are consistent with a slower transition zone and an overall faster shallow upper mantle, which is not compatible with a homogenous dry pyrolite composition. A gradual enrichment with depth in a garnet-rich component helps to reduce the observed discrepancies. A hydrated transition zone would help to lower the velocities in the transition zone, but it does not explain the seismic structure above it.

The thermal state and composition of Earth's upper mantle and transition zone dictate its dynamics from microscale (e.g., creep mechanisms and earthquakes) to macroscale (e.g., modality of mantle convection and plate tectonics). The knowledge of these fundamental physical parameters and their 3D variations in Earth's deep interior is indirect and relies entirely on the interpretation of geophysical data, based on insights from theoretical and experimental mineral physics. Among these data, seismological observations constitute a main source of information. Seismic waves record information about the elastic (and anelastic) structure of Earth. Long-period seismic data provide the most comprehensive global constraints on upper mantle shear velocity structure. Fundamental-mode surface waves are mostly sensitive to the uppermost mantle structure, and including overtones provides resolution in the transition zone.

Lateral temperature variations in the Earth have been inferred since the first tomographic studies began to reconstruct 3D seismic velocity structure (1, 2). However, despite the ever-improving resolution of seismic velocity models and a general agreement of different models on at least the large-scale structure (e.g., refs. 36), interpretation is still challenging. The few quantitative interpretation efforts made to date have shown that much of the seismic heterogeneity in the uppermost mantle is probably thermal in origin (7, 8), but it is likely that there is a compositional contribution under cratons (9, 10) and that fluids influence the very low velocities in the mantle wedge above subduction zones (8). By adding constraints from gravity or geoid data, joint inversions for seismic velocity and density have been performed because of the better chance to isolate thermal and compositional effects (1113). Whereas temperature sensitivity dominates for seismic velocities, especially at shallow depths (7, 14, 15), density also strongly varies with composition.

Besides the trade-off between temperature and composition, an important issue for seismic interpretation concerns the fact that a seismic model is not required to have a meaning in terms of these physical parameters.

Typically, the interpretation of a given seismic data set is performed in two steps. First, a seismic model is constructed that fits the data satisfactorily. Second, the model is interpreted in terms of temperature and composition based on the knowledge of the elastic and anelastic properties of mantle minerals plus constraints on plausible composition and temperature ranges from geochemistry [i.e., the signature of outcropping rocks (orogenic peridotites) and mantle inclusions (xenoliths and xenocrysts)], heat flow measurements at the surface, and conditions for melting mantle materials. The second part of the process is commonly regarded as critical for interpretation. Together with the trade-off between temperature and composition, other factors should be considered. Temperature-dependent anelasticity implies nonlinearity of the partial derivatives with temperature in the upper mantle (1). Consequently, absolute velocities are required for a correct interpretation. The presence of mineralogical phase transitions further complicates the interpretation of the transition zone structure (14, 16). Finally, the large uncertainties that still characterize the elastic and, even more, anelastic parameters derived by mineral physics (see refs. 14, 17, and 18) hamper a precise interpretation. Instead, it is often forgotten that the seismic models are already an interpretation of the data. The seismic models are not unique, and they depend on the parametrization, distribution, quality, and type of data used. More importantly, a given best-fit seismic model does not necessarily correspond to a physical model. For example, the previous work of Cammarano et al. (14) has pointed out the difficulties of interpreting seismic reference models [Preliminary Reference Earth Model (PREM) (19) and AK135 (20)] in terms of temperature and composition. The same problems are propagated to 3D tomography models that are obtained by perturbing a starting average model.

To overcome this problem, it is indeed important to test a physical hypothesis directly against the seismic data (as in refs. 21 and 22). The general idea is to use the mineral physics information already at an early stage of the seismic data processing, hence imposing a physical significance to the seismic tomography model.

Here we invert long-period seismic waveforms, which are mainly sensitive to shear velocity, with respect to a physical reference model instead of a standard 1D seismic reference model (e.g., PREM). Velocity variations will thus correspond to thermal (or compositional) variations and a consistent 3D density and V P structure may be determined. We start by assuming purely thermal variations because, in the upper mantle, composition plays a secondary role. We compare the thermal features constrained by the seismic data against expected thermal variations and estimated geotherms in various tectonic regions. The average velocity model extracted from the physically constrained tomography provides insights into the nature of the transition zone. To assess the uncertainties in the thermal model and estimate how well average velocity and velocity gradient with depth are constrained by long-period seismic data in the upper mantle, we perform a series of inversions starting with various models.

Seismic Procedure

The Starting Physical Model.

Because seismic inversion is based on perturbation theory, to be used as a reference model, the starting physical model must have an overall good fit to seismic data, comparable to the fit of radially symmetric seismic models (e.g., PREM). We tested our inversion by using some of the Physical Reference (PREF) models from Cammarano et al. (22). These are isotropic models with the same thermal and compositional structure, respectively 60-million-year-old oceanic geotherm in the lithosphere and a 1,300°C adiabat below that and pyrolite composition. The set of elastic and anelastic parameters for the mantle minerals were obtained from a Monte Carlo search within estimated bounds from ref. 14. These bounds encompass the total range of the experimental values existing at the time of the compilation (3).

Only the models that fit seismic travel times (P and S arrivals) and fundamental mode eigenfrequencies as well as the seismic reference models have been selected. Despite a variety of possible combinations of mineral physics parameters found (see ref. 21), the PREF models are very similar seismically [supporting information (SI) Fig. 4], and they are very different from the seismic reference models.

For example, all of the models have a low velocity zone in the uppermost mantle, where temperature effects prevail over pressure effects. The models do not have any sharp discontinuity at 220 km as in PREM; they have a larger jump than seismic reference models ≈410 km, associated with the olivine–wadsleyite phase transition, and a faster transition zone. All these features are governed by the given composition and thermal structure plus the strong constraint on the integrated average of upper mantle velocities given by the teleseimic travel times. A slower transition zone, balanced by an overall faster upper mantle above, would explain the travel time data better, but the mineral physics constraints do not allow that (SI Fig. 4) (22).

Because of the sparse coverage of travel times at far-regional distances (epicentral distances <30°), it was not possible to determine in that study (22) if such transition zone structure was required at a global scale. With the long-period overtone data used here, we were able to better constrain the average structure of the transition zone and draw conclusions on its nature.

We selected three PREF models that represent well the differences in seismic and mineral physics properties within the family of PREF models. The use of different starting PREF models, characterized by their own sensitivities of seismic velocities to temperature, helped assess the effects of mineral physics uncertainties on the physically constrained inversion.

At the same time, this method tested the assumption that our outcome does not depend on the starting model. To further investigate this point, we also inverted the seismic data starting with a model that has a much smaller jump at 410 km and different velocity gradients above and below this discontinuity (see Fig. 3, model B).

Data and Method.

The data we used included three-component, long-period fundamental waveforms and higher-order mode surface waveforms collected for global tomography efforts over the years at the University of California, Berkeley (5, 6, 23, 24). Events with moment magnitude larger than six and stations at epicentral distance between 15° and 165° were selected. The original seismograms have been deconvolved for the instrument response and filtered between 60 s and a variable maximum period, typically between 220 s and 1 h, which was chosen according to the event magnitude. Wave-packets for both fundamental and overtones were then extracted from the seismograms based on the comparison of the observed trace with the PREM synthetics (5). The PREM represents very well the average seismic structure of Earth. The bounds of the selection criteria are large enough to make the procedure appropriate also for the alternative average upper mantle structure based on a physical model. Minor and major arc paths were selected for our inversion. In total, the data set amounts to 39,829 fundamental waveforms and 59,831 overtones. The windowing scheme not only reduced the computation time of the inversion, but it allowed the application of appropriate weights to different phases. In our case, we assigned a double weight to the overtones compared with the dominant fundamental modes. Noise and path redundancy were also considered in the weighting scheme as described by Li and Romanowicz (6).

The seismic waveform tomographic method is based on the Nonlinear Asymptotic Coupling Theory method (25), which is a normal-mode perturbation approach that accounts for coupling both along and across dispersion branches. The 2D broadband sensitivity kernels computed with this methodology reproduce the body-wave character of the seismic wave propagation and are therefore able to accurately model the overtones.

Inversion.

We started the inversion with a low-resolution version of a recent anisotropic model (SAW642AN; ref. 5), but we replaced the background seismic reference model with a physical reference model (i.e., PREF). SAW642AN is parametrized in spherical splines with 642 equally spaced horizontal nodes and 16 vertical nodes from the top of the mantle to the core–mantle boundary. We reparametrized the model to 162 horizontal nodes, which corresponds to a spatial resolution of ≈3,500 km, and we kept the original vertical parametrization. At this resolution (corresponding to ∼12° in a spherical harmonics parametrization), global models from different groups are very similar. We inverted for the isotropic part of the model down to 1,000 km, and we kept fixed the radially anisotropic part of the model, represented by ξ (V SH 2/V SV 2). Radial anisotropy is required to simultaneously fit spheroidal and toroidal fundamental modes. We corrected for crustal structure by assuming the model CRUST2.0 (26). The inversion procedure consists of several iterations. For each step, we solved the inverse problem by using the classical least-squares approach (27). Regularization schemes involve scaling the equations according to the weighting scheme discussed and then damping. The number of iterations required to reach a variance reduction similar to the starting model SAW642AN range from four to 11 iterations for the different PREF models tested. The selected models do not fit both spheroidal and toroidal modes equally well. Because we do not vary the radial anisotropic structure, we map this discrepancy into the isotropic part and therefore into temperature. More discussion about the trade-offs with anisotropy structure follows (see Absolute temperatures below).

Although we did not invert for seismic attenuation, it is important to account for the effects of anelasticity on seismic velocities to retrieve a correct thermal interpretation. Pressure- and temperature-dependent anelasticity, defined as for the PREF models (14, 21) is included in our thermal models.

Results

Thermal Model.

Lateral variations.

The lateral variations in temperature (Fig. 1) confirm the primary role of temperature in determining the seismic structure of the upper mantle and transition zone. The thermal variations that explain the seismic data are indeed consistent with the range of temperatures expected from geodynamics modeling (28). However, the secondary compositional effect clearly emerges when one looks at the very low (i.e., as low as 460 K at 100 km and 500 K at 150 km) values of temperature obtained in the first 250 km beneath cratons. If more depleted composition (e.g., harzburgite) is considered for those regions, we expect a shift upward of the lowest temperatures of ≈200° (14, 29). Unlike in velocity models, the average temperature at a given depth is not balanced around the mean value but is shifted toward the hot regions (Fig. 1) because of the temperature-dependent anelasticity, which is included in our models. Accounting for its effects on seismic velocities helps to retrieve a more realistic thermal interpretation. In particular, because of different thermal structures beneath oceans and continents, attenuation is expected to have a very strong lateral heterogeneity, which is consistent with seismic attenuation global models (24). In our derived attenuation model, the lateral variations depend only on temperature; therefore, they resemble exactly the thermal variations shown in Fig. 1.

Fig. 1.

Lateral temperature variations inferred from the physically constrained waveform inversion. Positive anomalies in the uppermost mantle correspond to oceanic ridges, and negative anomalies correspond to cratonic areas. Negative thermal anomalies related to subduction zones appear below. The minimum and maximum values of anomaly compared with the mean temperature value is given for each depth.


The lateral thermal variations obtained by starting the inversion with the other two PREF models are very similar, indicating that the additional uncertainties in the temperature partial derivatives do not significantly affect the interpretation. Specifically, these uncertainties affect it much less than regularization schemes required by the inversion process.

Absolute temperatures.

The upper mantle geotherms constrained by long-period seismic data are able to reproduce the range of expected geotherms beneath oceans and cratons (Fig. 2). We compared the seismic geotherms with geotherms for different oceanic ages, which are based on the plate cooling model (30), and purely conductive continental geotherms, which were computed at steady-state and are based on surface heat flow and radiogenic heat production in the crust (31). Cratonic geotherms span a larger range than the oceanic ones, probably mostly due to the compositional differences between them, suggesting that it is necessary to include compositional variations in these regions. We sometimes found extremely low temperatures in cratons (Fig. 2) that are likely due to the secondary compositional component as well. Also in this case, the mineral physics uncertainties do not significantly affect the interpretation. Variations in the inferred geotherms obtained from different PREF starting models are generally small, and gradients with depth are similar. We found temperature variation as large as ≈75°C.

Fig. 2.

Seismic geotherms averaged over cratonic and oceanic regions. Selected cratonic regions are represented by solid areas (a), and respective geotherms are in solid lines (b). Oceanic regions are represented with dashed lines (a) and solid lines (b).


The capacity of inferring absolute temperatures and thermal gradients, however, depends on how precisely we can determine absolute seismic velocities and gradients. To this end, we should point out that our results rely on assumed crustal and anisotropy structure. Changes in the crust may affect the uppermost mantle. In addition, accounting for lateral variations in the crustal kernels may also modify the upper mantle structure down to 300 km (32). Even more important are the trade-offs between the isotropic and anisotropic parts of the model. The assessment of the trade-offs is not a simple task and is not addressed here thoroughly (see ref. 33). In general, because of the sensitivity of the data used, the trade-offs with anisotropy hamper a precise estimation of absolute isotropic velocities in the upper 300 km of the mantle. As a consequence, the thermal interpretation may change locally if a different anisotropic structure is assumed. For example, variations as large as hundreds of degrees and small (but dynamically important) changes in thermal gradients may result for the seismic geotherms shown in Fig. 2 when using another anisotropic model (we tested a previous version of the SAW642AN model; unpublished data).

In the transition zone, anisotropic effects have much less influence on the results. A purely thermal interpretation with our pyrolitic models requires unrealistic thermal profiles. This is discussed in a separate section.

The Average Model.

In Fig. 3, we show the inverted average 1D V S structure obtained from the 3D tomographic model. The average for an extremely different starting model also is shown. The other PREF models are similar to the one presented; hence, they have similar inverted structure (V S averages for different depth ranges are shown in Table 1). For a given anisotropy and crustal structure, absolute velocities and gradients are well constrained from long-period data, except in the vicinity of the 410 km discontinuity (compare the two inverted models in Fig. 3; see also Table 1). Trade-offs with anisotropy that are important for the upper mantle do not significantly affect the velocities and gradients below 300 km. We found that even a strongly anisotropic radial model, with variation from isotropy from +2% to −2% along the transition zone, does not greatly change the isotropic features of the inverted models. The regional velocity models beneath the same oceanic and cratonic regions as in Fig. 2 a are shown in thin lines in SI Fig. 5. The velocity structure at <300 km is very similar in these different tectonic regions, indicating the global validity of the inverted average model in the transition zone.

Fig. 3.

Inverted average structures for a PREF model (red) and a model with different gradient and smaller jump at 410 km (blue; model B in Table 1). Solid lines are the starting models; dashed lines are the inverted models. Circles on dashed lines indicate depths for which averages from the 3D model have been computed. Velocity models beneath the oceanic (thin dashed lines) and cratonic (thin solid lines) regions of Fig. 2 a are shown in the background in SI Fig. 5.


View this table:
Table 1.

Upper mantle V S averages (in m/s)


The inverted model tends toward PREM, which is confirmed to be a very good average seismic model, but there are striking differences due to the different starting parametrization. For example, because PREF does not have a 220-km discontinuity, as PREM does, the final model does not have this feature either. Instead, a velocity gradient higher than the starting PREF model is required by the seismic data around this depth (Fig. 3). In addition, long-period seismic waveforms are not directly sensitive to the mantle discontinuities. Therefore, the gradients near a discontinuity are not well constrained by those data, which can be clearly seen when we compare the gradients near the 410-km discontinuity between the PREF inverted model and model B with the small jump at 410 km (see Fig. 3).

Compared with the starting PREF model, the inverted model is characterized by a faster upper mantle, a slower transition zone, and different depth gradients throughout the upper mantle. We found high absolute shear velocities at ≈300 km. Around this depth, the velocity gradient is also much higher than the starting model (Fig. 3). The V S gradient between 500 and 600 km is slightly lower than PREM and more similar to the starting velocity gradient of the adiabatic pyrolite model (Fig. 3).

Overall, these features help to improve the fit to mode data for both toroidal and spheroidal branches (see SI Fig. 6). The PREF model(s) fit the fundamental modes satisfactorily, but the overtones, which are more sensitive to transition zone structure, are poorly fit (SI Fig. 6). This discrepancy supports the results of Cammarano et al. (22), who suggested the need for a slower transition zone for improving the travel-time data fit and extended it to the global scale. Note that both the starting PREF and the inverted model are isotropic. If the radial anisotropy parameters ξ and η (as defined in PREM) (19) are added assuming PREM values, the overall fit further improves (SI Fig. 6).

The Nature of the Transition Zone.

The seismic features obtained from our physically constrained inversion are different from the ones in the common seismic reference models. Absolute velocities and gradients with depth are well constrained throughout the upper mantle, except near discontinuities. Despite the large uncertainties of the seismic properties of mantle minerals (see ref. 22), our observations are robust enough to draw important conclusions on the nature of the upper mantle.

The velocity gradient with depth between 250 and 350 km (see Fig. 3) cannot be explained with any realistic thermal structure, assuming than composition does not change with depth. The velocity gradients of the PREF models within this depth-range are always similar (see SI Fig. 4), mainly because a relatively narrow range was chosen for the shear pressure derivatives of olivine (G′ was allowed to vary between 1.3 and 1.5; see table 1 of ref. 14). To reproduce this gradient with a 1,300°C adiabatic pyrolite, however, the values required of G′ are too high. If we vary only the olivine-G′ (all other values are from table 1 of ref. 14), a value of 2.5 is required. Thermal gradients are much less effective in reproducing this gradient. If we assume an isotherm (e.g., 1,700 K) instead of an adiabat, we still obtain a very high G′ (≈2.35). Although further experiments on the shear properties are warranted, such high values are very far from the experimental values available until now (e.g., refs. 34 and 35).

Another important aspect of the seismic features is the characteristic of the seismic jump at the 410-km discontinuity. The velocity jump imposed in the starting model at the olivine–wadsleyite transition dictates the gradients around it (Fig. 3). Pyrolitic models, which are ≈60 vol % olivine, have a larger jump than seismic reference models for this discontinuity [ΔV P, 6–10%; ΔV S, 6–12% (21); vs. 2.5% and 3.5%, respectively, for PREM and 3.5% and 4% for AK135]. If we assume in our starting model the large jump inferred from mineral physics, long-period seismic waveforms do require the changes in gradients just above and below the 410-km discontinuity shown in Fig. 3 (inverted PREF). A thermal interpretation of this structure is not feasible. By using our preferred pyrolitic models, we obtain temperatures of ≈1,250 K at 350 km and ≈1,750 K at ≈500 km and geodynamically unrealistic thermal gradients throughout the upper mantle.

Finally, the velocities in the transition zone are lower than for our preferred pyrolitic models. Note, however, that uncertainties in mineral physics are large for seismic properties of the transition zone (see Discussion).

We conclude that dry pyrolite cannot be reconciled with seismic data. This conclusion has a global character, as the similarity of the velocity profiles at <300 km beneath oceans and cratons show (SI Fig. 5).

Discussion

Uncertainties in Mineral Physics Data.

Although progress has been made in experimental mineral physics and theoretical computations of elastic properties at high pressure and temperature, there are still many uncertainties. In particular, shear properties of upper mantle minerals are difficult to measure at appropriate pressure and temperature ranges. For example, bulk and shear temperature-derivatives of the magnesium end-member of wadsleyite are still debated (∂K S/∂T and ∂G/∂T are −0.016 and −0.012 GPa/K, respectively, in ref. 37 and −0.012 and −0.017 GPa/K, respectively, in ref. 38). The result is that absolute seismic velocities of mineralogical models for a given composition and thermal structure are still very different (see refs. 14, 16, 39, and 40). However, velocity gradients with depth along a mantle adiabat are very similar in pyrolitic models because of the similar behavior of mantle minerals when they are compressed adiabatically. A third-order Birch–Murnagham equation of state with similar pressure derivatives (4–5 for K S and 1.5–2.5 for G) is usually appropriate to fit mineral physics data within the upper mantle pressure range (see, for instance, compilations in refs. 14, 16, and 3941). The elastic properties of olivine and wadsleyite at high pressure have been studied extensively and are known well enough to constrain the jump at the phase transition (see ref. 35). Note that a larger jump than PREM at the olivine–wadsleyite transition is consistently found in all of the available pyrolitic models. Considering all of the uncertainties of the elastic data, even the large uncertainties surrounding the shear pressure and temperature derivatives, the ΔV S jump at the phase transition has been found not <6% for pyrolite along a 1,300°C mantle adiabat (21). [Note that PREF models have a ΔV S jump closer to the average, ≈8.5% (see SI Fig. 4).] We found that a model with the same average V S and a ΔV S jump of 6%, with related slower transition zone and faster structure above it, is still not compatible with a purely thermal interpretation of the seismic data. The best way to further reduce the jump, as indicated by Duffy et al. (35), should be to have a much higher shear temperature-derivative of wadsleyite relative to olivine. Further studies are required to investigate this point, although the measured values on a large range of silicates (see ref. 35 and the available data on wadsleyite, e.g., refs. 8 and 38) suggest that such a case is very unlikely.

If Not Dry Pyrolite, Then What?

A piclogite composition, characterized by more garnet and pyroxenes and less olivine, has been proposed to provide a better fit to seismic 1D models in the transition zone (6). The original piclogite was only 16 vol % olivine (42). Based on later estimates of garnet elasticity at high pressure, a more olivine-rich mineral assemblage has been suggested (≈30% in ref. 43 and ≈40% in ref. 44). The piclogite model has been proposed in the frame of a mantle evolution that leads to a chemical differentiation of the shallow upper mantle from the transition zone (43). A chemical layering of the mantle is supported by several geochemical arguments. For example, studies of basalts worldwide indicate the need for different sources that have remained separate over a long time (>1 billion years) (45). The presence of such a chemical boundary implies layered mantle convection or an inefficient mixing mechanism. But seismic tomography provides clear evidence for subduction reaching far into the lower mantle (4), and modeling indicates that whole-mantle convection is sufficiently vigorous to mix large-scale heterogeneities (2). Several suggestions have been made to explain the discrepancy between geophysical and geochemical evidence [see review by van Keken et al. (46)]. Although the nature of mantle heterogeneities is still debated, there is a general consensus on excluding the 410-km discontinuity as a possible chemical boundary layer [with the exception of the water-filter hypothesis (47), which we discuss below in this section]. Furthermore, if the whole upper mantle, and not only the transition zone, is assumed to be piclogitic, the seismic velocities would be systematically lower than pyrolite. To reach similar average V S as pyrolite, the temperatures should be lower overall. Assuming adiabaticity, we estimate a 150 K cooler potential temperature for the upper mantle. This estimate is in conflict with the requirement of having a higher degree of partial melting than pyrolite to form basalts (48).

Small- and moderate-scale upper mantle compositional heterogeneities (102 to 105 m), which are beyond the resolution of the seismic data used here, have been inferred from geochemical studies (see ref. 49 and references therein). It has been suggested that the continuous reintroduction of heterogeneity in the upper mantle by subduction and partial melting may determine a nonequilibrium state of the mantle, the nature of which could be more complex than what is possible to infer from available modeling of mantle mixing (49). A mixture of enriched and depleted lithologies (e.g., eclogite plus harzburgite) instead of a unique mineral assemblage (e.g., pyrolite) may occur in the upper mantle. Such a hypothesis, however, does not accommodate the discrepancies we observe if the average chemical composition of the upper mantle is pyrolite. It has been proposed that small amounts of garnet pyroxenite in the upper mantle (≈5%) can explain the garnet signature in middle ocean ridge basalt (50). This percentage is a lower limit, because the melt derived from a peridotitic source with such proportions of pyroxenite will be largely dominated by the pyroxenite. A gradual enrichment in the olivine-free (eclogite or pyroxenite) component in the transition zone and right above it would help to reduce the discrepancies we observe. A garnet-rich composition would help increase seismic velocities and reduce the jump at the olivine–wadsleyite transition. For example, shear velocities for pyroxenite composition reach values of ≈4.8 km/s at a depth of 350 km (51), which is consistent with the absolute V S values required by long-period data (see Fig. 3). At this time, it is not clear whether such compositional gradation is dynamically possible. Note, for example, that the middle ocean ridge basalt (i.e., basalt) composition below ≈280 km is up to 5% denser than an olivine-rich lithology (>40 vol %), although it becomes lighter than that composition at ≈660 km (3). Recent thermochemical models (52), however, obtained such a compositional gradient. In terms of average composition, the hypothesis of a garnet-rich component in an olivine-rich matrix is similar to a piclogitic transition zone, but it does not prescribe the existence of a chemical boundary at 410 km, and it is more flexible for producing the observed geochemical heterogeneity (49).

Recently, it has been found that a large amount of water can be trapped in the nominally anhydrous minerals of the upper mantle. In particular, wadsleyite and ringwoodite can contain up to 3 wt % water (53) at pressures corresponding to transition-zone depths. This capacity tends to decrease as temperature increases, but it has been estimated that ≈1 wt % can be accommodated into the wadsleyite and ringwoodite structure at ≈1,400°C (54). Based on stability and melting relations of mantle phases (55), a hydrous transition zone that acts as a filter has been hypothesized by Bercovici and Karato (47) to reconcile geochemical and geophysical observations. The elastic properties of hydrous minerals are expected to change with the concentration of hydrogen (56), and, indeed, some of the recent experimental studies showed this (see, for example, refs. 36 and 57). Seismic velocities may be further reduced because of the dissipative effects of enhanced anelasticity (3). We note, however, that we found high velocities and gradients in the upper mantle above the 410-km discontinuity. Even considering the large uncertainties of elastic data, it is difficult to explain these features without invoking a compositional variation also above 410 km.

Conclusions

We show the feasibility of inverting seismic waveform data directly for physical parameters (temperature and composition) in the upper mantle. Lateral variations in temperature explain most of the data. Secondary compositional features, associated with cratons, emerge from the inversion as well. Absolute velocities in the shallow upper mantle depend on the crustal and anisotropic structure adopted. Hence, a detailed interpretation in terms of absolute temperatures and their gradients with depth rely on these features as well. The seismic geotherms we obtained are consistent, however, with expected geotherms in various tectonic settings. Including temperature-dependent anelasticity into the physical model is a key point for correctly interpreting the seismic data.

Below 300 km, average velocities and gradients with depth are well constrained by the long period data used, except near the mantle discontinuities. We found that seismic data globally require a slower transition zone and an overall faster shallow upper mantle than reference pyrolitic models. Despite the large uncertainties in mineral physics data, the observations are not compatible with a thermal interpretation assuming dry pyrolite. A gradual enrichment in a garnet-rich component throughout the upper mantle may help reconcile the observed discrepancies. A hydrated transition would help to lower the transition zone shear velocities, as required, but the high velocities and gradients above the wadsleyite stability field are not consistent with this hypothesis.

The interdisciplinary approach used here is general, and it can be applied to any seismic data set. Besides facilitating the interpretation of seismic data, it permits us to relate several parameters (that are not directly constrained by the seismic data used) to the same physical structure. For instance, the physical model is characterized by a given set of partial derivatives with respect to temperature (and composition) for seismic velocities and density. Even if the seismic data used are mostly sensitive only to shear velocity, the inverted 3D V S model can be directly related to 3D thermal structure (neglecting as a first guess the compositional signature), but also V P-, density-, and temperature-dependent attenuation structures can be predicted. By using this method, the same physical hypothesis can be tested against independent geophysical data. For example, dynamic topography or geoid anomalies can be used in principle to adjust the physical structure based on seismic data.

Other seismic data that complement the resolution of long-period waveform data can be included in the inversion. In particular, seismic data that are sensitive to mantle discontinuities, such as PP and SS precursors, and short-period reflected, refracted, and converted phases (as in ref. 58) would provide a better characterization of the seismic structure around the discontinuities. More investigation on shear properties at appropriate pressure and temperature conditions for both anhydrous and hydrous phases of the transition zone minerals would help to confirm or refine our conclusions.

In conclusion, we argue that a truly interdisciplinary approach that merges seismic observations, mineral physics data, and geodynamics constraints provides a valuable tool for obtaining a consistent physical model of the Earth's mantle.

Acknowledgments

We thank Federica Marone, Saskia Goes, Lars Stixrude, and Carolina Lithgow-Bertelloni for fruitful discussions; Federica Marone for constantly helping with the seismic codes; Lars Stixrude and Carolina Lithgow-Bertelloni for providing their mechanical–mixture model for pyrolite composition; and two anonymous reviewers for comments that improved the manuscript. This material is based on work partially supported by National Science Foundation Grants 0456335 and 0336951. This is contribution no. 07-04 of the Berkeley Seismological Laboratory, University of California.

Footnotes

  • *To whom correspondence should be addressed. E-mail: fabio{at}seismo.berkeley.edu
  • Author contributions: F.C. and B.R. designed research; F.C. and B.R. performed research; F.C. and B.R. analyzed data; and F.C. and B.R. wrote the paper.

  • The authors declare no conflict of interest.

  • This article is a PNAS Direct Submission.

  • This article contains supporting information online at www.pnas.org/cgi/content/full/0608075104/DC1.

  • Abbreviations:
    PREM,
    Preliminary Reference Earth Model;
    PREF,
    Physical Reference.

References