Subducted lithosphere beneath the Kuriles from migration of PP precursors

Earth and Planetary Science Letters 311 (2011) 101–111

Contents lists available at SciVerse ScienceDirect

Earth and Planetary Science Letters journal homepage: www.elsevier.com/locate/epsl

Subducted lithosphere beneath the Kuriles from migration of PP precursors Nicholas Schmerr a,⁎, Christine Thomas b a b

Department of Terrestrial Magnetism, Carnegie Institution of Washington, 5241 Broad Branch Rd., Washington, DC 20015 USA Institut für Geophysik, Westfälische Wilhelms-Universität, Corrensstraße 24, 48149 Münster, Germany

a r t i c l e

i n f o

Article history: Received 11 May 2011 Received in revised form 22 August 2011 Accepted 2 September 2011 Available online 4 October 2011 Editor: P. Shearer Keywords: Kurile Subduction Zone upper mantle discontinuities PP precursors P660P transition zone

a b s t r a c t We seismically image both thermal and chemical heterogeneity of the mantle beneath the Kurile subduction zone using P-wave energy reflected from the underside of discontinuities, arriving as precursors to the seismic phase PP. We take advantage of new broadband seismic data provided by the High Lava Plains Seismic Experiment and EarthScope's USArray, collecting a dataset of 31 high-quality Sumatran earthquakes sampling beneath the Kuriles. We employ high-resolution array analysis techniques, including migration and vespagrams, to identify precursory arrivals and study lateral variations in discontinuity depth, sharpness, and impedance of the mantle transition zone. We find the 410 km boundary is at 395 km near the subducting Kurile slab, though the boundary is 410–425 km deep elsewhere. In regions away from subduction, we do not detect a laterally continuous underside reflection of P-waves from the 660 km discontinuity. However, in the vicinity of the subducting Kurile slab, we detect robust P660P reflections from interfaces near 620–670 km depth, signifying an increase in the impedance contrast at 660 km depth. We also detect deeper reflectors, down to 720 km depth, beneath the Kurile slab in a localized area. Cold, aluminum-depleted harzburgitic lithosphere residing at the base of the transition zone best explains the local enhancement of the 660 km discontinuity P-wave impedance contrast. Our new discontinuity measurements support the hypothesis of cold, depleted lithosphere stagnating at the 660 km discontinuity beneath the Kuriles subduction zone, and imply the 660 km boundary can locally impede mantle flow and produce chemical heterogeneity within the transition zone. © 2011 Elsevier B.V. All rights reserved.

1. Introduction The Kurile subduction zone of the northwestern Pacific has been a key region for investigating the nature of convection in the mantle and establishing the dynamics of subduction zones. The Kurile island arc extends from the tip of the Kamchatka peninsula and terminates in a complex interaction with the Japanese island arc at Hokkaido Island. The arc is formed by the oblique convergence of relatively old (N100 Ma) (Müller et al., 1997) Pacific lithosphere subducting beneath the Sea of Okhotsk at ~8.2 cm/yr (DeMets, 1992). Seismicity associated with the descending slab beneath the Kurile Islands extends into the mantle transition zone (MTZ) to a depth of ~670 km, and numerous tomography and travel time studies have found evidence for high velocity materials at depth (e.g., Jordan, 1976; van der Hilst et al., 1993). Much of the early work in the Kuriles focused on determining the fate of the subducted lithosphere utilizing deep seismicity associated with the slab (Fedotov, 1965). Measurements of near source S, P, and depth phase travel time anomalies indicated the presence of a high velocity (N5%) anomaly in the vicinity of the Benioff zone (Fukao, 1977). Travel time sphere residual analyses by (Jordan, 1977) found ⁎ Corresponding author at: NASA Goddard Space Flight Center, Planetary Geodynamics Laboratory, Code 698, Greenbelt, MD 20771, USA. E-mail addresses: [email protected] (N. Schmerr), [email protected] (C. Thomas). 0012-821X/$ – see front matter © 2011 Elsevier B.V. All rights reserved. doi:10.1016/j.epsl.2011.09.002

that the Kurile slab penetrated into the lower mantle to at least 900 km depth, though later studies using the same technique indicated that the high velocities of the slab did not extend into the lower mantle (Fischer et al., 1988; Gaherty et al., 1991; Suetsugu, 1989). Early 3-dimensional tomography models of shear wave structure beneath the Kuriles found evidence for high velocities located just above the 660 km discontinuity, leading to the suggestion that the slab stalls or stagnates at the 660 km discontinuity (Fukao et al., 1992; Fukao et al., 2001; Fukao et al., 2009). A number of detailed tomography models for P and S waves (e.g., Ding and Grand, 1994; Li et al., 2008; Megnin and Romanowicz, 2000; Miller and Kennett, 2006; Ritsema et al., 2011; van der Hilst et al., 1993) image a high velocity anomaly above the 660 km discontinuity to the south near Hokkaido, stagnating on top the 660 km discontinuity in the down-dip direction, while to the north near Kamchatka this anomaly penetrates directly into the lower mantle. Tomographic studies find shear wave velocity anomalies in the slab dominate over bulk sound velocities (Gorbatov and Kennett, 2003), suggesting that the imaged velocity heterogeneity is produced by the temperature anomaly of the subducting lithosphere over the compositional effects of the slab, though uncertainties in this relationship remain unresolved. Another way to probe the structure beneath the Kuriles is through seismic imaging of the depth of the upper mantle discontinuities (e.g., Shearer, 2000). The two major upper mantle seismic discontinuities that delineate the MTZ arise from solid-to-solid phase changes in

102

N. Schmerr, C. Thomas / Earth and Planetary Science Letters 311 (2011) 101–111

the mineral olivine that are sensitive to both mantle temperature and composition (Bina and Helffrich, 1994). The exothermic transformation of olivine (ol) to wadsleyite (wd) generates the 410 km discontinuity (410) (Katsura and Ito, 1989), and the endothermic dissociation of ringwoodite (rw) into Mg-Perovskite (pv) and magnesiowüstite (mw) forms the 660 km discontinuity (660) (Ito and Takahashi, 1989). Upper mantle discontinuity depth is perturbed by lateral changes in mantle thermal and chemical properties, such as in the vicinity of cold subducted lithosphere in the MTZ. This owes to the depth of each discontinuity being controlled by the Clapeyron slope of the phase transitions, or the change in pressure of a phase transition for a given change in temperature. The Clapeyron slope of the 410 is positive, and is negative at the 660, resulting in a thickening of the MTZ in the presence of a vertically continuous cold mantle anomaly, and thinning of the MTZ in the presence of vertically continuous hot mantle anomaly (Lebedev et al., 2003). Lateral variations in mantle composition, such as the presence of hydrogen within a subducted slab, can also alter the depth of the discontinuities and produce topography on the boundaries (Karato, 2006). Several other seismic discontinuities have also been detected at MTZ depths. The transformation of wadsleyite to ringwoodite gives rise to a smaller discontinuity near 520 km depth (Shearer, 1990). In some regions, a second discontinuity near 560 km depth related to the transformation of garnet to Ca-Perovskite is observed (Deuss and Woodhouse, 2001). In the vicinity of subduction zones, a low velocity zone is often observed above the 410 km discontinuity, related to the presence of partial melt above this boundary (Song et al., 2004). The formation of ilmenite in colder regions of the mantle and the completion of the garnet phase transitions at 750 km depth can produce additional discontinuities near 660 depth (Weidner and Wang, 1998). Owing to the large variety of seismic sources and stations situated around the Kurile subduction zone, there have been numerous investigations utilizing reflected, refracted, and converted seismic waves interacting with the MTZ discontinuities. In the southern Kuriles, triplicated seismic waves detect a high velocity transition zone below 500 km depth, and a deep 690 km discontinuity, consistent with a stagnant slab in the MTZ (Tajima and Grand, 1995). The northwestern Pacific is extremely well-sampled by long-period underside reflected shear waves, which occur as precursory energy to the seismic phase SS and are sensitive to structure half way between the source and receiver (Flanagan and Shearer, 1998; Gu et al., 1998; Shearer and Masters, 1992). These seismic phases are 1–10% of the arriving SS amplitude, and require the stacking of datasets of hundreds to thousands of records to be observed robustly (Shearer, 1990). The SS precursory phases detect a large-scale depression of the 660 km discontinuity beneath the Kuriles that is associated with the regions of higher velocity and inferred slab stagnancy in the tomographic models (Cao et al., 2010; Shearer, 1991; Shearer and Masters, 1992). A comparable approach to studying discontinuity structure can be made with underside reflected P wave energy arriving as a precursor to the seismic phase PP (King et al., 1975). Studies of the PP precursors have largely focused on the underside reflection from the 410 (P410P), as the underside reflection from the 660 (P660P) is difficult to detect owing to extremely low amplitudes (Estabrook and Kind, 1996; Neele and Snieder, 1992). In several recent studies, P660P has been detected in different regions (Deuss et al., 2006; Thomas and Billen, 2009), suggesting that local heterogeneities can enhance the amplitude of this arrival. Past array studies of the Kuriles with short-period PP precursors find an elevated 410 km discontinuity (Rost and Weber, 2002), but did not detect the 660 boundary. Here we take advantage of the deployment of the broadband seismometers in the EarthScope USArray 1 and the recent Program for Array Seismic Studies of the Continental Lithosphere (PASSCAL)

