THE CRUST AND UPPERMOST MANTLE S-WAVE VELOCITY STRUCTURE BENEATH

We have studied the crust and uppermost mantle S-wave velocity structure beneath Japan Islands by using the teleseismic waveform data above Mw 6.0 recorded from January 2006 to February 2017 by 15 Hi-Net stations and the joint inversion technique of Pand S-wave receiver functions based on the Bayes theory. The results show that, beneath Japan Islands, the crust appear the characteristic of thinner in the south and east, and thicker in the north and west. The thinnest and thickest crust in the study region locate at station JSD (26km) and JGF(44km), respectively. The horizontal distribution of S-wave velocity in the study region are relatively complicated within upper and middle crust depth, while from the lower crust to the uppermost mantle depth, the velocity distribution is relatively uniform. The large earthquakes (above Mw 6.0) mainly took place at the edge of the high and low velocity zones.


INTRODUCTION
The Japan Islands locate among four active lithospheric plates (shown in Figure1). The subduction and collision of these plates has caused many destructive earthquakes in the islands or the surrounding region [1,2]. Many seismic researches have been conducted in order to detect the crust and upper mantle structure beneath Japan Islands. Using seismic body wave travel time tomography method and the seismic data recorded by more and more densely distribution of seismic stations in Japan, some scientists have obtained detailed P wave and S wave velocity perturbation images within upper mantle depth which has greatly improved the understanding about the deep structure and the tectonics beneath Japan and the relationship between large earthquakes and the deep velocity structure [3][4][5][6][7][8]. To the purpose of exploring the deep tectonic causes of 2011 Tohoku earthquake, some scientists studied the seismic attenuation of Northeast of Japan arc [9]. By the combination of seismic waveform and strong motion data and GPS observations, scholars also studied the rupture process of 2011 Sariku-Oki earthquake [7]. Some of researchers analyzed the generating mechanism of 2011 Tohoku earthquake and the induced tsunami by high resolution seismic images of the northeast Japan forearc [10].
Comparing with seismic travel time tomography technique, teleseismic receiver function method takes advantage of velocity discontinuities identification and obtaining absolute velocity information within crust and upper mantle depth (about 100km depth). The 410km and 660km discontinuities depth beneath Japan subduction zone have been studied by using P wave receiver function migration and summation technique [11,12]. They proposed a new receiver function migration method based on waveequation post stack, and applied it to study the subsurface structures of the Japan subduction zone, and obtained the slab image and the dipping angle of Pacific plate subduction down to Japan [13]. However, by using the receiver function migration method, only the interfaces lateral variation within crust and upper mantle depth could be obtained. P-and S-wave receiver function joint inversion method has only been proposed and developed in recent years, by which the layered P-wave and S-wave velocity structure models within crust and uppermost mantle depth (about 300km) beneath some study region could be obtained [14][15][16]. In this paper, we would like to use P-and S-wave receiver function joint inversion technique based on Bayes theory to analyze the crust and uppermost mantle S-wave velocity structure, which is much different in many aspects comparing with global inversion method [14][15][16][17].
Furthermore, only S-wave velocity structure within 100km depth is determined by the observed P-wave and S-wave receiver functions, which implies fewer target variables, less inversion uncertainties.

DATA AND METHOD
The seismic waveforms used in this study are recorded by 15 stations of Hi-net (Figure 1), which are operated by the Japanese national universities, Japan Meteorological Agency (JMA) and National Research Institute for Earth Science and Disaster Prevention (NIED) for the High -sensitivity seismic network. Totally 292 records of teleseismic events with Mw≥6.5 from January 2006 to February 2017 with epicentral distances in the range of 30-90 degrees for P-wave receiver functions (PRF) and 60-80 degrees for S-wave receiver functions (SRF), are collected. We selected records with high signal-to-noise ratio and clear onset of P-waves and S-waves. The waveforms were rotated from the north-east-vertical (N-E-Z) to the radial-transverse-vertical (R-T-Z) coordinate using the back-azimuth. For the PRFs, the three-component records were then cut in the time window of 5s prior to and 20s after the P-arrival. We then constructed the receiver functions by deconvolving the vertical component from radial component using an iterative approach [18]. For the SRFs, the time window is selected from -25s to 25s and also constructed by the iterative approach. A low-pass Gaussian filter was applied to smooth the PRFs to 1.5Hz and the SRFs to 0.5Hz. We then did the summation of PRFs and SRFs from different azimuth for each station, which could suppress the laterally heterogeneity of the media beneath each station and could also enhance the amplitudes of conversion waves. The PRF and SRF waveform data after summation of each station then were transformed into spectral domain by FFT algorithm. So far, the data for joint inversion were ready as our inversion would be processed in spectral domain, which would speed up the calculation and reduce the computation It is necessary to point out that we only invert the S-wave velocity instead of both P and S, which means the more reliable results could be expected for the numbers of target variables were reduced. The detailed processes of joint inversion would not be introduced here (For the details, the reader is referred to the other paper by Wang and Liu [17]) Some examples for the joint inversion in Japan Islands have been given ( Figure 2). Three different and possible initial models (PREM, IASP91, and CRUST 1.0) (Left panel of each subgraph in Figure 2) were attempted for the inversion of each station to test the inversion reliability and stability. From the final models (Middle panels), much similar results could be observed from different initial models (Figure 2 (a) and (b)), but the final model's consistency for the third station JSD is not so good as the upper two stations (ASAJ and JMM). We infer to two reasons, the first of which is that the station JSD located on the sea (Figure 1), the crust structure of which is much different with the stations located on the continent. The second reason might be the forward calculation in the inversion process was based on the global continental model.  Figure 1). In each subgraph, the left panel represents Initial models (PREM (red), IASP91 (black) and CRUST 1.0 (green)) for the inversion. The middle panel shows final models after inversion. The line color is the same as shown in the left panel. Moho depth is denoted by black arrow. The waveform fitting of P-and S-wave receiver functions of each station for the joint inversion are shown in the right panel. Black lines represent Observed receiver functions while color lines represent predicted theoretical receiver functions after joint inversion from different initial models. P-to-s and S-to-p conversion waves from Moho are marked out by black arrows corresponding to the Moho depth in the middle panel.

