Scattering strength at active volcanoes in Japan as inferred from the peak ratio analysis of teleseismic P waves

Small-scale seismic velocity heterogeneity has been studied through the calculation of peak amplitude ratio as a means to quantify the strength of seismic wave scattering at volcanoes in Japan. This ratio is defined as the ratio of the maximum (peak) P wave energy in the transverse component seismogram envelope over that of the three-component sum seismogram envelope (transverse + radial + vertical). According to the previous study using Japan’s Hi-net seismometer network, the peak ratio is observed to be larger near the (active) quaternary volcanoes. However, these Hi-net stations are not positioned on the volcanoes themselves. This study systematically examines the peak ratios at 47 active volcanoes across Japan, using seismometers operated by the Japan Meteorological Agency (JMA). Analyses were performed at four frequency bands: 0.5–1, 1–2, 2–4, and 4–8 Hz. We found that the JMA stations yield higher peak ratios than the Hi-net stations. Their differences are statistically significant at the 99.9% confidence level in all frequency bands. We also examined the differences between the ground surface and borehole stations of the JMA network. The former shows larger peak ratios, and for most frequency bands, the differences are also statistically significant at the 99.9% confidence level. This suggests an intensification of small-scale medium heterogeneities especially at shallow depths at active volcanoes, and that scattering might have been enhanced at the very shallow parts.


Introduction
Small-scale heterogeneities in seismic velocity exist in the crust and lithosphere, and they can invoke the scattering of seismic waves (Sato et al. 2012). In Japan, the far-reaching and densely spaced Hi-net seismometer network allows for the broad applications of seismic scattering studies. For example, Carcolé and Sato (2010) estimated the spatial distribution of the intrinsic and scattering Q factors by applying the multiple lapse time window analysis (MLTWA) to the S wave-coda energies of local earthquakes. In addition, Takahashi et al. (2009) estimated the 3D distribution of the small-scale heterogeneity's power spectral density function (PSDF) in the Tohoku and Hokkaido regions. Based on the Markov approximation of the parabolic wave equation, they analysed the peak delay times of high frequency S wave from local microseismic events at multiple frequency bands. The estimated heterogeneity is spatially well correlated with the distribution of quaternary volcanoes.
One alternative to quantify the strength of seismic scattering and thus heterogeneity is the computation of peak amplitude ratio (hereafter called peak ratio). This approach is based on the partition of P wave energies into the transverse component due to scattering. In theory, obliquely incoming direct P waves could only generate shaking in the vertical and radial directions due to the longitudinal nature of the particle movement in horizontally layered media. On the other hand, scattered P wave arrivals off three-dimensional medium heterogeneities can be detected in all three components. Therefore, transverse component seismogram contains information on the scattered waves, and can be used to measure the scattering strength in the medium underneath the seismic station.
The peak ratio is defined as the ratio of the maximum P wave energy between the transverse component envelope and the 3-component sum envelope. According to Sato (2006), assuming the Gaussian type random medium, peak ratio can be approximated as 1.81ε 2 h/a where ε describes the root-mean-squared fractional fluctuation in seismic velocity which controls the amplitude of velocity heterogeneity within the medium, a is the correlation distance which controls the size and shape of patches of velocity heterogeneity, and h is the thickness of the heterogeneous medium. In addition, theoretical models (e.g., Sato 2006) and numerical simulations (e.g., Takemura et al. 2015) predict a positive correlation between peak ratio and the peak delay time, representing the lag time from the onset to the maximum amplitude. This follows from the notion of envelope broadening, from which a longer delay suggests that more energy is being partitioned into the transverse component, hence an increase in the peak ratio.
The peak ratio approach has been applied to evaluate the scattering strength. On the global scale, Kubanza et al. (2007) observed that the peak ratios on the stable continents are small in comparison to those at the active tectonic regions. On the regional scale, Nishimura (2012) analysed the envelopes of teleseismic P waves recorded by the Hi-net borehole stations in Japan, and found that tectonically active regions, for example, those in the proximity of quaternary volcanoes or active faults, showed higher peak ratios. That is presumably due to an increase in heterogeneous elements following the increase in tectonic activity.
The study of Nishimura (2012) touched on the premise that high peak ratio can be expected for stations close to quaternary volcanoes. However, Hi-net seismometers are not installed on the actual volcanic edifice. Therefore, to specifically examine the peak ratio estimation in volcanic settings, we shall make use of seismometers that the Japan Meteorological Agency (JMA) operates to monitor 47 active volcanoes in Japan. These stations are mostly deployed within 10 km of the active craters. Thus, by utilising them, we can categorically quantify the strength of seismic wave scattering at the quaternary volcanoes. One advantage of the peak ratio approach is that it does not require the observation of local earthquakes in the study area (Nishimura 2012). In other words, it enables us to perform a systematic survey of the peak ratio at volcanoes across Japan, using the same set of teleseismic events.