1

http://www.earthscope.org

High Lava Plains Seismic Experiment (HLP) (Carlson et al., 2005) to investigate detailed discontinuity structure beneath the Kuriles with PP precursors. We utilize an array migration technique that allows us to isolate the arrivals of underside reflected P energy from the noise and other interfering seismic phases, and focus in particular on mapping the character of the 660 km discontinuity. Our seismic sources allow us to investigate structure along a profile across the southern Kuriles and further test the hypothesis of a cold, stagnant slab in the lower MTZ with a method that complements recent tomography results from this region.

2. Dataset To study discontinuity structure, we collect a dataset of broadband velocity seismograms that sample beneath the Kurile subduction zone. Data were downloaded from the Incorporated Institutions for Seismology Data Management Center. 2 The sources used originated in the Sumatran subduction zone recorded between January 2005 to September 2009 by the EarthScope USArray and the HLP Seismic Experiment (Fig. 1, Table 1). The HLP stations formed the backbone of our array and were supplemented by USArray stations deployed within 750 km from the center of the HLP. Events were required to have a Mw ≥ 5.8 to ensure sufficient seismic energy for our analysis. The instrument response was deconvolved from each seismogram and we studied the underside reflections of PP on the vertical component. The seismograms of each event were visually inspected to determine the quality of the body wave arrivals; high quality events did not require any filtering to identify the PP arrivals from the background noise. In addition to these high quality events, we kept earthquakes that possessed clear PP arrivals when low-pass filtered at a corner of 10 s. The dataset initially consisted of 169 events and 14,047 seismograms, after visual inspection, the dataset consisted of 31 high-quality events with 2180 seismograms, with an average of 116 records per event (minimum 22, maximum 231). For our analysis, we marked a reference time for the PP arrivals using the travel time of the highest amplitude arrival of the first swing of the PP waveform. We measured a signal to noise ratio for every seismogram by comparing the enveloped maximum amplitude in a 50 second wide window centered on the pick of the PP phase to the enveloped maximum amplitude in a 100 second window ending 50 s before the predicted arrival time for the Pdiff phase. To investigate frequency dependence, we generated several subsets of the data with band-pass filters with upper and lower corners at 1–10 s, 3–25 s, 6–50 s, 10– 100 s, and 15–100 s.

3. Method To study the upper mantle discontinuities beneath the Kuriles, we utilize several array methods to amplify the PP precursory arrivals out of the noise levels: velocity spectral analysis (vespagram), slowness backazimuth analysis, and migration. To test the accuracy of our methodology, we also generate reflectivity synthetic seismograms (Fuchs and Müller, 1971) for each event depth and source focal mechanism and analyze the synthetics alongside the data. We start by generating 4th root vespagrams (Davies et al., 1971; Muirhead and Datt, 1976; Rost and Thomas, 2002) for each event to verify that the arrival times and slowness of the precursors matched predictions from ak135 (Kennett et al., 1995). We only further analyze events with a signal to noise ratio ≥ 3.0 for both PP and P410P in the vespagrams (Fig. 2). The vespagrams allow us to determine whether the precursory energy is well separated from the interfering energy of PKiKP, topside reflections, and other non-underside reflection seismic phases. 2

http://www.iris.edu

N. Schmerr, C. Thomas / Earth and Planetary Science Letters 311 (2011) 101–111

a)

103

b)

15

45

5

40

-5

235 90

100

240

245

110

c) 60

Siberia Kamchatka Sea of Okhotsk

d)

50

PP

41 66 0 0

Ku

rile

s

P410P

CMB Hokkaido 40

ICB

Japan

140

150

Fig. 1. Seismic dataset parameters and geography of the Kurile subduction zone. a) Sumatran subduction zone source region with events used in this study (stars). b) Station geometry of the HLP seismic array (black triangles) and location of the EarthScope Transportable Array stations used in this study (inverted triangles). c) Kurile subduction zone study region and location of underside reflections (circles). d) Example raypaths of PP and P410P. Underside reflections occur at the center of the seismic raypath. CMB—core mantle boundary, ICB—inner core boundary. Slab depth contours are from (Gudmundsson and Sambridge, 1998).

When the arrivals were not separated well, the event was not considered further. We retrieve the depths of the upper mantle discontinuities using a migration method adapted from (Thomas and Billen, 2009). The migration method allows the back projection of precursory energy to the underside reflection point and lessens the size of the Fresnel zone thereby enhancing resolution (e.g., Rost and Thomas, 2009). To migrate, we generate a 40° by 40° grid in 1° increments centered on the reflection point halfway between the average array center and source location. The grid is incremented every 5 km between 0 and 900 km depth. We then calculate travel times between each grid point and array station, as well as between each grid point and the source location using a raytracing code through ak135. Each seismic trace is shifted by the resulting delay time and stacked. We select the maximum amplitude in a 3 second window centered on the theoretical arrival time of the

precursor, and utilize a bootstrap resampling algorithm with 300 random replacement resamples to evaluate the 95% confidence amplitude in each stack (Efron and Tibshirani, 1986). The resulting migrated energy recreates the X-like shape of the PP Fresnel zone, owing to the minimum–maximum travel time of underside reflected waves (Choy and Richards, 1975). The stacked energy from the migrated section will spread along one isochrone if only one source and receiver is used. However, using a seismic array, the spreading of energy will be less and the largest amplitude will be focused at the reflection point, enhancing resolution. The backazimuth range for the data used in this study is on the order of 7°, which also lessens the size of the effective Fresnel zone. Past work has demonstrated array enhancement of resolution, albeit using a smaller array than the one used here (Rost and Thomas, 2009). Fig. 3 includes examples of synthetics where focusing of

N. Schmerr, C. Thomas / Earth and Planetary Science Letters 311 (2011) 101–111

Table 1 Event parameters and associated reflector depths. Lat

Lon

Depth Bouncepoint Lat

Reflectors

Lon

d410 d660 dX1 km

dd/mm/yyyy hh:mm

degN

degE

km

degN

degE

24/07/2005 25/04/2006 27/06/2006 17/07/2006 27/07/2006 11/08/2006 01/12/2006 07/04/2007 08/08/2007 13/09/2007 13/09/2007 14/09/2007 19/09/2007 20/09/2007 24/10/2007 22/01/2008 20/02/2008 25/02/2008 25/02/2008 25/02/2008 29/03/2008 27/06/2008 08/08/2008 22/11/2008 30/12/2008 15/04/2009 15/04/2009 13/07/2009 16/08/2009 18/08/2009 30/09/2009

7.92 1.99 6.50 − 9.42 1.71 2.40 3.39 2.92 − 5.93 − 2.13 − 3.17 − 4.07 − 2.75 − 2.00 − 3.40 1.01 2.77 − 2.49 − 2.33 − 2.24 2.86 11.01 − 3.92 − 4.35 − 4.30 − 3.09 − 3.12 − 9.13 − 1.48 − 0.91 − 0.72