RESULTS
By using the PRF and SRF joint inversion method mentioned above, we obtained 15 stations results, including PRFs and SRFs profile, and S-wave velocity model of each station. Three profiles (Figure 1) were designed in order to analyze the horizontal variation characteristic of S-wave velocity beneath Japan Islands. Obvious lateral variations of Moho depth could be observed from profile A-A' (Figure 3), which is along near south-to-north direction from Philippine Sea plate pass through east of Japan Island to the northeast of Japan (Figure 1). The Moho depth increase from South to North along profile A-A', which could be supported by both PRF profile (Figure 3 (a)) and S-wave velocity profile (Figure 3(c)). The shallow Moho depth in this profile is about 30km, located at Philippine Sea plate, while the station JEW located at northeast Japan has the deepest Moho, which is about 40km. We couldn't catch much deep structure information from SRF profile (Figure 3 (b)), but it would provide effective constrains for the joint inversion. The sedimentary layers beneath two stations (JHJ2 and JEW) could be observed in the PRF profile (special waveforms nearby 0 sec) and S-wave velocity model Obvious low velocity values at 0km depth could be observed in Figure  3 (a) and (c), which are marked out by red arrows.
Another profile B-B' along near south-to-north direction passes through the west part of Japan Islands. The remarkable lateral fluctuation of Moho along this profile could be easily observed from both PRFs (Figure 4 (a)) and S velocity models (Figure 4 (c)). Taking the station JSD as the center, to the south, The Moho is about 42-44km, and 36-42km to the north. The Moho depth is only 26km beneath the station JSD located in the sea.  Figure 1).
The profile C-C' (Figure 5) passes through the center of Japan Islands along southwest-to-northeast direction (Figure 1). The Moho appears much deeper and with good con- Figure 5: The same as Figure 3, but for profile C-C' (location shown in Figure 1). sistence (42-44km) to the west of station MAJ, and shallower (30-34km) to the east of MAJ. An obvious positive phase (red color) at about 3-4 sec could be identified between station JWT to station MAJ from the PRF profile ( Figure 5 (a)), which might be the conversion wave from Conrad interface. This phase could also be observed in west of Japan Islands (Figure 4 (a)), but couldn't be found in the east of the Islands. (Figure 3 (a)).
To further analyze the results, the S-wave velocity images at different depths from 0km to 100km with 5km or 10km interval were drawn by the two-dimensional in Figure 6: The obtained S-wave velocity distribution at different depth (0 to 25km with 5km interval).
Red and blue colors represent low and high velocities, respectively. The absolute velocity value is shown at right. The purple triangles denote stations. The yellow circles, blue circles and the red star represent epicenters of earthquakes above Mw6.0, Mw6.5 and Mw7.0, respectively. trepidation technique (Figure 6-8). The S-wave velocity changes gradually with depth increasing, which could be observed from 0-25km images ( Figure 6). The relatively uniform image at 0km (Figure 6 (a)) becomes much more complicated to the depth of 5-25km ( Figure 6 (bf)). We could also find that the epicenter of large earthquakes (up to Mw 6.0) mainly locate at the edge of high S-wave velocity bodies and the low ones. From the image of 25km (Figure 6 (f)), at which depth the Tohoku Mw 9.0 earth Figure 7: The same as Figure 6, but for the depth of 30 to 55km with 5km interval.  Figure 6, but for the depth of 60 to 100km with 10km interval. quake, 2011 took place, we couldn't get the exact velocity structure beneath the epicenter, but we could observe that to the southwest of the epicenter the S-wave velocity shows special characteristic of lowhigh-low.
This feature extends to 35km depth (Figure 7 (a), (b)). From 40km to 100km, the other simple characteristic could be concluded that the high velocity shown in the large part of the northeast of Japan Islands, and the low velocity appear in a relatively small region of the southwest of the Islands. The exceptions appear in the images of 50km (Figure 7(e)) and 70km (Figure 8 (b)), which show more complicated lateral changes. The upper boundary of the subducting Pacific slab appear no close relationship with the velocity distribution (Figure 8 (a), (c), (e)).

DISCUSSION
In this paper, we studied the crust and uppermost mantle S-wave velocity structure beneath Japan Islands by using the teleseismic waveform data above Mw 6.0 recorded from January 2006 to February 2017 by 15 Hi-Net stations and the joint inversion technique of P-and Swave receiver functions based on the Bayes theory. The results show that, beneath Japan Islands, the crust depth appear the characteristic of thinner in the south and east, and thicker in the north and west. The thinnest and thickest crust in this region locate at station JSD (26km) and JGF(44km), respectively. The horizontal distribution of S-wave velocity in the study region are relatively complicated within upper and middle crust depth, while from the lower crust to the uppermost mantle depth, the velocity distribution is relatively uniform. The large earthquakes (above Mw 6.0) mainly took place at the edge of the high and low velocity zones.
This paper aims at testing the new joint inversion technique in Japan Islands, so only the data recorded by 15 stations of Hi-Net have been used. More detailed crust and upper mantle structure information could be expected, if more data recorded by dense seismic array were processed by using this method. For the complication of inversion problems, all the inversion techniques must face the problem of results non-uniqueness. The joint inversion technique used in this paper could be developed by adding surface wave or ambient noise dada as the inversion constrains.