Data and analyses
We use 99 teleseismic events recorded by the JMA network between 2010 and 2019. In accordance to the previous study, we choose only earthquakes deeper than 300 km so that contamination by the depth phases could be avoided (Nishimura 2012). The choosing of only deep events ensures that hypocentral distances are suitably far enough for the analysis window to exclude direct S wave arrivals. The magnitude of the events ranges between 5 and 6.5, with epicentral distance ranging between 0° and 60°. The back-azimuths are mostly distributed in the NE-NW and the SE-SW directions. See Additional file 1: Fig.  S1 for the distribution of these earthquakes.
We follow the procedures developed in Kubanza et al. (2007) to compute the mean-squared (MS) envelopes at individual stations and calculate their peak ratios. After correcting for the instrumental response, we bandpass filter the data into four frequency bands: 0.5-1, 1-2, 2-4, and 4-8 Hz. We rotate the data into the vertical-radialtransverse (ZRT) components, according to the backazimuths for the stations and earthquakes coordinates. Source-receiver pairs with signal-to-noise ratio (S/N) of less than 10 are discarded. We compute only the envelopes of the remaining pairs by squaring the seismogram traces. Since stations record events of various distances and magnitudes, the squared seismograms are normalised in accordance to Eq. 3 of Kubanza et al. (2007). For stations that record more than 10 earthquakes, we stack their seismogram envelopes to produce the MS envelopes. We calculate the peak ratio by dividing the peak of the transverse component MS envelope over that of the three-component sum MS envelope.
Two examples of MS envelopes at two stations are shown in Fig. 1. These two stations are shown because and I peak sum , respectively, from which the peak amplitude ratio is calculated. Time is given in relation to the theoretical P wave arrival according to the ak135 model. Theoretical S wave arrival is outside the given time window both are borehole stations located close to the Sakurajima volcano in the region of Kyushu, but they can be differentiated by the network that they belong to. Station SFT2 is one of the JMA stations that monitor the Sakurajima volcano, whereas station AIRH is the closest Hi-net station to Sakurajima. For stations AIRH and SFT2, the peak ratios are estimated to be 0.217 and 0.294, respectively, for the 2-4 Hz frequency band in the examples.
We also calculate the transverse peak delay time which we define as the time between the direct P arrival time, and the transverse MS envelope peak time. For simplicity, the reference P wave onset is determined based on the ak135 earth-model which gives "a significantly better fit to a broad range of phases than is provided by the iasp91 and sp6 models" (Kennett et al. 1995), and is widely used for example, as the global reference model for 3D teleseismic tomography studies of the Murray Basin in southeast Australia (Rawlinson et al. 2006). In this study, data were normalised and stacked for a 45 s time-window commencing from 10 s prior to the theoretical arrival time. Thus, the time axis in Fig. 1 is aligned and shown in relation to this theoretical P arrival time.