92.19 16 97.00 21 92.79 28 108.32 21 97.15 20 96.35 22 99.08 204 95.70 30 107.68 291 99.63 22 101.52 53 101.17 23 100.89 35 100.14 30 101.02 21 97.44 20 95.96 26 99.97 25 99.89 25 99.81 25 95.30 20 91.82 17 101.08 35 101.26 24 101.22 20 100.41 6 100.47 22 119.32 67 99.49 20 97.95 10 99.87 81

56.18 49.97 55.74 35.21 49.84 50.93 50.27 52.09 38.71 46.10 44.24 43.96 44.89 45.87 44.06 49.86 52.03 45.83 46.04 46.10 52.51 59.92 44.05 45.23 45.14 46.95 46.75 32.73 46.56 47.90 47.17

145.43 425 – 143.72 415 700 139.07 410 – 156.05 410 – 143.78 420 – 142.75 420 630d 145.51 410 – 141.53 415 640d 154.52 400 – 145.99 440n – 148.09 410n – 147.57 405 – 147.43 – 660 146.59 410 680 147.66 405 690 143.47 420 630d 141.80 410 710 146.43 410 670 146.35 400 675 146.23 395 670 140.93 420 620d 137.66 420 – 147.82 410 690 147.21 420 – 147.42 410 – 145.79 425 – 145.96 415 164.66 400 – 145.99 420 – 144.26 410 650 145.39 405 –

15:42 18:26 18:07 15:45 11:16 20:54 03:58 09:51 17:04 03:35 16:09 06:01 07:27 08:31 02:50 17:14 08:08 08:36 18:06 21:02 17:30 11:40 06:37 16:01 19:49 17:47 20:01 10:52 07:38 17:50 10:16

km

km 740

750

Epicentral Distance (deg)

Event

a)

PP

135

Pdiff

p410pPdiff

P410P PKiKP P660P PKPp410p

p660pPdiff

PKPp660p

130

125

120

775

−300

−200

−100

0

Time (s)

b)

750

730 830

n: neg. polarity; –: reflector not observed; d: depth phase.

stacked energy at the saddle point of the precursor reflection is observed. A problem, however, could be a precursor phase that has traveled out of plane. In these cases we measure the deviation from the theoretical reflection point using slowness-back azimuth analysis (Rost and Thomas, 2002; Rost and Weber, 2002) and if the reflection point is within the Fresnel zone of a wave that would travel along the great circle path, we utilized the measurement (Fig. 4). Events with poor focusing of the precursors in the migration and/or slowness backazimuth analysis were discarded. To select the depth of the discontinuity, we take an amplitude profile at the calculated reflection point for the source and array center (Fig. 3). This migration technique removes the effects of seismic energy at other slownesses and enhances the arrivals from the precursors (Thomas and Billen, 2009). 4. Results We compare the information from the vespagrams and migration profiles to determine the depth and amplitude of PP precursory arrivals from the 410 and 660 km discontinuities. In the vespagrams, the PP precursors arrive in a slowness window of 6.2–6.9 s/deg, though often the actual arrivals are spread out over slightly larger slowness ranges (Fig. 2). Examples of well-defined precursory arrivals in our migrations and vespagrams are shown for the events in Figs. 2–5. We investigate discontinuity structure at several band-pass filters of the dataset to test for frequency-dependent effects. Higher frequencies are contaminated by more noise, but often still contain useful information about short-scale discontinuity structure. We use the bootstrap-resampling algorithm to determine whether migrated amplitudes at higher frequencies are statistically greater than the

Slowness (s/deg)

104

10

P660P

p410pPdiff

5 Pdiff

P410P

PP

p660pPdiff

PKPp410p

PKPp660p

PKiKP

0

−300

−200

−100

0

Time (s) Fig. 2. Waveforms and vespagram from earthquake occurring on 15/04/2009-20:01 for a robust P410P and non-detection of P660P. a) Seismograms recorded at the HLP and TA arrays, aligned on the maximum amplitude of the PP arrival. The data are band-pass filtered with corners at 10 and 100 s. Dotted lines show the predicted travel times for select seismic phases. b) 4th root vespagram of the event, positive amplitudes are black, and negative amplitudes are white. Crosses indicate the predicted slowness and travel time for the labeled seismic phases. In both panels, theoretical travel times and slownesses for seismic phases are from ak135.

background noise in the migration. A standard practice in the analysis of underside reflections is to use amplitudes that fall above the 2σ confidence interval (e.g., Flanagan and Shearer, 1998; Gu et al., 1998). The validity of this approach is demonstrated in Fig. 5, where results for different filters of the data are shown. The time domain corners of the band-pass filter are given in the top right of each panel. Though more noise is present in the higher frequency results, the bootstrap algorithm reveals two significant peaks at 400 and 675 km depth in all three band-pass filters. We obtain consistent results for data with a high frequency corner set at 0.1 Hz (10 s), though several events produce clear arrivals in the 6–50 second filter. The filters with corners at 1–10 s and 3–25 s are too contaminated by noise for further interpretation. The frequency analysis reveals several statistically robust peaks and troughs; these are related to the waveform of the seismic arrivals, such as the depth phase of P410P and negative sidelobe of PP (Fig. 5). Synthetic seismograms generated for each set of event parameters and station geometries facilitate the identification of these other arrivals. We only further analyze events with at least a strong P410P in the vespagram, clear focusing in the slowness backazimuth analysis, and a well-defined arrival in the migration. The resulting 31 events are detailed in Table 1 and form a nearly linear profile across the southern Kurile subduction zone (Fig. 1) with a kink at 40°N mirroring the bend of the Sumatran subduction zone source region. The depth of each discontinuity is determined by finding the maximum amplitude near the predicted depth for each boundary in the migrated sections. The accuracy of the retrieved discontinuity depth is controlled by the fineness of the migrated grid, though the resolution of small-scale features is limited by the migration of long wavelength precursory seismic energy. We found that a grid spacing of 5 km allowed us to determine the peak amplitude

N. Schmerr, C. Thomas / Earth and Planetary Science Letters 311 (2011) 101–111

a) 0

105

b)

ak135

data

c)

ak135

0 km

data

Depth (km)

250

400 km 400

400

675

675

500

675 km

750

0

10

20

30

40 0

10

20

30

40

0

10

20

Distance (deg)

30

40 0

10

20

30

40

Distance (deg)

Fig. 3. Migration results for the event 25/02/2008-18:06 and ak135 reflectivity synthetics, with a prominent P660P detection. a) Raypath parallel cross sections through the migrated data and synthetic amplitudes. Positive amplitudes are blue, negative are red. The orientation of the cross sections and lateral extent of the migration grid is given by the small globe. At the center of each panel, we show the migration amplitudes at the predicted reflection point. The black and gray shading are respectively the positive and negative 95% confident amplitudes from the bootstrap resampling. b) Raypath perpendicular cross sections; the details are the same as in part (a). c) Horizontal slices through the migrated energy of the event. In all panels, the maximum amplitude is normalized to one to enhance the visibility of the migrated amplitudes.

for each precursor and provides a reasonable tradeoff between the approximate sensitivity of long period P-waves (Chaljub and Tarantola, 1997) and efficient computation of the migrated grids. The cross section in Fig. 6 presents migration results demonstrating the lateral variation in the discontinuity depth beneath the Kuriles. Owing to the source and array geometry of our dataset, the majority of our observations are densely clustered in a region spanning 1000 km from 43°N to 53°N. Sampling further to the south and north of the subduction zone is sparser but still provides information about the depths of the upper mantle discontinuities every 500 km. The depths of the upper mantle discontinuities must be corrected for variations in crustal thickness and lateral mantle heterogeneity as both affect the resulting depths. We use 1-D raytracing through CRUST 2.0 (Bassin et al., 2000) and MITP08 (Li et al., 2008) to compute correction values for the entire path of PP and the precursors. The resulting corrections are relatively small (b2 s) compared to the travel time heterogeneity observed in the depths of the discontinuities. In the vespagrams and migrated profiles, we observe a clear reflection from the 410 km discontinuity in 30 of the studied events (Table 1). Beneath the Kuriles, the 410 km discontinuity occurs at a depth of 410 ± 15 km, with the shallowest measurements in the vicinity of the subducting Pacific lithosphere (Fig. 6). A total of 15 events produce reflections from the 660 km discontinuity, and with

Backazimuth (deg)

t= -315 s

Amax =0.029

t= -135 s

Amax =0.004

t= -125 s

the majority of 660 detections in the vicinity of the subducting Pacific lithosphere. Topography on the 660 km discontinuity is complex, and in several regions, we detect multiple reflectors in the depth range of 620–830 km. We do not detect the 660 km discontinuity in events sampling beneath the Pacific Ocean. However, from 43°N to 48°N, corresponding to the mantle lithosphere directly beneath the subducted slab, we detect a laterally continuous reflection at 660– 700 km depth (Fig. 6), and a second intermittent and weak reflector at N800 km depth. In the MTZ region from 48°N to 53°N, associated with elevated P-wave velocities in MITP08 (Li et al., 2008), multiple reflectors become apparent, one at 620–650 km depth, and much deeper reflector ranging from 700 to 800 km depth. Further to the north (towards Siberia), where the slab is no longer present in the MITP08 tomography, we detect no reflection from the 660 km discontinuity. We also measure the amplitude ratios of the precursory energy relative to the reference PP phase from the vespagrams (Fig. 7) where precursors are well separated from interfering phases such as PKiKP. The amplitude ratios provide information on the sharpness and strength of the impedance contrast across each discontinuity (Chambers et al., 2005; Shearer and Flanagan, 1999). Amplitude ratios must also be corrected for the effects of geometrical spreading of the seismic wavefront, source radiation patterns, and anelasticity Amax =0.024

Amax =0.009

t= -80 s

Amax =0.224

t= 0 s

Amax

360

P660P

Pdiff

PKiKP

PP

P410P

300

4

8

4

8

4

8

4

8

4

8

0

Slowness (s/deg) Fig. 4. Slowness backazimuth analysis of event 25/02/2008-18:06 for select seismic phases, including P660P. Amplitudes are normalized to one; the reference amplitude (Amax) is shown in the upper right of each panel. The theoretical slowness and backazimuth of each seismic phase is marked with a white cross; relative travel times and slownesses are calculated from ak135.

106

N. Schmerr, C. Thomas / Earth and Planetary Science Letters 311 (2011) 101–111

a)

15-100s

10

b) PP

P660P

P410P

PP P675P

P400P pP400P

5

900

700

600

500

400

300

200

100

0

300

200

100

0

300

200

100

0

Depth (km)

10-100s

10

Slowness (s/deg)

800

5

900

800

700

600

6-50s

10

500

400

Depth (km)

5

900 −300

−200

−100

0

800

700

600

500

400

Depth (km)

Time(s) Fig. 5. Vespagrams and migration amplitudes for event 25/02/2008-18:06 with a prominent P660P detection. a) Vespagrams for three band-pass filters of the dataset. The band-pass filter used is shown in the upper right corner of each panel. Positive amplitudes are black, negative are white. b) Migration profiles of energy at the band-pass filters of the data in (a). Arrows indicate the transition zone discontinuities. The profiles are located at the predicted underside reflection point for the center of the array and source location. Black and gray shading are respectively the positive and negative 95% confident amplitudes from the bootstrap resampling.

effects. We correct for these effects by multiplying each measured amplitude ratio by a correction factor computed from ak135 reflectivity synthetics, using the PREM (Dziewonski and Anderson, 1981) values of attenuation for the mantle. We correct our amplitude ratios to a reference distance of 125°; the correction factor is then the dividend between the synthetic PdP/PP ratio at 125° and at the epicentral distance of the source and array. The resulting scaling factor is used to correct the measured amplitude ratio from the dataset for the effects of differential attenuation, source radiation, and geometrical spreading. We do not correct for lateral variations in attenuation structure in the mantle; large lateral variations in attenuation of the PP waveform in the upper mantle would produce correlated amplitude ratios between the 410 and 660, which we do not observe. The average amplitude ratio for P410P/PP is 0.0375 ± 0.018, close to the predicted value of 0.0365 for the reflection coefficient at a reference epicentral distance of 125° in ak135. The average amplitude ratio for observations of P660P/PP is 0.0246 ± 0.021, nearly half that of the predicted value of 0.0452 from ak135. However, in the region underlying the slab, the P660P/PP amplitude ratios become comparable to those produced by ak135. The deepest reflectors (N700 km) have amplitude ratios that are below 0.02 and the reflections are weaker to the north. 5. Discussion Past investigations of discontinuity structure utilizing underside reflected S-waves in the Kuriles region detect a transition zone 250– 260 km thick (Flanagan and Shearer, 1998; Gu et al., 1998; Shearer,

1991; Shearer and Masters, 1992). This thickening is the result of a major depression on the 660 km discontinuity throughout the region, consistent with the presence of cold, subducted lithosphere in the MTZ. Our results support this hypothesis and allow us to further characterize the properties of the subducted lithosphere in the MTZ. Most notable is the presence of P660P, a seismic phase absent from global investigations of upper mantle structure. Long-period stacks of underside reflected body waves reveal well-developed P410P and S410S arrivals, however, only the S660S is routinely observed and P660P phase is missing in most analyses (Deuss, 2009). Only recently have several studies revealed robust intermittent detections of P660P (Deuss, 2009; Deuss et al., 2006; Thomas and Billen, 2009). Several past investigations of P660P amplitude ratios explain the difficulty in detecting this seismic phase, finding that the density of velocity jumps at 660 km depth are about 50% of their PREM values (Estabrook and Kind, 1996; Shearer and Flanagan, 1999). Our results indicate a region with enhanced reflectivity at the 660 km discontinuity. The reflection coefficient for an upward traveling seismic wave is dependent upon the incidence angle and impedance contrast across the discontinuity. A high impedance contrast across a discontinuity will increase the reflection coefficient. However, the reflection coefficient is also dependent on the incidence angle of an incoming wave. We define incidence angle as the angle between the upgoing raypath and a vertical vector pointing towards the center of the Earth. For a given impedance contrast and increasing epicentral distance, the incidence angle of an underside reflection increases, approaching the horizontal for deep reflectors; for P-waves this will lower the reflection

N. Schmerr, C. Thomas / Earth and Planetary Science Letters 311 (2011) 101–111

107

a) A

A

China

0

Sea of Okhotsk

Kuriles

50

45

A’

Pacific Ocean

A’ Depth (km)

200

400

pP660P

600

800 65

60

55

40

35

30

Latitude ( )

b) 0

0

0

0

0

100

100

100

100

100

200

200

200

200

200

300

300

300

300

500

400

500

400

500

600

600

600

600

700

700

700

700

700

800

800

800

800

800

900 27/06/2006-18:07

900 29/03/2008-17:30

900 18/08/2009-17:50

900 25/02/2008-18:06

600

Depth (km)

500

400

Depth (km)

400

Depth (km)

500

Depth (km)

400

Depth (km)

pP410P

300

pP660P

900 17/07/2006-15:45

Fig. 6. Detailed upper mantle discontinuity structure beneath the Kurile subduction zone and surrounding regions. a) Cross section through MITP08 tomography (Li et al., 2008) showing the discontinuity depths measured from the migration of each event. Triangles indicate the 410 km discontinuity, inverted triangles the 660 km discontinuity, and stars a second reflector near 650–800 km depth. Green circles are earthquake hypocenters within 100 km of the cross section. The dotted vertical line indicates the location of the kink in A–A′. b) Migration results for selected events. Details of the migration profiles are similar to Figs. 3 and 5. The horizontal arrows indicate robust precursory arrivals in each event.