Results
From 204 JMA seismic stations deployed around active volcanoes in Japan, peak ratios were computed at 134, 140, 136, and 109 stations under the frequency bands of 0.5-1, 1-2, 2-4, and 4-8 Hz, respectively. These stations are positioned on active volcanoes that are mainly distributed along the volcanic front. Similarly, the lack of stations in the western part of Japan (e.g., Kinki, Chugoku, and Shikoku regions) reflects the scarcity of active volcanoes in these areas. Figure 2 shows the distribution of the peak ratio at different regions for the frequency band of 2-4 Hz averaged for each volcano. The colour scheme used in Fig. 2 mimics the scale and range used in Nishimura (2012) to allow for a more direct comparison with the Hi-net observation. In this mid frequency band, volcanoes in central Japan and Kyushu generally show medium peak ratios roughly around 0.34. In Hokkaido, we see stretches of high peak ratio estimation around 0.40 and on the contrary, stretches of small peak ratio estimation in the Tohoku region around 0.24. Overall though, peak ratios at volcanoes are consistently of high values-especially considering that peak ratios at the Hinet stations include ratios < 0.12 quite considerably (Fig. 2 of Nishimura 2012 along with Additional file 1: Fig. S2)which suggests strong seismic scattering at the active volcanoes.
The uncertainties in our peak ratio estimations at the JMA stations were evaluated using the delete-1 Jack-knife approach. With this approach, event recordings at individual stations are disregarded one by one, to obtain n estimates of peak ratio where n is the number of events. The standard deviation of this set of n estimates can then be calculated and presented as the coefficient of variation in percentage to express the uncertainty. The maximum errors for the frequency bands of 0.5-1, 1-2, 2-4, and 4-8 Hz are 11.7%, 14.2%, 15.5%, and 17.2%, respectively. However, the averages of these errors are much smaller at 3.9%, 3.3%, 3.5%, and 4.0%, respectively. In other words, the typical errors found in our results are considerably less than 5 percent. Figure 3 shows the histograms of peak ratio distributions at the JMA stations across the four frequency bands. The distribution is normalised as the probability density function. For each bin, the probability density function is calculated by dividing its bin count over the multiplication of the total bin count and the discrete bin interval (= 0.05), making the area under the whole histogram equates to one, i.e., p i w = 100% where p i is the probability density at the i th bin and w is the discrete bin width. Therefore, the probability at a given bin is obtained by multiplying the probability density with the discrete bin width (= 0.05).
The peak ratio averages plus-and-minus one standard deviation for the JMA stations are 0.32 ± 0.12 at the frequency band of 0.5-1 Hz, 0.35 ± 0.13 at 1-2 Hz, 0.36 ± 0.13 at 2-4 Hz, and 0.39 ± 0.12 at the 4-8 Hz frequency band (Table 1; yellow histograms in Fig. 3). On the other hand, the peak ratio averages plus-and-minus one standard deviation for all Hi-net stations at the same frequency bands are 0.14 ± 0.11, 0.14 ± 0.10, 0.16 ± 0.11, and 0.20 ± 0.14, respectively (Nishimura 2012). The averages for the JMA stations are higher than those at Hi-net by 129%, 150%, 125%, and 95% for the 0.5-1, 1-2, 2-4, and 4-8 Hz, respectively. That is, the peak ratios at the JMA stations are about twice as large as the Hi-net stations. On the other hand, their standard deviations are relatively stable, around 0.1, signifying a clear distinction in peak ratio distribution between the data sets. To verify this contrast, we perform the t test which is a statistical tool to compare whether two sample means are statistically different from one another (e.g., Brase and Brase 2013). We confirm that for all frequency bands, the differences in peak ratio between the JMA and Hi-net stations are statistically significant with the confidence level of 99.9%. Furthermore, we compare the peak ratio of the JMA stations only with Hi-net stations categorised in Nishimura (2012) as belonging to the quaternary volcano group, as well as with the quaternary volcano plus active fault group. The quaternary volcano group consists of stations in the proximity (within 20 km) of at least one volcano. The quaternary volcano plus active fault group is exclusive to the previous group, containing stations in the proximity of both a volcano and an active fault. Against the former group, the average differences become 100%, 150%, 112%, and 63% at the four frequency bands, respectively, with the t test results affirming that these differences are statistically significant at the 99.9% confidence level. Against the latter group, the average differences become 78%, 119%, 89%, and 70%, respectively, with the t test results signifying differences statistically significant at the 99.9% confidence level too. Therefore, we still see considerable differences in peak ratio between the JMA stations and the Hi-net subgroups with more similar tectonic settings.
The JMA stations can be further differentiated into the ground surface and borehole stations. We define seismometers to be borehole stations if they are installed at a depth of greater than 45 m from the ground surface. Most borehole stations are installed at depths of roughly 100 m. In Fig. 3, we also include the histograms of the peak ratios at the four frequency bands for the ground surface (red histograms) and borehole (blue histograms) stations. The normalisation allows us to directly compare At all of the frequency bands, the ground surface stations show larger peak ratios than the borehole stations with average differences at the 0.5-1, 1-2, 2-4, and 4-8 Hz frequency bands of 0.044, 0.136, 0.169, and 0.085, respectively (Table 1; Fig. 3). With reference to the borehole stations, this is equivalent to the increase in peak ratio of about 15%, 52%, 66%, and 25%, respectively. We conducted similar t tests to the ground surface and borehole station groups across all frequency bands. The t test results indicate that apart from the 0.5-1 Hz frequency band which returns a confidence level of 95%, the ground surface and borehole stations peak ratios are significantly different from one another at the 99.9% confidence level.
Finally, we examined the transverse peak delay time at the ground surface and borehole stations across the different frequencies (Table 1; Additional file 1: Fig. S8). The averages in the transverse peak delay time for the surface stations are smaller (came earlier) than those for the borehole stations at the lowest frequency bands by 1.583 s or 15%. This is in contrast to the average increase in peak ratio from borehole to surface stations at the same band (Table 1; Fig. 3). On the other hand, at higher frequency bands, the averages in the transverse peak delay time increase by 0.340 s or 4% (1-2 Hz), 1.021 s or Fig. 3 Histograms showing the probability density functions of peak amplitude ratios for the four frequency bands: 0.5-1, 1-2, 2-4, and 4-8 Hz. Yellow histograms show the distributions of peak ratios when all JMA stations are considered. Red and blue histograms separate the peak ratio distributions into JMA ground surface and borehole stations, respectively. The former is showing mostly larger values, and both are greater than the reference numbers for the Hi-net stations. The grey bands indicate the peak ratio mean plus-and-minus one standard deviation for all Hi-net data according to Table 1 of Nishimura (2012) 15% (2-4 Hz), and 0.653 s or 10% (4-8 Hz) from borehole to ground surface stations. More notably, within each frequency band, we do not see a clear disparity differentiating the Hi-net transverse peak delay times and those for the ground surface and borehole stations (Table 1; Additional file 1: Fig. S8), unlike in the case of peak ratio (Table 1; Fig. 3) where we can more easily distinguish the three in a rather orderly fashion (that is, JMA surface > JMA borehole > Hi-net peak ratio).