coefficient owing to more energy converting to S-waves at the interface. The opposite behavior is present in S-waves, as the incidence angle increases, the amplitude of the underside reflection will also increase. Thus deep reflectors with weak impedance contrasts are difficult to detect with PP precursors—the impedance contrast must be sufficiently high to offset the effects of the increasing incidence angle. Constraining the impedance contrast requires knowledge of the sharpness of the discontinuity, a parameter directly related to the shape of the phase loop of the associated mineralogical phase change (Stixrude, 1997). The presence of a broad gradient diffuses seismic energy and reduces the amplitude of the underside reflections. High frequency investigations of the sharpness of the 660 km discontinuity around the globe with precursors to P′P′ indicate the width of the post-spinel phase transition is extremely sharp (b5 km) (Benz and Vidale, 1993; Xu et al., 2003). The sensitivity of precursory amplitudes to gradient width grows at longer periods, thus large changes in the width of the gradient will produce frequency dependent changes in the precursory amplitudes. To test the effect of gradient width, we generated reflectivity synthetic seismograms for a 10-

second dominant period P660P using seismic models with increasing gradient widths at 660 km depth. These tests indicate a gradient width N 35 km is needed to reduce P660P amplitudes by more than 15%, a value inconsistent with the P′P′ observations. We further investigate frequency dependence in our dataset using 5 different band pass filter windows with upper and lower corners at 1–10 s, 3–25 s, 6–50 s, 10–100 s, and 15–100 s. The amplitudes of P660P and P410P in the 1–10 and 3–25 second windows were unconstrained owing to large amounts of noise present in the seismograms. However, there was no systematic variation in precursory amplitudes at the longer period band pass filters. Thus, the presence of large gradients would reduce the observed P660P amplitudes making them more difficult to detect, a result inconsistent with the localized enhancement of P660P amplitudes observed beneath the Kuriles, and the observations of a globally sharp interface imaged by P′P′ studies. (Xu et al., 2003) Several mechanisms are potentially responsible for locally enhancing the impedance contrast across the 660 km discontinuity beneath the Kuriles subduction zone. Seismic tomography reveals

108

N. Schmerr, C. Thomas / Earth and Planetary Science Letters 311 (2011) 101–111

Siberia

Sea of Okhotsk

Kuriles

410 Depth (km)

300

400 PdP/PP Amplitude 500 0.05 0.04

660 Depth (km)

600

0.03 pP660P

0.02 0.01

700

800 55

50

45

Latitude ( ) Fig. 7. Detailed discontinuity amplitude ratios and depths for robust P410P and P660P detections. The cross section corresponds to the densely sampled region beneath the Kuriles subduction zone in Fig. 6. Global average discontinuity depths are shown as dotted lines, the circled values are P660P depth phase detections.

a)

A

Pyrolitic Mantle (~5 wt% Al2O3)

vp

MORB (>15 wt% Al2O3)

Harzburgite (~1 wt% Al2O3)

A’

Pyrolitic Mantle (~5 wt% Al2O3)

ilm Mg-pv

660

Al-content (wt %)

WW98-1700K-5%Al

WW98-1700K-3%Al

WW98-1900K-3%Al

WW98-1900K-5%Al

A

A’

15 12 9 6 3 0

2000 1500 1000

A

Temperature (K)

c)

rw Mg-pv+mw

WW98-1900K-5%Al

b)

630

A’

0.05

PdP/PP

ak135

WW98 0.00

P660P Detections Fig. 8. Interpreted thermal and chemical structure beneath the Kuriles subduction zone. a) Predicted velocity structure for a slab stagnating in the transition zone, (Weidner and Wang, 1998) and 660 km discontinuity topography from variations in thermal and chemical structure at the boundary. b) Hypothetical variations in Al-content (gray) and temperature (black) at 660 km depth. c) Reflectivity synthetic amplitude ratios (solid black line) for P660P/PP, and P660P/PP amplitude ratio for the ak135 model (dotted line). These are compared to the P660P/ PP amplitude ratios measured from the data sampling beneath the Kuriles, arranged from left to right by decreasing latitude (gray circles).

N. Schmerr, C. Thomas / Earth and Planetary Science Letters 311 (2011) 101–111

high velocity anomalies present in the MTZ beneath the Kuriles region coincident with our detections of P660P (e.g., Grand, 1994; Li et al., 2008; Miller et al., 2006; van der Hilst, 1995; Widiyantoro and van der Hilst, 1997). This high velocity anomaly extends across the 660 km discontinuity into the lower mantle over a region approximately 500 km wide. We show a schematic interpretation of this structure in Fig. 8. Seismic tomography, travel time measurements, past investigations of discontinuity structure, and our result suggest that the slab stagnates near the Hokkaido bend. In this region we detect a weak P660P arrival at a depth of 680–750 km. This is consistent with the Clapeyron slope of the post-spinel phase transition, as lowered temperatures of the slab increase the depth of the 660 km discontinuity, however, the lowered temperatures alone do not explain an increase in the impedance contrast at the phase transition, allowing us to detect P660P. In addition to the thermal heterogeneity introduced by the downgoing slab, there will be a significant change in the composition of materials at 660 km depth, in particular the aluminum-content. Transformations in the non-spinel components of the mantle incorporating Al can contribute significantly to the change in elasticity across the 660 km discontinuity (Akaogi et al., 2002; Wang et al., 2004). The slab consists of a 10–15 km thick layer of mid-ocean ridge basalt (MORB), and 80–90 km of depleted lithosphere and underlying entrained mantle. The MORB veneer has N15 wt.% Al2O3 contained within the garnet phases (Weidner and Wang, 1998). The aluminum is extracted from the melting of the underlying mantle at the mid-ocean ridge and concentrated into the crust, leaving an Aldepleted harzburgitic lithosphere. Changes in densities and seismic velocities for theoretical mineral assemblages present at various mantle temperatures and Al-content were calculated by Weidner and Wang (1998) (hereafter referred to as WW98). We compare our observed discontinuity structure with predicted amplitude ratios calculated using reflectivity synthetic seismograms generated for seismic velocity profiles calculated in WW98 (Fig. 8). In all cases, the amplitudes predicted by WW98 are significantly less than ak135, but similar to our observed amplitude ratios and we are able to reproduce the pattern of amplitudes ratios in Fig. 8c. It is important to stress we are not attempting to replicate the exact amplitudes of WW98, but rather the patterns predicted by the mineral physics calculations for changes in seismic structure. These perturbations in the temperature and chemical content of the mantle produce observable changes in impedance contrast and reflectivity. For ambient, higher aluminum, pyrolitic mantle, the seismic structure calculated by WW98 (Fig. 8a, model for 1900 K and 5% Al), predicts relatively low P660P/PP amplitudes. A low amplitude arrival would correspond to regions away from subduction where we do not detect P660P. In the areas of lowered mantle temperatures (1700 K, and varying aluminum content), the WW98 model predicts two separate discontinuities near 660 km depth, a slightly shallower boundary from the transformation of ilmenite to Mg-perovskite, and the deeper postspinel phase transition. Synthetic waveforms produced for these models indicate that the energy from this second reflector is nearly coincident with the P660P arrival, broadening the precursory pulse but not necessarily detectable as a separate arrival. WW98 predicts a P660P reflector of higher amplitude in colder regions than the reference high aluminum pyrolitic model. For ambient mantle temperatures (1900 K) and Al-depletion (3%), the seismic structure predicted by WW98 produces the highest synthetic P660P/PP amplitude ratios of any of the models. Aluminum depletion reduces the amount of garnet present in the model, increasing the bulk modulus contrast across the phase transition and enhancing the impedance contrast at the 660 km discontinuity (Hirose et al., 1999; Hirose et al., 2005). To the south and west of the weak, deep P660P detections associated with the slab (between 43 and 48 °), the 660 km discontinuity is present as a single discontinuity and not directly associated with the

109

high velocities of the subducting slab, even though it is still relatively deep (660–690 km). It is in this region that we observe consistently high P660P amplitudes (Fig. 7). We propose this localized enhancement arises from the subduction of Al-depleted harzburgitic lithosphere across the 660 km discontinuity. In both regions, the appearance of P660P is consistent with oceanic lithosphere stagnating within the transition zone. WW98 also explored the seismic structure predicted for relatively high mantle temperatures. High mantle temperatures lead to the breakdown of rw to majorite garnet, perturbing the depth and sharpness of the 660 km discontinuity (Deuss, 2009; Deuss et al., 2006; Houser and Williams, 2010). At temperatures N2100 K, the rw phase decomposes into majorite garnet, switching the sign of the Clapeyron slope at the 660 km discontinuity (Hirose, 2002), causing it to deepen in the presence of high temperatures. The seismic structure predicted in WW98 for high temperature (2100 K) produces synthetic P660P amplitudes lower than those for the pyrolitic model, and thus difficult to detect in the Earth. It is interesting to note that tomographic imaging of the subducting Pacific plate suggests a low velocity anomaly lying directly beneath the Kurile slab, leading to the suggestion that high temperature plume material is entrained into the transition zone (Honda et al., 2007). This region corresponds with our observations of a deep 660 beneath the subducting slab at latitudes 43 to 48°. However, as suggested by the synthetic modeling of WW98, this feature is not predicted to sharpen the 660 km discontinuity (e.g., Deuss et al., 2006). Instead, the majorite garnet to Mg-pv phase transition occurs over a much wider depth interval than the post-spinel phase transition, creating a broad gradient and a weak P660P arrival, inconsistent with our observation of strong P660P arrivals within this region. In the vicinity of the high velocity slab in the MTZ, we observe multiple weak P660P reflectors (Figs. 6, 7) near 620–630 km depth and near 680–750 km depth. It is likely that the shallower reflectors are the depth phases of the PP precursors (i.e., pP660P) that stack up as coherent precursory arrivals. In each case we have to calculate the time delay of the depth phase in relation to the main arrival and thus discriminate between cases of two reflectors and a reflector and depth phase (see Fig. 6 for an example). If a strong pP660P is visible, there should also be a strong pP410P phase visible as well. In several bins with multiple P660P observations we also observe an arrival after P410P suggesting the shallower reflector in these events is pP660P (indicated in Figs. 6 and 7). Additional arrivals near 520 km depth are energy from PKiKP mapped onto the migrated profile. Examination of the associated vespagram and synthetic migration allows us to readily identify this energy, and is the motivation for not interpreting any arrivals near 520 km depth. In addition to the vespagrams, the backazimuth slowness analysis further identifies contaminating PKiKP arrivals. With all these complexities taken into consideration, there still appears to be a deep P660P reflector associated with the downgoing Kurile slab. We also consider that folding and buckling of the slab at the increased viscosity of the lower mantle would produce alignment and fabrics within the subducted lithosphere that may locally enhance the impedance contrast at 660 km depth. Lattice preferred orientations of ringwoodite or shape preferred alignments of crustal materials would produce horizontal discontinuities within or in the vicinity of the slab. However, measurements of shear wave splitting beneath the Kuriles finds that the majority of anisotropy observed in the slab is confined to the uppermost mantle (Fischer and Yang, 1994), and it is unclear if the MTZ mineral assemblage is capable of producing the required anisotropy for generating underside reflections (Tommasi et al., 2004). We observe significantly less complex topography on the 410 km discontinuity. Beneath the Sea of Okhotsk, tomography detects lowered mantle velocities, and we find the 410 km discontinuity near 410– 425 km. Near the Kurile slab, the 410 is elevated (395 km) consistent

110

N. Schmerr, C. Thomas / Earth and Planetary Science Letters 311 (2011) 101–111

with the presence of cold subducted lithosphere entering the MTZ. In two events located near the subducting slab we find evidence for a polarity reversal of the 410 km discontinuity (Table 1). A negative polarity P410P would require an impedance drop with increasing depth. One possible explanation is the presence of a metastable wedge of olivine within the subducted lithosphere (e.g.,Castle and Creager, 1998), the higher velocity and density of wadsleyite overlying a wedge of metastable olivine within the slab would produce the conditions necessary for a negative P410P reflection. However, events sampling near these observations do not detect the strong amplitude and depth variations associated with a metastable wedge, and the Fresnel zone of the PP precursors would make the detection of such a feature difficult to resolve (Chaljub and Tarantola, 1997). The P410P/PP amplitude ratios to the northwest of the Kurile slab are lower than the ratios observed under the Pacific (Fig. 7). Mantle transition zone mineralogy is capable of storing several weight percent of H2O (Karato, 2006), and subduction processes are an ideal mechanism for locally enhancing the hydrogen content of the MTZ (Ohtani et al., 2004; Richard et al., 2006). Increased hydrogen content reduces the pressure of the ol-wd phase transition, lowers seismic velocity, and broadens the phase transition stability field (Ohtani and Litasov, 2006). The presence of hydrogen in wadsleyite will reduce the reflection coefficient across the 410 and bring the discontinuity to shallower depths. Alternatively, buoyant hydrated material residing at the top of the transition zone and melting processes would locally enhance the impedance contrast near the slab and produce deeper reflections (Bercovici and Karato, 2003; Schmerr and Garnero, 2007). We do not detect strong variations in 410 topography or amplitude ratios, suggesting that lowered temperatures are capable of explaining our results, though we cannot rule out small-scale compositional heterogeneity at scale-lengths below our resolution (e.g., Zheng et al., 2007). 6. Conclusions The detailed topography on the 410 and 660 km discontinuities beneath the Kuriles subduction zone was imaged using array techniques applied to the PP precursors. The unprecedented dense station spacing provided by the HLP seismic array and USArray Transportable Array enabled this approach. Our detailed array analysis showed P410P reflections from the olivine to wadsleyite phase transition at 410–425 km depth, perturbed to 395 km depth in the vicinity of the subducting Pacific lithosphere. The array analysis also revealed P660P reflections in this area from reflectors at 620–800 km depth, associated with the ilmenite to perovskite and ringwoodite to Mgperovskite + magnesiowüstite phase transitions. This region of robust, deep P660P reflections is coincident with the lowered velocities associated with the subducting Kurile slab, and represents the first detailed analysis of P660P arrivals from beneath the Kuriles. The detection of elevated P660P/PP amplitude ratios indicates the presence of depleted harzburgitic lithosphere at the base of the transition zone; mineral physical experiments and synthetic seismic modeling shows depleted materials enhance the 660 km discontinuity impedance contrast, producing P660P underside reflections. P660P is not detected in regions situated away from subduction, further supporting a chemical enhancement of the discontinuity impedance contrast. The observed thermal and chemical heterogeneity visible in both 410 and 660 reflections is consistent with cold, depleted lithosphere in the MTZ, stagnating at the 660 km discontinuity beneath the Kuriles subduction zone, implying post-spinel boundary can locally impede mantle flow and produce chemical heterogeneity within the transition zone. Acknowledgments We are grateful for the constructive comments and suggestions provided by the editor Peter Shearer and two anonymous reviewers that helped to greatly improve the quality of the manuscript. Data

were collected with the Standing Order for Data (Owens et al., 2004) software. We thank the HLP project for providing access to their data. Data from the TA network were made freely available as part of the EarthScope USArray facility supported by the National Science Foundation. The facilities of the IRIS Data Management System, and specifically the IRIS Data Management Center, were used for access to waveform, metadata or products required in this study. The IRIS-DMS is funded through the National Science Foundation and specifically the GEO Directorate through the Instrumentation and Facilities Program of the National Science Foundation. Data analyses were performed with the TauP (Crotwell et al., 1999), Seismic Handler (Stammler, 1993), and Seismic Analysis Code (Goldstein et al., 2003) toolkits. Figures were generated using GMT (Wessel and Smith, 1998). N.S. was supported by a Department of Terrestrial Magnetism Postdoctoral Fellowship, and C. T. was supported under grant DFG1530/2-1.