Conclusions and discussion
Assuming the Gaussian-type infinite random media following Sato (2006), peak ratio can be approximated as 1.81ε 2 h/a . Kubanza et al. (2007) coined the term ε 2 h/a as the quantity of randomness to evaluate the strength of heterogeneity. They estimated this quantity on the global scale across different frequency bands between 0.5 and 4 Hz, and showed that values in Japan are found to be typically between 0.024 and 0.15. Across the same frequency bands, the peak ratios in this study are typically around 0.3 and 0.4, which in terms of ε 2 h/a are estimated to be around 0.17 to 0.22 for volcanoes, which is larger than the Kubanza et al. (2007) estimation.
If we take ε 2 h/a as 0.2 with thickness h of 20-100 km, ε 2 /a becomes in the order of 10 -2 to 10 -3 km −1 at volcanoes in this study. Referring to Table 2 of Sato (2019) which compiles ε and a values from study areas featuring strong heterogeneity, this estimation of ε 2 /a is comparable to the Nikko region in Japan ( ε 2 /a = 7.94 × 10 −3 km −1 ) which is located on the volcanic front, and the Pyrenees mountain range in Europe ( ε 2 /a = 6.73 × 10 −3 km −1 ). On the other hand, it is larger than the presented range for the oceanic plates of Japan ( ε 2 /a = 4 × 10 −5 to 8 × 10 −4 km −1 ).
In the previous section, we have seen how the estimated peak ratios at the JMA stations are significantly larger than those of the Hi-net network. Unlike the Hinet stations, the JMA stations are positioned on the volcanic edifices, so there seems to be a correlation between scattering strength and the proximity to the volcanic centre. To examine this further, we plot the peak ratio against station's distance to the volcanic crater or summit according to the 4th edition National Catalogue of the Active Volcanoes in Japan (Japan Meteorological Agency 2013) at the 2-4 Hz band (Fig. 4). Even though we do see small peak ratios close to the summit as well, for borehole stations, larger peak ratios of > 0.3 are confined at shorter distances within 5 km or so. Furthermore, as we go further from the summit, the peak ratio ostensibly decreases. The continuation of this trend emerged when the surrounding Hi-net stations' peak ratios are also considered. Therefore, the observed stronger scattering can be interpreted as the intensification of small-scale medium heterogeneity as we get closer to the volcanic centre, for example, as a result of magmatism or hydrothermal activities. We could also see a similar decline for the ground surface stations in Fig. 4, but it is less prominent than for the borehole stations.
We observed stronger scattering (larger peak ratios) at the JMA ground surface stations in comparison to their borehole counterparts (Figs. 3, 4). This may reflect the stronger heterogeneity at the near surface that may consist of pyroclastic materials and/or lava flow deposits. Similar accord is also suggested by Hirose et al. (2019) who used passive ambient noise cross-correlation functions to estimate the scattering parameters of Rayleigh waves at Sakurajima volcano. The technique has been validated to be consistent with estimations made using active shot experiments, which are often applied at active volcanoes (e.g., Asama volcano, Yamamoto and Sato 2010). The strong scattering parameters (short mean free path) inferred by these approaches suggest intensive medium heterogeneities in the shallow parts of active volcanoes. This is because such methods mostly sample the shallow depths, and their results indicate stronger scattering (shorter mean free path) than what is inferred Table 1 Mean and standard deviation of peak amplitude ratios and transverse peak delay times N is the number of stations available. Some stations are unclassified on the JMA website. Hi-net numbers are from Table 1 of Nishimura (2012) We divide the data set into three categories (all JMA, ground surface, and borehole stations), as well as into different frequency bands Hi-net all data 0.5-1 0.14 ± 0.11 10.7 ± 6.7 635 1-2 0.14 ± 0.10 8.6 ± 5.2 725 2-4 0.16 ± 0.11 7.0 ± 4.7 712 4-8 0.20 ± 0.14 6.5 ± 4.6 574 for the typical Earth's lithosphere (e.g., Sato et al. 2012). Our observation of the peak ratio contrast between the ground surface and borehole stations builds upon this by suggesting that scattering may have been enhanced at the very shallow depths, in the first few hundred metres or so following the installation depths of borehole stations. This observation also compels the careful examination of regional peak ratios by considering the seismometer installation depth. As we have established, assuming the Gaussiantype media, the peak ratio is proportional to ε 2 and a . Kubanza et al. (2007) following the forward scattering model of Sato (2006) also suggested that for the case of the vertical incidence plane P wavelet, the peak delay of the horizontal component intensity spectral density (MS envelope) is scaled by the characteristic time t M , which in the case of the Gaussian-type media increases in proportion to the ratio of ε 2 and a as well. In addition, Takemura et al. (2015) numerically performed the finite difference method (FDM) of wave propagation simulation in a von Karman type random medium, and they observed that scattering off the small-scale heterogeneities can progressively accumulate as the wave travels, which would partition more energies into the transverse component. Therefore, one would think of a positive correlation between the peak ratio and the transverse peak delay time. However, interestingly this relationship is not observed in our results (Table 1, as well as Fig. 3 and Additional file 1: Fig. S8). And so, there might be additional factors at play that were not included in the Sato (2006) model in regard to the observed scattering strength. For example, the effect of the free surface (e.g., Emoto et al. 2010) might be an influencing factor. Another factor might be the irregular topography of the volcanoes which may trigger P-to-S and body-to-surface wave conversions (e.g., Rayleigh wave) that may also excite the transverse component. We also need to consider systematic errors in the selection of envelope peaks especially for the transverse component. The shape of the transverse envelope often develops more gradually unlike the clear and early peak of the three-component sum envelope. In addition to the fluctuating nature of the Fig. 4 Peak amplitude ratio plotted against the estimated distance to the volcanic centre according to the 4th edition National Catalogue of the Active Volcanoes in Japan (Japan Meteorological Agency 2013). The plots are for the JMA stations at the frequency band of 2-4 Hz. JMA stations can be differentiated into ground surface (red) and borehole (blue) stations. Stations with green symbols are not classified on the JMA website. Also shown are peak ratios at Hi-net borehole stations that are within 20 km of volcanoes (Nishimura 2012). Overall, the JMA ground surface stations show larger peak ratios than their borehole counterparts. For the JMA borehole stations, we see an ostensible decay in peak ratio with distance. This pattern is continued by the Hi-net stations envelope, this could result in more than one local peaks with similar amplitudes but significantly different time delays. Accordingly, the effect to the peak ratio estimation is minimal, but it may introduce ambiguity to the peak delay time determination. A clearer understanding bridging the observed and theoretical relationship between peak ratio and peak delay time can perhaps help alleviate the above idiosyncrasies and be an avenue for future research.
Additional file 1: Figure S1. Distribution of teleseismic sources used in this study. Figure S2. Peak amplitude ratio at the JMA seismometers averaged for individual volcanoes at all frequency bands evaluated in this study. The colour scale is the same as the one used in Nishimura (2012). Figure S3. The distribution of peak amplitude ratio calculated for the JMA volcanic seismometer network averaged for each volcano at the 0.5-1 Hz frequency band across different regions: (a) Tohoku (b) Hokkaido (c) Kyushu (d) Kanto-Chubu. The colour scheme used mimics the one from Nishimura (2012). Figure S4. The distribution of peak amplitude ratio calculated for the JMA volcanic seismometer network averaged for each volcano at the 1-2 Hz frequency band across different regions: (a) Tohoku (b) Hokkaido (c) Kyushu (d) Kanto-Chubu. The colour scheme used mimics the one from Nishimura (2012). Figure S5. The distribution of peak amplitude ratio calculated for the JMA volcanic seismometer network averaged for each volcano at the 4-8 Hz frequency band across different regions: (a) Tohoku (b) Hokkaido (c) Kyushu (d) Kanto-Chubu. The colour scheme used mimics the one from Nishimura (2012). Figure S6. Station peak amplitude ratios plotted against their estimated distance to the volcanic centre according to the 4th edition National Catalogue of the Active Volcanoes in Japan (Japan Meteorological Agency, 2013) at all four frequency bands. JMA stations can be differentiated into ground surface (red) and borehole (blue) stations. Stations with green symbols are not classified on the JMA website. Also shown are peak ratios at Hi-net borehole stations that are within 20 km of volcanoes (Nishimura, 2012). Overall, the JMA ground surface stations show larger ratios than their borehole counterpart. For the JMA borehole stations especially at the lower frequencies, we see an ostensible decay in peak ratio with distance. This pattern is continued by the Hi-net stations. Figure S7. Peak amplitude ratios according to the volcano where the stations are positioned. Stations are further differentiated into ground surface (red) and borehole (blue) stations. Stations with green symbols are not classified on the JMA website. Results are presented for the four frequency bands used in this study. From left to right, volcanoes are sorted from South to North. Grey denotes the 95% confidence interval for the peak ratio of quaternary volcano group cited in Nishimura (2012). For the 4-8 Hz band, there are no data for the Asosan and Yakedake volcanoes. Figure S8. Probability density functions of transverse peak delay times at four frequency bands: 0.5-1, 1-2, 2-4, 4-8 Hz. Yellow histograms show the distributions of transverse peak delay times when all JMA stations are considered. Red and blue histograms separate the transverse peak delay time distributions into JMA ground surface and borehole stations, respectively. The grey bands indicate the transverse peak delay time mean plus-and-minus one standard deviation for all Hi-net data according to Table 1 of Nishimura (2012).