References Akaogi, M., Tanaka, A., Ito, E., 2002. Garnet–ilmenite–perovskite transitions in the system Mg4Si4O12–Mg3Al2Si3O12 at high pressures and high temperatures: phase equilibria, calorimetry and implications for mantle structure. Phys. Earth Planet. Inter. 132, 303–324. Bassin, C., Laske, G., Masters, G., 2000. The current limits of resolution for surface wave tomography in North America. Eos Transactions 81, Fall Meeting Suppl., San Francisco, CA, p. F897. Benz, H.M., Vidale, J.E., 1993. Sharpness of upper-mantle discontinuities determined from high-frequency reflections. Nature 365, 147–150. Bercovici, D., Karato, S., 2003. Whole-mantle convection and the transition-zone water filter. Nature 425, 39–44. Bina, C.R., Helffrich, G., 1994. Phase-transition clapeyron slopes and transition zone seismic discontinuity topography. J. Geophys. Res. Solid Earth 99, 15853–15860. Cao, Q., Wang, P., van der Hilst, R.D., de Hoop, M.V., Shim, S.H., 2010. Imaging the upper mantle transition zone with a generalized Radon transform of SS precursors. Phys. Earth Planet. Inter. 180, 80–91. Carlson, R.W., James, D.E., Fouch, M.J., Grove, T.L., Hart, W.K., Grunder, A.L., Duncan, R.A., Keller, G.R., Harder, S.H., Kincaid, C.R., 2005. On the cause of voluminous magmatism in the northwestern United States Abstracts with Programs Geol. Soc. Am. 37, 125. Castle, J.C., Creager, K.C., 1998. NW Pacific slab rheology, the seismicity cutoff, and the olivine to spinel phase change. Earth Planets Space 50, 977–985. Chaljub, E., Tarantola, A., 1997. Sensitivity of SS precursors to topography on the uppermantle 660-km discontinuity. Geophys. Res. Lett. 24, 2613–2616. Chambers, K., Deuss, A., Woodhouse, J.H., 2005. Reflectivity of the 410-km discontinuity from PP and SS precursors. J. Geophys. Res. Solid Earth 110. doi:10.1029/ 2004JB003345. Choy, G.L., Richards, P.G., 1975. Pulse distortion and Hilbert transformation in multiply reflected and refracted body waves. Bull. Seismol. Soc. Am. 65, 55–70. Crotwell, H.P., Owens, T.J., Ritsema, J., 1999. The TauP toolkit: flexible seismic travel-time and ray-path utilities. Seismol. Res. Lett. 70, 154–160. Davies, D., Kelly, E.J., Filson, J.R., 1971. Vespa process for analysis of seismic signals. Nature Phys. Sci. 232, 8–13. DeMets, C., 1992. Oblique convergence and deformation along the Kuril and Japan trenches. J. Geophys. Res. Solid Earth 97, 17615–17625. Deuss, A., 2009. Global observations of mantle discontinuities using SS and PP precursors. Surv. Geophys. 30, 301–326. Deuss, A., Woodhouse, J., 2001. Seismic observations of splitting of the mid-transition zone discontinuity in Earth's mantle. Science 294, 354–357. Deuss, A., Redfern, S.A.T., Chambers, K., Woodhouse, J.H., 2006. The nature of the 660kilometer discontinuity in Earth's mantle from global seismic observations of PP precursors. Science 311, 198–201. Ding, X.Y., Grand, S.P., 1994. Seismic structure of the deep Kurile subduction zone. J. Geophys. Res. Solid Earth 99, 23767–23786. Dziewonski, A.M., Anderson, D.L., 1981. Preliminary reference earth model. Phys. Earth Planet. Inter. 25, 297–356. Efron, B., Tibshirani, R., 1986. Bootstrap methods for standard errors, confidence intervals, and other measures of statistical accuracy. Stat. Sci. 1, 54–75. Estabrook, C.H., Kind, R., 1996. The nature of the 660-kilometer upper-mantle seismic discontinuity from precursors to the PP phase. Science 274, 1179–1182. Fedotov, S.A., 1965. Upper mantle properties of southern part of Kuril island arc according to detailed seismological investigation data. Tectonophysics 2, 219–225. Fischer, K.M., Yang, X.P., 1994. Anisotropy in Kuril-Kamchatka subduction zone structure. Geophys. Res. Lett. 21, 5–8. Fischer, K.M., Jordan, T.H., Creager, K.C., 1988. Seismic constraints on the morphology of deep slabs. J. Geophys. Res. Solid Earth Planets 93, 4773–4783. Flanagan, M.P., Shearer, P.M., 1998. Global mapping of topography on transition zone velocity discontinuities by stacking SS precursors. J. Geophys. Res. Solid Earth 103, 2673–2692. Fuchs, K., Müller, G., 1971. Computation of synthetic seismograms with the reflectivity method and comparison with observations. Geophys. J. R. Astron. Soc. 23, 417–433.

N. Schmerr, C. Thomas / Earth and Planetary Science Letters 311 (2011) 101–111 Fukao, Y., 1977. Upper mantle P-structure on ocean side of Japan-Kurile arc. Geophys. J. R. Astron. Soc. 50, 621–642. Fukao, Y., Obayashi, M., Inoue, H., Nenbai, M., 1992. Subducting slabs stagnant in the mantle transition zone. J. Geophys. Res. Solid Earth 97, 4809–4822. Fukao, Y., Widiyantoro, S., Obayashi, M., 2001. Stagnant slabs in the upper and lower mantle transition region. Rev. Geophys. 39, 291–323. Fukao, Y., Obayashi, M., Nakakuki, T., Deep Slab Project, G, 2009. Stagnant slab: a review. Ann. Rev. Earth Planet. Sci. 37, 19–46. Gaherty, J.B., Lay, T., Vidale, J.E., 1991. Investigation of deep slab structure using longperiod S-waves. J. Geophys. Res. Solid Earth 96, 16349–16367. Goldstein, P., Dodge, D., Firpo, M., Minner, L., 2003. SAC2000: signal processing and analysis tools for seismologists and engineers. The IASPEI International Handbook of Earthquake and Engineering Seismology. Academic Press, London. Gorbatov, A., Kennett, B.L.N., 2003. Joint bulk-sound and shear tomography for Western Pacific subduction zones. Earth Planet. Sci. Lett. 210, 527–543. Grand, S.P., 1994. Mantle shear structure beneath the America and surrounding oceans. J. Geophys. Res. Solid Earth 99, 11591–11621. Gu, Y., Dziewonski, A.M., Agee, C.B., 1998. Global de-correlation of the topography of transition zone discontinuities. Earth Planet. Sci. Lett. 157, 57–67. Gudmundsson, O., Sambridge, M., 1998. A regionalized upper mantle (RUM) seismic model. J. Geophys. Res. Solid Earth 103, 7121–7136. Hirose, K., 2002. Phase transitions in pyrolitic mantle around 670-km depth: Implications for upwelling of plumes from the lower mantle. J. Geophys. Res. Solid Earth 107, 2078 doi:2010.1029/2001JB000597. Hirose, K., Fei, Y.W., Ma, Y.Z., Mao, H.K., 1999. The fate of subducted basaltic crust in the Earth's lower mantle. Nature 397, 53–56. Hirose, K., Takafuji, N., Sata, N., Ohishi, Y., 2005. Phase transition and density of subducted MORB crust in the lower mantle. Earth Planet. Sci. Lett. 237, 239–251. Honda, S., Morishige, M., Orihashi, Y., 2007. Sinking hot anomaly trapped at the 410 km discontinuity near the Honshu subduction zone. Jpn. Earth Planet. Sci. Lett. 261, 565–577. Houser, C., Williams, Q., 2010. Reconciling Pacific 410 and 660 km discontinuity topography, transition zone shear velocity patterns, and mantle phase transitions. Earth Planet. Sci. Lett. 296, 255–266. Ito, E., Takahashi, E., 1989. Postspinel transformations in the system Mg2SiO4–Fe2SiO4 and some geophysical implications. J. Geophys. Res. Solid Earth Planets 94, 10637–10646. Jordan, T.H., 1976. Lithospheric slab penetration into lower mantle. TransactionsAmerican Geophysical Union, 57. 962-962. Jordan, T.H., 1977. Lithospheric slab penetration into lower mantle beneath Sea of Okhotsk. J. Geophys. Z. Geophys. 43, 473–496. Karato, S., 2006. Remote sensing of hydrogen in Earth's mantle. In: Keppler, H., Smyth, J.R. (Eds.), Water in Nominally Anhydrous Minerals. : Reviews in Mineralogy & Geochemistry, 62. Geochemical Society, St. Louis, MO, pp. 343–375. Katsura, T., Ito, E., 1989. The system Mg2SiO4–Fe2SiO4 at high-pressures and temperatures —precise determination of stabilities of olivine, modified spinel, and spinel. J. Geophys. Res. Solid Earth Planets 94, 15663–15670. Kennett, B.L.N., Engdahl, E.R., Buland, R., 1995. Constraints on seismic velocities in the earth from travel-times. Geophys. J. Int. 122, 108–124. King, D.W., Haddon, R.A.W., Husebye, E.S., 1975. Precursors to PP. Phys. Earth Planet. Inter. 10, 103–127. Lebedev, S., Chevrot, S., van der Hilst, R.D., 2003. Correlation between the shear-speed structure and thickness of the mantle transition zone. Phys. Earth Planet. Inter. 136, 25–40. Li, C., van der Hilst, R.D., Engdahl, E.R., Burdick, S., 2008. A new global model for P wave speed variations in Earth's mantle. Geochem. Geophys. Geosyst. 9, Q05018 doi:05010.01029/02007GC001806. Megnin, C., Romanowicz, B., 2000. The three-dimensional shear velocity structure of the mantle from the inversion of body, surface and higher-mode waveforms. Geophys. J. Int. 143, 709–728. Miller, M.S., Kennett, B.L.N., 2006. Evolution of mantle structure beneath the northwest Pacific: evidence from seismic tomography and paleogeographic reconstructions. Tectonics 25, TC4002 doi:4010.1029/2005TC001909. Miller, M.S., Kennett, B.L.N., Gorbatov, A., 2006. Morphology of the distorted subducted Pacific slab beneath the Hokkaido corner. Jpn. Phys. Earth Planet. Inter. 156, 1–11. Muirhead, K.J., Datt, R., 1976. N-Th root process applied to seismic array data. Geophys. J. R. Astron. Soc. 47, 197–210. Müller, R.D., Roest, W.R., Royer, J.Y., Gahagan, L.M., Sclater, J.G., 1997. Digital isochrons of the world's ocean floor. J. Geophys. Res. Solid Earth 102, 3211–3214. Neele, F., Snieder, R., 1992. Topography of the 400-km discontinuity from observations of long-period P400p phases. Geophys. J. Int. 109, 670–682. Ohtani, E., Litasov, K.D., 2006. The effect of water on mantle phase transitions, water in nominally anhydrous minerals. Rev. Mineral. Geochem. 62, 397–419.

111

Ohtani, E., Litasov, K., Hosoya, T., Kubo, T., Kondo, T., 2004. Water transport into the deep mantle and formation of a hydrous transition zone. Phys. Earth Planet. Inter. 143–44, 255–269. Owens, T.J., Crotwell, H.P., Groves, C., Oliver-Paul, P., 2004. SOD: standing order for data. Seismol. Res. Lett. 75, 515–520. Richard, G., Bercovici, D., Karato, S.I., 2006. Slab dehydration in the Earth's mantle transition zone. Earth Planet. Sci. Lett. 251, 156–167. Ritsema, J., Deuss, A., van Heijst, H.J., Woodhouse, J.H., 2011. S40RTS: a degree-40 shear-velocity model for the mantle from new Rayleigh wave dispersion, teleseismic traveltime and normal-mode splitting function measurements. Geophys. J. Int. 184, 1223–1236. Rost, S., Thomas, C., 2002. Array seismology: methods and applications. Rev. Geophys. 40. doi:10.1007/S10712-10009-19070-10716. Rost, S., Thomas, C., 2009. Improving seismic resolution through array processing techniques. Surv. Geophys. 30, 271–299. Rost, S., Weber, M., 2002. The upper mantle transition zone discontinuities in the Pacific as determined by short-period array data. Earth Planet. Sci. Lett. 204, 347–361. Schmerr, N., Garnero, E., 2007. Upper mantle discontinuity topography from thermal and chemical heterogeneity. Science 318, 623–626. Shearer, P.M., 1990. Seismic imaging of upper-mantle structure with new evidence for a 520-km discontinuity. Nature 344, 121–126. Shearer, P.M., 1991. Constraints on upper mantle discontinuities from observations of long-period reflected and converted phases. J. Geophys. Res. Solid Earth 96, 18147–18182. Shearer, P.M., 2000. Upper Mantle Discontinuities. In: Karato, S., Forte, A.M., Liebermann, R.C., Masters, G., Stixrude, L. (Eds.), Earth's Deep Interior: Mineral Physics and Tomography From the Atomic to the Global Scale, Geophysical Monograph 117. AGU, Washington D.C., pp. 115–131. Shearer, P.M., Flanagan, M.P., 1999. Seismic velocity and density jumps across the 410and 660-kilometer discontinuities. Science 285, 1545–1548. Shearer, P.M., Masters, T.G., 1992. Global mapping of topography on the 660-km discontinuity. Nature 355, 791–796. Song, T.R.A., Helmberger, D.V., Grand, S.P., 2004. Low-velocity zone atop the 410-km seismic discontinuity in the northwestern United States. Nature 427, 530–533. Stammler, K., 1993. Seismichandler programmable multichannel data handler for interactive and automatic processing of seismological analyses. Comput. Geosci. 19, 135–140. Stixrude, L., 1997. Structure and sharpness of phase transitions and mantle discontinuities. J. Geophys. Res. Solid Earth 102, 14835–14852. Suetsugu, D., 1989. Lower mantle high-velocity zone beneath the Kurils as inferred from P-wave travel time and amplitude data. J. Phys. Earth 37, 265–295. Tajima, F., Grand, S.P., 1995. Evidence of high-velocity anomalies in the transition zone associated with southern Kurile subduction zone. Geophys. Res. Lett. 22, 3139–3142. Thomas, C., Billen, M.I., 2009. Mantle transition zone structure along a profile in the SW Pacific: thermal and compositional variations. Geophys. J. Int. 176, 113–125. Tommasi, A., Mainprice, D., Cordier, P., Thoraval, C., Couvy, H., 2004. Strain-induced seismic anisotropy of wadsleyite polycrystals and flow patterns in the mantle transition zone. J. Geophys. Res. Solid Earth 109, B12405 dot:12410.11029/ 12004jb003158. van der Hilst, R., 1995. Complex morphology of subducted lithosphere in the mantle beneath the Tonga Trench. Nature 374, 154–157. van der Hilst, R.D., Engdahl, E.R., Spakman, W., 1993. Tomographic inversion of P-data and PP-data for aspherical mantle structure below the northwest Pacific region. Geophys. J. Int. 115, 264–302. Wang, Y.B., Uchida, T., Zhang, J.Z., Rivers, M.L., Sutton, S.R., 2004. Thermal equation of state of akimotoite MgSiO3 and effects of the akimotoite–garnet transformation on seismic structure near the 660 km discontinuity. Phys. Earth Planet. Inter. 143, 57–80. Weidner, D.J., Wang, Y.B., 1998. Chemical- and clapeyron-induced buoyancy at the 660 km discontinuity. J. Geophys. Res. Solid Earth 103, 7431–7441. Wessel, P., Smith, W.H.F., 1998. New, improved version of Generic Mapping Tools released. EOS Trans. AGU 79, 579. Widiyantoro, S., van der Hilst, R., 1997. Mantle structure beneath Indonesia inferred from high-resolution tomographic imaging. Geophys. J. Int. 130, 167–182. Xu, F., Vidale, J.E., Earle, P.S., 2003. Survey of precursors to P ' P ': fine structure of mantle discontinuities. J. Geophys. Res. Solid Earth 108. Zheng, Y.C., Lay, T., Flanagan, M.P., Williams, Q., 2007. Pervasive seismic wave reflectivity and metasomatism of the Tonga mantle wedge. Science 316, 855–859.