# Anomalous Attenuation of High-Frequency Seismic Waves in Taiwan: Observation, Model and Interpretation

Marie Calvet, Ludovic Margerin, Shu-Huei Hung

► To cite this version:

Marie Calvet, Ludovic Margerin, Shu-Huei Hung. Anomalous Attenuation of High-Frequency Seismic Waves in Taiwan: Observation, Model and Interpretation. Journal of Geophysical Research: Solid Earth, 2023, 128 (3), pp.e2022JB025211. 10.1029/2022jb025211. hal-04008752

HAL Id: hal-04008752

https://hal.science/hal-04008752

Submitted on 28 Feb 2023

HAL is a multi-disciplinary open access
archive for the deposit and dissemination of sci-
entific research documents, whether they are pub-
lished or not. The documents may come from
teaching and research institutions in France or
abroad, or from public or private research centers.

L'archive ouverte pluridisciplinaire HAL, est
destinée au dépôt et à la diffusion de documents
scientifiques de niveau recherche, publiés ou non,
émanant des établissements d'enseignement et de
recherche français ou étrangers, des laboratoires
publics ou privés.

<figure><img src="image_1.png" /></figure>

Distributed under a Creative Commons Attribution 4.0 International License
JGR Solid Earth

RESEARCH ARTICLE

10.1029/2022JB025211

**Key Points:**

* We map crustal seismic attenuation parameters across Taiwan in the 1–8 Hz band with a spatial resolution of 30 km

* Scattering loss varies over more than one order of magnitude across Taiwan; variations of absorption are of the order of 30%

* Anomalous energy distribution in the coda observed in porous sedimentary and highly faulted areas is explained by the presence of fluids

**Supporting Information:**

Supporting Information may be found in the online version of this article.

**Correspondence to:**

M. Calvet,  
marie.calvet@irap.omp.eu

**Citation:**

Calvet, M., Margerin, L., & Hung, S.-H. (2023). Anomalous attenuation of high-frequency seismic waves in Taiwan: Observation, model and interpretation. *Journal of Geophysical Research: Solid Earth, 128*, e2022JB025211. https://doi.org/10.1029/2022JB025211

Received 25 JUL 2022
Accepted 2 FEB 2023

**Author Contributions:**

**Conceptualization:** M. Calvet, L. Margerin

**Formal analysis:** M. Calvet, L. Margerin

**Funding acquisition:** L. Margerin, S.-H. Hung

**Investigation:** M. Calvet

**Methodology:** M. Calvet, L. Margerin

**Project Administration:** L. Margerin

**Resources:** S.-H. Hung

**Software:** M. Calvet, L. Margerin

**Supervision:** S.-H. Hung

**Validation:** M. Calvet, L. Margerin, S.-H. Hung

**Writing – original draft:** M. Calvet, L. Margerin, S.-H. Hung

© 2023. The Authors.

This is an open access article under
the terms of the Creative Commons
Attribution License, which permits use,
distribution and reproduction in any
medium, provided the original work is
properly cited.

# Anomalous Attenuation of High-Frequency Seismic Waves in Taiwan: Observation, Model and Interpretation

<figure><img src="image_2.png" /></figure>

M. Calvet¹, L. Margerin¹, and S.-H. Hung²

¹CNRS, Institut de Recherche en Astrophysique et Planétologie, Université de Toulouse 3, Toulouse, France, ²Department of Geosciences, National Taiwan University, Taipei, Taiwan

**Abstract** High resolution maps of seismic attenuation parameters in Taiwan have been obtained from a modified “Multiple Lapse Time Window Analysis” (MLTWA). At most of the stations in porous sedimentary and highly faulted areas in Taiwan, the conventional modeling of MLTWA based on the scalar theory of radiative transfer in a half-space with isotropic scattering fails to explain the spatio-temporal distribution of the whole S-wave train. Using Monte Carlo simulations of wave transport, we demonstrate that this anomalous energy distribution in the coda may be modeled by multiple anisotropic scattering of seismic waves. In addition to the scattering quality factor $Q_{sc}$, we introduce a parameter $g$ (independent of $Q_{sc}$) which determines the angular redistribution of energy upon scattering (scattering anisotropy). We inverted the attenuation parameters $Q_{sc}^{-1}$, $Q_{l}^{-1}$ and $g$ in three frequency bands (1–2, 2–4, and 4–8 Hz). Overall, Taiwan is more attenuating than most orogens with a mean effective scattering loss $(Q_{sc}^{*})^{-1} = Q_{sc}^{-1}(1-g)$ and a mean intrinsic absorption $Q_{l}^{-1}$ of $2.5 \times 10^{-3}$ and $9 \times 10^{-3}$ at 1.5 Hz, respectively. Scattering loss $(Q_{sc}^{*})^{-1}$ varies over more than one order of magnitude across Taiwan while absorption variations reach approximately 30%. The more attenuating zones are the Coastal Range and the Coastal Plain where scattering dominates over absorption at low frequency, and inversely at high frequency. These regions are also characterized by strong backscattering ($g < -0.85$) at 1.5 Hz and rather high $V_P/V_S$ ratio. We propose that the observed strong back-scattering at low frequency, related to large impedance fluctuations in the crust, is induced by the presence of fluids.

**Plain Language Summary** Seismic attenuation is a key parameter which controls the amplitude decay of ground motions with the distance from an earthquake. Scattering and absorption are the two mechanisms controlling attenuation. The former is caused by the presence of small-scale geological heterogeneities in earth's crust. The later is caused by anelastic processes related to fluids. The analysis of the distribution in space and time of the seismic energy allows us to separate the contribution of the two mechanisms. We implement a new propagation model which takes into account the nature of heterogeneities. This model has been used to map seismic attenuation across Taiwan. Overall, this orogen is more attenuating than most others worldwide. Spatial variations in scattering attenuation reach about one order of magnitude across the Island. By contrast, absorption variations are of the order of 30%. High tectonic or volcanic activity is correlated with strong scattering and absorption, low seismic wavespeeds and high velocity ratios between seismic compressional and shear waves. In the same areas, our model also reveals that seismic energy is strongly backscattered, hence tends to diffuse out more slowly from the source region. The combination of these observations may be indicative of a high concentration of fluids/melts.

# 1. Introduction

The Taiwan orogen is located at the junction of two opposite-verging trenches: to the northeast the Philippine Sea plate (PSP) subducts northward beneath the Eurasian plate (EP) along the Ryukyu trench while to the south the EP subducts eastward beneath the PSP along the Manila trench (Figure 1). Its formation has been widely thought to result from the oblique arc-continent collision between the EP and the PSP (Teng, 1990), presently with a convergence rate of about 82 mm/yr in a direction close to SE-NW (Yu et al., 1997). The Longitudinal Valley (LV) in eastern Taiwan is generally marked as the collision suture between the two plates that accommodates nearly half the modern convergence rate (Figure 1). The orogen is still young and very active as manifested by high rates of uplift/erosion (Ching et al., 2011; Dadson et al., 2003) and frequent and dense occurrence of earthquakes (Y.-L. Chen et al., 2016).
<figure><img src="image_3.png" /><figcaption>**Figure 1.** Topographic and tectonic map of Taiwan which mainly results from the oblique convergence between the Philippine Sea Plate (PSP) and the Eurasian Plate (EP). The white arrow indicates the velocity vector of the PSP relative to the EP. The boundaries of plate convergence along the Manila Trench in the south and the Ryukyu Trench in the northeast of Taiwan are shown in white with arrowheads indicating opposite subduction polarities. The dashed line indicates a 2000-m isopach of the Cenozoic sediments from A. Lin et al. (2003). Surface traces of active faults from Shyu et al. (2016) are shown in red lines. The geology of Taiwan is characterized by five major lithotectonic units (delineated with thick dark lines): from west to east, the Coastal Plain (CP), the Western Foothills (WF), the Hsuehshan Range (HR), the Central Range (CR), and the Coastal Range (CoR). The deformation front (DF), also the boundary between the CP and WF denoted by a dashed white line, extending northward from the northern Manila trench offshore SW Taiwan to the northern coast of Taiwan marks the frontal limit of collision-related deformation. The longitudinal Valley (LV) corresponds to the collision suture zone between the EP and PSP. The Okinawa Trough (OT) behind the Ryukyu arc-trench system is a back-arc rifting basin which migrates southwestward to the Ilan Plain (IP) in northeast Taiwan. The map also shows the seismological stations used in the present study (left panel) and the distribution of earthquakes (right panel) which occurred between 1994 and 2016 with a focal depth shallower than 40 km and a local magnitude greater than 3.5 reported from the Central Weather Bureau of Taiwan.</figcaption></figure>

The geological style of Taiwan is essentially a concatenation of lithotectonically distinct terranes trending approximately north-northeast separated by major boundary faults. It is mainly divided into five units from west to east: the Coastal Plain (CP), the Western Foothills (WF), the Hsuehshan Range (HR), the Central Range (CR), and the Coastal Range (CoR) (Ho, 1988) (Figure 1). Prior to the arc-continent collision, all of these units except for the CoR belonged to the rifted Eurasian continental margin floored with flat-bedded, up to 10-km thick Cenozoic sedimentary strata unconformably overlying the pre-Cenozoic basement (Teng & Lin, 2004). The CP and the offshore Taiwan Strait sitting in the inner part of the continental margin are not involved in the orogeny, whereas the outer margin encompassing the WF, HR, and CR has been uplifted and deformed into the Taiwan orogenic belt. The deformation front bounds the CP and WF and represents the westernmost limit of the orogenic contractional zone. The CP, as the undeformed onshore part of the foreland basin, is covered with unconsolidated alluvial deposits eroded from the uplifting mountain ranges to the east. The WF comprises a west-verging fold-and-thrust belt loaded by the Cenozoic siliciclastic deposition (Teng, 1990). The HR and western part of the CR (also named the Backbone Range (BR)) make up the Taiwan Slate Belt composed of the Eocene-Miocene low-grade metasedimentary strata. The Eastern Central Range represents the oldest outcropping lithostratigraphic unit of Taiwan that consists of the pre-Tertiary, high-grade poly-metamorphic complex with an assemblage of meta-sedimentary and meta-magmatic rocks known as the Tananao Schist or Basement Complex (TC) (W.-S. Chen et al., 2017). The exotic CoR formed by the accreted Luzon arc and forearc basin of the PSP affinity is comprised of the Neogene andesitic volcaniclastic and flysch sedimentary rocks and accretionary
mélange complexes (Hsu, 1976). On the other hand, the northern Taiwan orogenic belt has transitioned from compressional collision into extensional collapse since the Plio-Pleistocene times (Teng, 1996). The lithospheric stretching has induced the post-collisional magmatism to produce the Northern Taiwan Volcanic Zone (NTVZ) including onshore Tatun Volcano Group (TVG) and Keelung Volcano Group (KVG) at the northern tip of Taiwan and several offshore volcanoes composed dominantly of andesitic lavas and pyroclastic rocks (K.-L. Wang et al., 1999) (Figure 1). Meanwhile, it has reactivated the ceased back-arc rifting of the Okinawa Trough (OT) and facilitated its southwestward extension onshore to the Ilan Plain (IP) in northeast Taiwan.

Over the past few decades, crustal elastic wave velocities of Taiwan have been intensively explored largely through tomographic inversion or forward modeling of massive *P* and *S* arrival times from local earthquakes and active sources (Huang et al., 2014; Y.-P. Lin et al., 2011). Though the resulting distributions of *P*-wave velocity (*V*<sub>P</sub>), *S*-wave velocity (*V*<sub>S</sub>), and their ratio (*V*<sub>P</sub>/*V*<sub>S</sub>) have illuminated our understanding of the subsurface structure and tectonic evolution of the Taiwan orogen as well as the regional characterization of seismic ground motion and risk, they still remain insufficient to explain all seismic waveform signals and to infer physical properties and bulk compositions of crustal substances. As evidenced by ubiquitous coda trailing the direct ballistic waves and its exponential decay of amplitude, the crust of Taiwan being created, deformed, and destroyed through diverse geological processes exhibits the highly heterogeneous and anelastic nature. It is well known that they both contribute to seismic wave attenuation but its character and physical origin have not yet been fully investigated.

Previous attenuation models of the crust and upper mantle in the Taiwan region are mostly constrained from fitting the observed displacement spectra of direct *P*- and *S*-waves of local earthquakes at a single station or amplitude ratios at paired stations as a function of frequency to those predicted by an assumed Brune-type ($\omega^{-2}$) source model (Brune, 1970). The corresponding corner frequency ($f_c$) indicative of the finite dimension of an earthquake at which its source spectrum starts to drop abruptly, and the attenuation operator, $t^*$, defined as the integral of the reciprocal of quality factor $Q(Q^{-1})$ along the path of travel, are then simultaneously determined for inferring the earthquake stress drop and dimension (Ko, Kuo, & Hung, 2012) and spatial distributions of $Q_P$ and $Q_S$ (Roth et al., 1999), respectively. Assuming frequency-independent quality factors, Y.-J. Wang et al. (2010) presented 3-D $Q_P$ and $Q_S$ models of the crust beneath the Taiwan orogenic belt from the tomographic inversion of $t^*$ values measured at high frequencies of 2–30 Hz for *P* and 1–20 Hz for *S*. However, due to the lack of crossing ray paths and thus poor tomographic resolution in the near surface, there is neither obvious correlation found between the lateral variations of *Q* and properties of rock outcrops nor sharp change of attenuation across major geological boundaries or structural lineaments. The earlier attenuation models from the inversion of fewer $t^*$ data estimated directly from the exponential decay of *P*- and *S*-wave amplitude spectra with frequency declared that the regions of lower $Q_P$ and $Q_S$ correlate with higher seismicity in the upper crust (<15 km depth), whereas it is opposite in the lower crust and upper mantle (K. J. Chen, 1998; K.-J. Chen et al., 1996). Cheng (2013) inverted the maximum acceleration of *S*-waves converted from seismic intensity reported by the Central Weather Bureau (CWB) of Taiwan for the $Q_S$ structure, indicating the anomalously low $Q_S$ regions extending from the upper crust to the lithospheric mantle at 75 km depth beneath the offshore volcanic island, Kueishantao, IP and NTVZ.

Despite the laboratory and field evidence for frequency-dependence of attenuation (Mavko et al., 2020), all the aforementioned studies assumed frequency-independent *Q*, overlooked the tradeoff problem between the estimated corner frequency ($f_c$) and $t^*$, neglected the possible scattering loss causing the amplitude decay of high-frequency body waves, and implicitly attributed all the resulting *Q* variations to the anelastic properties of the crust and upper mantle varied by temperature, fluid/melt content, population of fractures or cracks, and etc. Though these models reveal some common features, such as relatively lower $Q_P$ and $Q_S$ in the volcanic areas of northern Taiwan, the alluvial deposit plain of southwestern Taiwan, and the high heat flow area under the SE CR, there still exist discrepancies as to details among them. For example, both Y.-J. Wang et al. (2010) and Cheng (2013) observed the significant lateral variation of attenuation across the CLF in their $Q_S$ models. However, compared to $Q_S$ imaged in the foot-wall block of the CLF at depths $\le 10$ km, the former yielded relatively higher $Q_S$ anomalies in the hanging-wall block between the CLF and Chuchih fault (CF) separating the WF and HR, while the latter resulted in comparatively lower $Q_S$, stronger attenuation therein. Ko, Kuo, and Hung (2012) introduced a cluster event method that substituted a conventionally single event with a cluster of nearby events pairing with stations for robust determination of magnitude-dependent $f_c$ and path-dependent $t^*$. Ko, Kuo, Wang, et al. (2012) applied it along with two-station spectral ratio measurements to imaging the $Q_P$ variation across northeast Taiwan and the offshore southwestern Ryukyu subduction zone. Unlike the broadly
extended low $Q_p$ and $Q_s$ zones under north and northeast Taiwan observed in Cheng (2013), a sharp, narrow zone of high $Q_p$ anomaly which immediately abuts a low $Q_p$ region beneath Kueishantao was found at tens to 90 km depths.

On the other hand, Chung et al. (2009) mapped the lateral variations of the coda quality factor $Q_c$ at several narrow frequency bands between 1.5 and 18 Hz in Taiwan from local S-wave coda starting at twice the travel time of direct S. The results indicate higher $Q_c$ in the western plain and eastern coastal areas and lower in the central mountain range at the central frequencies of 1.5 and 3 Hz but an almost reverse trend as the frequency increases. However, as demonstrated in Calvet and Margerin (2013), the measured $Q_c$ strongly depending on the chosen lapse-time window which could lead to the biased variation of $Q_c$ with earthquake hypocentral distance has not been cautiously investigated in this study. Zhang and Papageorgiou (2010) made a first attempt to separate the contribution of intrinsic absorption ($Q_i$) and scattering loss ($Q_{sc}$) to total apparent S-wave attenuation ($Q_{sc}$) by adopting a multiple lapse time window analysis (MLTWA) (Fehler et al., 1992) and found the obtained $Q_c$ in accord with S-wave attenuation estimated by the coda normalization method. However, they only determined intrinsic and scattering $Q$ at few sparsely distributed stations under the assumption of isotropic scattering. The results remain too preliminary and uncertain to shed light on the geological structure and bulk crustal properties in Taiwan.

In this work, we shall pursue the efforts put into the mapping of seismic attenuation across Taiwan with a special emphasis on the crust. While our approach is based on the classical MLTWA, our observations highlight some unusual spatial distributions of S-wave energy measured in a time window comprising the ballistic S-wave and its early coda. To interpret this original behavior, we introduce a non-isotropic scattering model which offers a simple physical explanation for the observations. Hence, in addition to the classical quality factors we are led to consider an anisotropy parameter $g$ that quantifies the tendency of seismic waves to scatter backward ($g < 0$) or forward ($g > 0$). We examine both the frequency dependence and the spatial variations of the three attenuation parameters across Taiwan. In addition, we investigate the relation between the distribution of anomalous stations where non-isotropic is required to explain the data and salient geological and seismological features highlighted by previous investigations.

Our work has been made possible by the efforts of the CWB of Taiwan which, since the 1990s, has established the CWB Seismic Network (CWBSN) for monitoring regional seismic activities and routinely locating earthquakes hypocenters with manually-picked P and S arrival times. Beginning in 2012, the short-period, broadband, and strong-motion seismograph networks on more than 200 sites distributed in Taiwan and offshore islands are all equipped with high resolution 24-bit recorders and integrated to a unified data platform that has significantly increased the detection capability of local seismicity (Shin et al., 2013) and lateral resolution of crustal velocity heterogeneity across Taiwan (Huang et al., 2014). Owing to recent advances in theoretical formulation of seismic energy transport over a broad range of propagation regime and coda wave sensitivities to scattering and absorption quality factors in the anisotropically scattering crust (Margerin et al., 2016; Mayor et al., 2014), the vast amount of earthquake data set would open up new opportunities to comprehensively reassess the scattering and absorption properties of the crust in Taiwan.

## 2. Data and Observations

In this study, we make a first selection of earthquakes with a focal depth less than 40 km that occurred between 1994 and 2016. A total of 17,149 events with a local magnitude greater than 3.5 are extracted from a catalog provided by the CWB (Figure 1). The maximal magnitude in our data set is 7.07. Waveforms are collected on broadband, short-period and accelerometric networks operated by the CWBSN and BATS (Broadband Array in Taiwan for Seismology, Institute of Earth Science, Academia Sinica). A total of 236 stations are used (see Figure 1). We select data verifying the following criteria in the frequency band 2-4 Hz: (a) the hypocentral distance is smaller than 100 km (b) S-wave coda duration is greater than 75 s with a signal-to-noise ratio (SN) greater than 4; (c) the average intensity decreases monotonically in the coda window from the maximum of the envelop. This last criterion eliminates undesired signals such as aftershocks in the coda. Our final data set contains about 330,000 high-quality waveforms.

Next, the spatio-temporal distribution of the energy of the whole S-wave train is analyzed using the Multiple Lapse Time Window method developed by Fehler et al. (1992) and Hoshiba (1993). At every station,
we collect events with a hypocentral distance $r_j$ less than 60 km and a focal depth smaller than 15 km. The recorded waveforms are first filtered in three frequency bands ([1–2], [2–4], and [4–8] Hz) using bandpass Butterworth filters of order 3. Contrary to previous studies (Carcolé & Sato, 2010; Hoshiba, 1995), we were unable to exploit higher frequencies. The primary reason is the severe degradation of the signal-to-noise ratio in the S-wave coda, which entails the rejection of the vast majority of waveforms according to the aforementioned selection criterion (2).

Following the MLTWA procedure introduced by Fehler et al. (1992), we compute the instantaneous energy of the signals by summing the squares of the ground velocities on three components. Next, energy integrals denoted by $U_i$ ($i$ = 1, 2, 3) are calculated in three consecutive time windows of 15 s duration. The first window starts 1 s before the S-wave arrival time $t_S$ in order to properly take into account the whole energy of the ballistic wave. To avoid contamination by the P-wave coda, the signal is set to zero for times smaller than ($t_S$ − 1) s. To correct for the source magnitude and local site amplification, we compute the fourth energy integral, $U_4$, in a 15 s duration time window starting 60 s after the origin time of the earthquake. We subsequently normalize the energy in the windows $U_i$ ($i$ = 1, 2, 3) by the energy in window $U_4$, and correct for the geometrical spreading with a factor $4\pi r_j^2$. This last operation compensates for the divergence of $U_1$ as the hypocentral distance approaches 0. In spite of all reasonable care taken to select the data and remove source and site effects, our measurements are sometimes contaminated by outliers. To circumvent this difficulty, we employ a simple denoising technique to the data: (a) energy ratios are sorted as a function of hypocentral distance; (b) the trimmed mean of the energy ratios is computed in moving a window of length seven samples and a stride of one sample. In effect, the minimum and the maximal values are removed before the mean is computed. Because the window is relatively short, the procedure does not introduce too much spatial smearing. Furthermore, stations with less than 20 measurements are systematically discarded.

Figure 2 shows, in a semi-logarithmic scale, the normalized energy ratios, $4\pi r_j^2 U_i / U_4$, as a function of hypocentral distance at two stations located in distinct tectonic units in Taiwan: CHK in the CoR and STY in the CR. At low frequency (1–2 Hz), the energy distributions at CHK and STY are clearly distinct. STY exhibits a “classical” behavior as previously reported in a number of studies based on MLTWA from the literature (see, as an archetypal example of MLTWA, the measurements in Japan by Carcolé and Sato (2010)): (a) at small hypocentral distance $r$, the energy ratio between the first and the second window is greater than the energy ratio between the second and the third window, so that we have, in general, $\log_{10}(U_1/U_2) > \log_{10}(U_2/U_3)$; (b) The normalized energy ratio in the first window $U_1$ decreases with hypocentral distance, whereas the ratios $U_{2,3}$ show the opposite behavior. By contrast, CHK exhibits an “anomalous” behavior: (a) the energy is much more “uniformly” distributed in the coda: the energy ratio between the first and the second window is close to the energy ratio between the second and the third window, namely, $\log_{10}(U_1/U_2) \sim \log_{10}(U_2/U_3)$; (b) the three energy ratios $\log_{10}(4\pi r_j^2 U_i / U_4)$ show a clear curvature with a positive slope at small hypocentral distance. Most of the stations located in the WF and the CoR exhibit an “anomalous” energy distribution similar to CHK. At higher frequency (>4 Hz), we note that most stations exhibit the more “classical” energy distribution in the coda, indicating a change in the propagation regime.

### 3. Modeling the Anomalous Distribution of Seismic Energy

Conventional modeling of MLTWA is based on the scalar theory of radiative transfer (RT) in a half-space containing isotropic scatterers. In this model, the spatial distribution of seismic energy around frequency $f$ depends on the scattering and absorption quality factors $Q_{sc}(f)$ and $Q_i(f)$ only. To illustrate the unusual distribution of seismic energy at CHK, we show in Figure 3 the a posteriori fits of MLTWA measurements obtained with the conventional isotropic multiple-scattering model (gray lines). This model provides an excellent match to the data at station STY but fails at capturing the curvature of the spatial energy distribution near the source at CHK.

Several authors have previously argued that non-isotropic scattering is important to model the envelope shape of seismograms near the ballistic arrivals (Calvet & Margerin, 2013; Gaebler et al., 2015; Gusev & Abubakirov, 1996; Heller et al., 2022; Hoshiba, 1995; Zeng, 2017). Hence, it is natural to hypothesize that the assumption of isotropic scattering may be at the origin of the discrepancy between the observations and predictions of the RT model at station CHK. To confirm the plausibility of this argument, we have incorporated a non-isotropic scattering pattern into Monte-Carlo simulations of the multiple-scattering process, following previous works of Calvet and Margerin (2013). We note that for elastic waves, both forward and backward scattering are possible depend-
<figure><img src="image_4.png" /><figcaption>**Figure 2.** Examples of multiple lapse time window analysis measurements in three frequency bands [1–2], [2–4], and [4–8] Hz (from top to bottom) at two stations (CHK, STY) located in two distinct tectonic regions of Taiwan (map). Examples of records and enveloppes are shown in Figure S1 in Supporting Information S1. Blue, green and red dots correspond to the logarithm of the energy ratio in the first, second and third window, respectively, corrected for the geometrical spreading $4\pi r^2$, where *r* is the hypocentral distance. In the map, surface traces of active faults from Shyu et al. (2016) are shown in red lines and the major lithotectonic units are delineated with thick dark lines.</figcaption></figure>

<figure><img src="image_5.png" /><figcaption>**Figure 3.** Fit of multiple lapse time window analysis measurements (dots) at stations STY (left) and CHK (right) using isotropic (dashed line) or anisotropic scattering models (solid line). Scattering anisotropy (backscattering) clearly improves the agreement with data at station CHK.</figcaption></figure>
<figure><img src="image_6.png" /><figcaption>**Figure 4.** Effect of scattering anisotropy on the spatial distribution of energy contained in the early coda. The polar plot shows the scattering patterns for different values of the anisotropy parameter g, as indicated in inset, and the corresponding normalized energy ratios are shown on the right. Note that a logarithmic scale is used for the polar plot.</figcaption></figure>

ing on the nature of the heterogeneities and the frequency regime (Wu & Aki, 1985). A relatively simple and realistic analytical model of scattering non-isotropy which handles the two situations is provided by the so-called Henyey-Greenstein function (Henyey & Greenstein, 1941):

$$p(\mathbf{k}, \mathbf{k}'; g) = \frac{1 - g^2}{4\pi(1 - 2g\mathbf{k} \cdot \mathbf{k}' + g^2)^{3/2}} \quad (1)$$

In Equation 1, the unit vectors k' and k represent the directions of propagation before and after scattering, respectively. The parameter g is the mean value of the cosine of the scattering angle and bounded between -1 and 1. Negative (resp. positive) values of g correspond to preferentially backward (resp. forward) scattering, while isotropic scattering is obtained by setting g = 0. As shown by Margerin and Nolet (2003), the Henyey-Greenstein model may be understood as a limiting case of scattering in an inhomogeneous random medium characterized by a Von-Kármán power spectrum, when the Hurst exponent κ tends to 0. Physically, κ controls the roll-off of the power spectrum at large wavenumbers so that κ → 0 corresponds to random media that are rich in small-scale fluctuations. We refer the reader to Chapter 2 in Sato et al. (2012) for an in-depth discussion of heterogeneity power spectra.

The impact of non-isotropic scattering on the spatial distribution of seismic energy contained in the time window of ballistic waves + the first 15 s of the coda (referred to as the first time window hereafter) is illustrated in Figure 4. These simulations are carried out for scattering anisotropy that varies from predominantly forward (g = 0.85) to strongly backward (g = -0.95). The isotropic case (g = 0) is also shown for reference. The values of intrinsic and scattering attenuation, as quantified by inverse quality factors, are fixed at the values that are typical of station CHK as summarized in Table 2 (isotropic model). It is customary to quantify the scattering attenuation with a length scale l termed scattering mean free path and defined as:

$$l = \frac{\nu Q_{sc}}{2\pi f}, \quad (2)$$

with v the S-wave velocity. It represents the attenuation length of coherent waves or the typical distance traveled freely by the waves between two scattering events.

In Figure 4, we observe that the main effect of forward scattering is to increase uniformly the normalized energy in the first time window, as compared to the isotropic case. To understand at least qualitatively this effect, it is necessary to introduce yet another quantity known as the transport mean free path and defined as:

$$l^* = l/(1 - g). \quad (3)$$

The transport mean free path is a crucial length scale as it governs the energy level in the late coda. Indeed, from first principles of multiple-scattering theory, it can be shown that the transport of coda wave energy is governed by a diffusion equation at long lapse-time, with a diffusion constant given by $D = \nu l^*/3$ (in 3-D space) (see Sato (2019) for a review). Classical solutions of the diffusion equation in full-space indicate that the coda
excitation is proportional to $D^{-3/2}$. Furthermore it has been proposed that $l^*$ is a better length scale than $l$ to quantify the attenuation of ballistic waves, as it takes into account the fact that only the energy scattered at large angles contributes to the attenuation of direct waves observed on seismograms (Sato, 1982; Wu, 1982). Hence, for a given value of the scattering mean free path, forward scattering simultaneously increases the energy contained in the ballistic waves and reduces the coda level at long lapse-time. When combined, these two effects explain the strong increase of the normalized energy contained in the first time window in the case of forward scattering, as compared to the isotropic case.

Following the same line of reasoning as in the forward-scattering case, we expect backward scattering to simultaneously increase the attenuation of ballistic waves and the coda level at long lapse-time. As anticipated, we observe in Figure 4 that backscattering tends to reduce the normalized energy in the first time window. Furthermore, backscattering changes the curvature of the curve representing the energy in the first time window as a function of hypocentral distance from straight to concave. This effect is also clearly observed in the data at station CHK. Hence, strong backscattering of the waves appears as a good candidate to explain the anomalous seismic energy distribution at station CHK. This hypothesis is confirmed by the examination of Figure 3, which shows the bet-fitting non-isotropic models obtained by a least-squares inversion of the data where we let $Q_{sc}^{-1}$, $Q_{l}^{-1}$ and $g$ vary (see the next section for implementation details). On the one hand, at station STY, anisotropic scattering is not necessary to explain the data and the best-fitting value of $g$ is close to 0. On the other hand, at station CHK, the introduction of strong back-scattering ($g = -0.99$) improves dramatically the fit between the multiple-scattering model and the observations. This result supports the use of the Henyey-Greenstein model to represent the scattering anisotropy of seismic waves. In the next section we present an inversion method to simultaneously retrieve absorption, scattering attenuation and anisotropy in a MLTWA approach.

## 4. Inversion Method

Following Carcolé and Sato (2010), we implement a Levenberg-Marquardt algorithm to minimize a misfit function defined as the sum of the squared differences between the base 10 logarithms of observed and modeled normalized energy ratios:

$$ \chi^2 = \sum_{j} \sum_{i=1}^{l=3} \left( \log_{10} \left( \frac{4\pi r_j^2 U_i^m}{U_4^m} \right) - \log_{10} \left( \frac{4\pi r_j^2 U_i^d}{U_4^d} \right) \right)^2, \quad (4) $$

where the index $j$ runs over all earthquakes and the superscripts $m$ and $d$ refer to model-predicted and observed data, respectively. To guide the inversion, we compute the partial logarithmic derivatives $\{\partial_g, \partial_{g^{-1}}, \partial_{Q_{sc}^{-1}}\}$ of the normalized energy integrals $(4\pi r^2 U_i / U_4)$ using the approach of Takeuchi (2016) for the last two parameters. To compute the sensitivity with respect to $g$, we use the basic ideas of the differential Monte-Carlo method (Lux & Koblinger, 2018), which allows to compute the seismic energy in two media with slightly different anisotropy parameters $g$ (unperturbed) and $g + \delta g$ (perturbed) in a single simulation. It suffices to note that for a given sequence of $n$ scattering events, the energy $U_i$ detected in the coda in the perturbed and unperturbed media are related by:

$$ U_i(g + \delta g) = U_i(g) \prod_{l=1}^{l=n} \frac{p(\mathbf{k}_{l+1}, \mathbf{k}_l; g + \delta g)}{p(\mathbf{k}_{l+1}, \mathbf{k}_l; g)}. \quad (5) $$

In Equation 5, each factor of the product is equal to the probability to be scattered from direction $\mathbf{k}_l$ to direction $\mathbf{k}_{l+1}$ in the perturbed medium, divided by the probability of the exact same event in the unperturbed medium. Each factor should therefore be interpreted as a weight which restores a posteriori the correct frequency of occurrence of the scattering from direction $\mathbf{k}_l$ to direction $\mathbf{k}_{l+1}$ in the perturbed medium. Applying a Taylor expansion to Equation 5, we arrive at the following formula for the partial derivative of the energy with respect to the anisotropy parameter $g$:

$$ \partial_g U_i(g) = U_i(g) \prod_{l=1}^{l=n} \frac{\partial_g p(\mathbf{k}_{l+1}, \mathbf{k}_l; g)}{p(\mathbf{k}_{l+1}, \mathbf{k}_l; g)}, \quad (6) $$

where the logarithmic derivative of the Henyey-Greenstein function is given by:
<figure><img src="image_7.png" /><figcaption>**Figure 5.** Illustrations of multiple lapse time window analysis measurements and inversion results in the frequency band [1–2] Hz for two example stations, STY and CHK. Results for other stations are shown in Figures S2 to S14 in Supporting Information S1. Energy ratio measurements (dots), the best-fitting isotropic (dashed lines) and anisotropic (solid lines) models are shown on the left panels, where blue, green and red dots correspond to energy ratios in the first, second and third window, respectively. The station locations are indicated on the map furthest to the right: red lines are surface traces of active faults from Shyu et al. (2016) and the five major lithotectonic units are delineated with thick dark lines. The middle panels show the misfit as a function of the scattering anisotropy (g), while the right panels indicate $1/Q_i$ (black solid lines) and $1/Q_{sc}^*$ (gray solid lines) for the best models as a function of the scattering anisotropy (g) with standard errors delineated by dashed lines.</figcaption></figure>

$$ \frac{\partial_g p(\mathbf{k}_{l+1}, \mathbf{k}_l; g)}{p(\mathbf{k}_{l+1}, \mathbf{k}_l; g)} = \frac{2g}{g^2 - 1} + \frac{3(\mathbf{k}_l \cdot \mathbf{k}_{l+1} - g)}{(1 - 2g\mathbf{k}_l \cdot \mathbf{k}_{l+1} + g^2)}. \quad (7) $$

It should be noted that in practice Equation 6 has to be averaged over a large number of realizations to obtain a sufficiently accurate estimate of both the energy and its partial derivative. Empirically, we found that the simulation of $4 \times 10^6$ independent trajectories was sufficient for our purposes.

Since the Levenberg-Marquardt algorithm is a local optimization method, it is important that the starting model be sufficiently close to the best-fitting model in the parameter space, in order to avoid a local minimum trap. We therefore adopt a two-step inversion scheme to estimate the attenuation parameters $Q_i^{-1}$, $Q_{sc}^{-1}$ and g at every station. First, $Q_{sc}^{-1}$ and $Q_i^{-1}$ are estimated for 19 fixed values of g ranging from -0.95 to 0.85 with an increment $\Delta g = 0.1$. Second, starting from the best solution obtained at step 1, the three parameters ($Q_i$, $Q_{sc}$, g) are optimized with the aid of the Levenberg-Marquardt algorithm, thanks to the implementation provided by the GNU Scientific Library (Galassi et al., 2002). The first inversion step also allows us to study the possible trade-off between the inverse quality factors and the anisotropy parameter.

Figure 5 illustrates step 1 of our inversion approach and the results obtained at stations STY and CHK, which display distinct energy distributions in the coda. A posteriori fits of the MLTWA measurements (dots) for the best anisotropic (black solid lines) and isotropic (gray dashed lines) models are shown on the left. For the simulations, the S-wave velocity is deduced from the experimental traveltime curves constructed at each station from the seismic bulletin provided by the CWT. Table 1 provides the best-fitting values of $Q_i^{-1}$, $Q_{sc}^{-1}$ and g for these

two stations. The misfit function, the effective scattering attenuation ($1/Q_{sc}^*$) and the absorption ($Q_i^{-1}$) (obtained after step 1) are also plotted as a function of g together with the typical uncertainties. Following our discussion on the transport mean free path, the effective scattering loss is defined as $1/Q_{sc}^* = Q_{sc}^{-1}(1-g)$. This quantity is thought to be a better estimate of the scattering attenuation of direct waves than $Q_{sc}^{-1}$ since it takes into account the fact that the energy scattered around the forward direction is not entirely lost. From the uncertainties shown in Figure 5, it is clear that absorption $Q_i^{-1}$ is better constrained than the scattering parameters $Q_{sc}^{-1}$ and g. We also remark that $Q_{sc}^{-1}$ is strongly correlated with g, in contrast with $Q_i^{-1}$.

<table><caption>**Table 1**
Best Isotropic and Anisotropic Models at STY and CHK</caption><thead><tr><th rowspan="2">Station</th><th colspan="3">Isotropic model</th><th colspan="3">Anisotropic model</th></tr><tr><th>1/Q<sub>i</sub></th><th>1/Q<sub>sc</sub></th><th>g</th><th>1/Q<sub>i</sub></th><th>1/Q<sub>sc</sub></th><th>g</th></tr></thead><tbody><tr><td>STY</td><td>0.00849</td><td>0.00829</td><td>0</td><td>0.00852</td><td>0.00905</td><td>0.01</td></tr><tr><td>CHK</td><td>0.00826</td><td>0.01204</td><td>0</td><td>0.01102</td><td>0.01658</td><td>-0.99</td></tr></tbody></table>
<figure><img src="image_8.png" /><figcaption>**Figure 6.** Comparison of scattering and intrinsic attenuation ($1/Q_{sc}^*$ and $1/Q_i$) between isotropic and anisotropic models in three frequency bands [1–2], [2–4], [4–8] Hz. The scattering anisotropy *g* is given by the color scale (negative values are in red and positive ones in blue). The size of the circles indicates the misfit ratio between isotropic and anisotropic models ($\chi_{iso}/\chi_{aniso}$).</figcaption></figure>

Station STY is characterized by a “classical” distribution of the seismic energy. The misfit function decreases with *g* in the range [−0.95, 0.0] and increases rather slowly for positive *g*. We do not observe a significant difference between the a posteriori fits for anisotropic models with *g* → 1 and a posteriori fit for the isotropic multiple-scattering models. In other words, all models with *g* > −0.4 can explain equally well the MLTWA measurements at STY. For *g* > −0.4, we observe on the one hand that $Q_i^{-1}$ is almost constant, regardless of the level of anisotropy. On the other hand $Q_{sc}^{-1}$ increases with *g* indicating that strong forward-scattering is compatible with the data but a higher scattering level is required.

Unlike STY, station CHK shows an “anomalous” distribution of the seismic energy, characterized by the downward curvature of the normalized energy as a function of hypocentral distance near *r* = 0 km. The corresponding misfit function increases with *g* over the whole range [−0.95, 0.85]. A posteriori fits for models with strong back-scattering ($g \approx -1$) are significantly improved, in particular at short hypocentral distances for the first time window. To summarize, an anomalous distribution of the energy in the coda is better explained by a multiple-scattering model with strong back-scattering ($g < -0.8$), while a “classical” energy distribution in the coda can be equally well explained by isotropic or forward-scattering models. In the latter case, we note an increase of scattering attenuation $Q_{sc}^{-1}$ with the anisotropy parameter *g*. The equivalence of forward-scattering and isotropic scattering models in the case of a “normal” energy distribution was previously pointed out in the context of sedimentary basins by Gaebler et al. (2015).

Results at all available stations in the three frequency bands [1, 2], [2, 4], and [4, 8] Hz are summarized in Table 4. In Figure 6, we compare intrinsic and scattering attenuation ($Q_i^{-1}$ and $Q_{sc}^*$) for isotropic ($g = 0$) and anisotropic ($g \neq 0$) models. We recall that $Q_{sc} = Q_{sc}^*$ for isotropic scattering ($g = 0$). We also compute the misfit ratio, $\chi_{iso}^2/\chi_{aniso}^2$, to compare quantitatively the goodness of the fit to the data for anisotropic and isotropic models.
<figure><img src="image_9.png" /><figcaption>**Figure 7.** Total scattering maps ($1/Q_{sc}^*$) in three frequency bands: results of multiple lapse time window analysis are plotted at the location of the station. Values higher (resp. lower) than the average one are shown in red (resp. blue). The scale indicates the maximum and the minimum values in each frequency band. Pink lines are surface traces of active faults from Shyu et al. (2016) and the five major lithotectonic units are delineated with thick dark lines (see Figure 1 for more details).</figcaption></figure>

At low frequency (1–4 Hz), the data fits of anisotropic models are systematically better ($\chi_{iso}/\chi_{aniso} > 1.5$) when the optimal value of $g$ is strongly negative (typically $-1 < g < -0.8$). In that case, the isotropic models systematically underestimate the absorption (by a factor ~1.5) and the effective scattering attenuation (by a factor 2–3) compared to the best anisotropic model with negative $g$. We also observe that at most of the stations where the best-fitting model has positive $g$, the improvement of the fit to the data is modest compared to isotropic models. The scattering attenuation ($1/Q_{sc}$) found in isotropic models is generally close to the effective scattering attenuation ($1/Q_{sc}^*$) deduced from anisotropic scattering models. This result agrees with the findings of Gaebler et al. (2015). We also find that the absorption retrieved from the data does not vary significantly between isotropic and anisotropic models with $g > 0$.

At high frequency ($f = 6$ Hz), there is a significant reduction of the number of stations at which the inversion can be performed, due to the overall strong level of attenuation. As a consequence of the degraded data coverage, the inversion of the scattering anisotropy is more delicate. Nevertheless, it is worth noting in Figure 6 that models with positive $g$ exhibit stronger effective attenuation than their isotropic counterparts ($1/Q_{sc}^*(\text{anisotropic}) > 1/Q_{sc}(\text{isotropic})$). This behavior is also visible to a lesser extent at lower frequencies for models with $g > 0$ and weak attenuation. In these models, the transport mean free path can be of the same order or even larger than the typical hypocentral distance. In that case, energy transport is likely sensitive to details of the scattering that are not encapsulated in the parameter $g$. Indeed, complete information on the scattering pattern requires the knowledge of all the moments of the cosine of the scattering angle. In summary, in the weak scattering regime, the effective scattering attenuation $1/Q_{sc}^*$ is not sufficient to describe the transport of seismic energy. Sufficiently detailed knowledge of the scattering pattern appears to be required. The inequivalence of isotropic and anisotropic models also applies to the strong backscattering regime. In the next section, the inversion results are further explored by analyzing the spatial variations of attenuation and scattering anisotropy across Taiwan.

## 5. Results: Attenuation Maps

Figures 7–10 present maps of effective scattering loss ($Q_{sc}^{-1}$), absorption ($Q_i^{-1}$), effective albedo ($B_0^* = 1/Q_{sc}^*/(1/Q_{sc}^* + 1/Q_i)$), that is, the relative contribution of scattering versus absorption, and scattering anisotropy ($g$) in the [1, 2], [2, 4], and [4, 8] Hz frequency bands. Results of MLTWA inversion (see Table 4) are assigned to the location of each station as previously done by Carcolé and Sato (2010). This simple mapping approach can be justified on the basis of the sensitivity analysis of Mayor et al. (2014) and Margerin
et al. (2016). These authors have demonstrated that the sensitivity of the intensity in the coda to attenuation parameters is greater at the location of the station and the source. As MLTWA is a single station approach, there is a redundancy of paths sampling the vicinity of the seismic station. Our results are therefore expected to be representative of the attenuation properties in the vicinity of each station. We will provide in this section a description of the most striking features observed on the various maps. In the case of absorption $Q_i^{-1}$ and scattering $Q_{sc}^{*-1}$, the color scale has been specifically adapted to better highlight the lateral variations: white color corresponds to the average attenuation for Taiwan (Table 2), while dark red and dark blue correspond to the minimum and maximum values obtained in each frequency band, respectively.

## 5.1. Scattering

Compared to values reported worldwide (Sato, 2019), the mean scattering level is globally high in Taiwan, particularly at low frequency (see Table 2). Furthermore, we note that the mean value of scattering attenuation $\overline{Q_{sc}^{-1}}$ decreases significantly with frequency, over about one order of magnitude between 1.5 and 6 Hz (Table 2). This strong frequency dependence is also observed on single-station estimates (Table 4). Figure 7 exhibits huge lateral variations of scattering, about 2 orders (resp. 1 order) of magnitude at low (resp. high) frequency. Despite the relatively small geographical area, scattering fluctuations in Taiwan are one order of magnitude greater than those observed in Japan or in Korea (Carcolé & Sato, 2010; Rachman et al., 2015). In Table 3, similar tectonic units were grouped together to compare their average level of seismic attenuation. In the 1–4 Hz frequency band, scattering is consistently strong around active faults in the CP and northern WF, as well as in alluvial deposits of the IP and Pingtung Plain (south of the WF). The LV and the CoR, composed of exotic volcanic arc rocks and associated mélange materials, are also characterized by strong scattering. In comparison, the HR, CR, and central WF exhibit relatively low scattering strength. In these regions, we find that the scattering is, respectively, 4 times and 3 times less intense in the 1–2 and 2–4 Hz frequency bands, as compared to neighboring tectonically active or volcanic regions (Table 3). In the highest frequency range (4–8 Hz), scattering is on the whole much weaker. The HR (except in the vicinity of the geological boundary fault on the eastern edge) and south-central WF are still marked by the lower values of scattering, while the eastern flank of the HR, CR and CoR are globally the most attenuating areas. Even if the coverage is sparse, active faults are consistently marked by rather low scattering at high frequency.

## 5.2. Absorption

Similar to scattering attenuation, the mean level of absorption is rather high in Taiwan and slightly decreases with frequency from $\overline{Q_i^{-1}} = 9.37 \times 10^{-3}$ at 1.5 Hz to $\overline{Q_i^{-1}} = 3.22 \times 10^{-3}$ at 6 Hz (Table 2). Figure 8 illustrates the rather strong lateral variations of absorption, with an amplitude of a few tens of % at a scale of a few tens of kilometers. Yet, it is worth pointing out that these variations are modest compared to those shown on the scattering attenuation maps. The spatial distribution of absorption depends on the frequency and the amplitude of the variations slightly decreases with frequency (Table 2). As shown on the 1–2 Hz map, stations located in the LV and CoR affiliated with the accreted Luzon arc and forearc of the PSP affinity, composed mainly of andesitic volcanoclastics, flysch deposits, and ophiolite-bearing mélange, show higher absorption than those located in the CR and HR that constitute the metamorphic slate belt and schist complex (see also Table 3). Sedimentary materials located in highly tectonized area in the CP and WF, that is, in the vicinity of active faults, exhibit on average stronger absorption than sedimentary materials in the IP or near the Peikan Basement High. On the 2–4 Hz map, we observe grossly the same features except for the IP where absorption is now almost at the same level as in the CoR/LV. The absorption signature of active faults in the WF is also less clear. Even if the spatial coverage is much more sparse at high frequency ($f > 4$ Hz), the spatial distribution of high-frequency absorption does not seem to be controlled by the surface geology: we rather observe an east-west gradient.

## 5.3. Albedo

Figure 9 shows the relative contribution of scattering and absorption, also termed albedo ($B_0^*$), in the three frequency bands [1, 2], [2, 4], and [4, 8] Hz. The albedo is an important quantity for seismic hazard assessment. For the same level of total attenuation, a small value of albedo implies that attenuation will reduce the amplitude of the whole signal. To the contrary, a large value of albedo is expected to promote the transfer of energy
<figure><img src="image_10.png" /><figcaption>**Figure 8.** Absorption maps ($1/Q_i$) in three frequency bands: results of multiple lapse time window analysis are plotted at the location of the station. Values higher (resp. smaller) than the average one are shown in red (resp. blue). The scales indicate the maximum and the minimum values in each frequency band. Pink lines are surface traces of active faults from Shyu et al. (2016) and the five major lithotectonic units are delineated with thick dark lines (see Figure 1 for more details).</figcaption></figure>

from direct to coda waves. In this case the reduction of the amplitude of the direct waves can be associated to an increase in the duration of the ground motion. On average, $B_0^*$ decreases as the frequency increases. We observe that scattering is universally dominant ($B_0^* > 0.5$) at 1–2 Hz in areas where scattering is strong around active faults distributed in the CP, WF, LV, CoR, and IP while scattering and absorption are equivalent ($B_0^* \sim 0.5$) in the HR and CR. At higher frequency ($f > 4$ Hz), absorption becomes globally dominant ($B_0^* < 0.5$) in the whole Taiwanese crust (see also Tables 2 and 3).

<figure><img src="image_11.png" /><figcaption>**Figure 9.** Albedo maps ($B_0^* = Q_{sc}^{-1}/(Q_i^{-1} + Q_{sc}^{-1})$) in three frequency bands: results of multiple lapse time window analysis are plotted at the location of the station. Red (resp. blue) colors indicate that scattering (resp. absorption) dominates over absorption. Pink lines are surface traces of active faults from Shyu et al. (2016) and the five major lithotectonic units are delineated with thick dark lines (see Figure 1 for more details).</figcaption></figure>
<figure><img src="image_12.png" /><figcaption>**Figure 10.** Scattering anisotropy (g) maps in three frequency bands: results of multiple lapse time window analysis are plotted at the location of the station. Models with strong backscattering (resp. strong forward scattering) are shown in red (resp. blue). Isotropic models are in white. The size of the circles correspond to the misfit ratio between isotropic and anisotropic models ($\chi_{iso}/\chi_{aniso}$). Pink lines are surface traces of active faults from Shyu et al. (2016) and the five major lithotectonic units are delineated with thick dark lines (see Figure 1 for more details).</figcaption></figure>

## 5.4. Scattering Anisotropy

Figure 10 shows the spatial variations of scattering anisotropy (g) across Taiwan in the three frequency bands ([1, 2], [2, 4], and [4, 8] Hz). When $-0.5 < g < 0.5$, scattering is mostly isotropic and the corresponding values of g are shown in white. Red (resp. blue) colors correspond to predominantly backward (resp. forward)-scattering. The size of the colored disks depends on the ratio between the misfit function ($\chi^2$) for an isotropic and an anisotropic model. $\chi_{iso}^2/\chi_{aniso}^2 > 1.5$ indicates that an anisotropic model significantly improves the fit to the data. As previously discussed, anisotropic models with strong back-scattering ($g < -0.8$) better explain the data than isotropic ones in certain areas. At low frequency (1–2 Hz), g is globally negative in the crust of Taiwan (Table 2). The zones of strong back-scattering ($g < -0.8$) are mainly located close to active faults in CP, WF and IP. Stations in LV and CoR are also globally characterized by negative values of g. By contrast, crystalline massifs in CR and HR exhibit higher values of scattering anisotropy ($g > -0.5$). On average g slightly increases with frequency (see also Table 2), but the spatial distribution of g is almost similar in the two first frequency bands. At higher frequency ($f > 4$ Hz), the spatial distribution of g is not spatially coherent and does not show clear relations with tectonic structures or geological units.

To summarize, we observe that volcanic materials in LV are characterized by strong scattering and absorption. This property of volcanos is broadly compatible with the literature (Carcolé & Sato, 2010; Del Pezzo et al., 2001;

<table><caption>**Table 2**</caption><caption>Mean, Minimum and Maximum of Absorption ($Q_i^{-1}$), Effective Scattering ($Q_{sc}^{-1}$) and Scattering Anisotropy (g) in Taiwan in Three Frequency Bands</caption><thead><tr><th>Frequency</th><th>$Q_i^{-1}$</th><th>Min</th><th>Max</th><th>$Q_{sc}^{-1}$</th><th>Min</th><th>Max</th><th>$\bar{g}$</th><th>Min</th><th>Max</th></tr></thead><tbody><tr><td>1–2 Hz</td><td>0.00937</td><td>0.00537</td><td>0.01262</td><td>0.02531</td><td>0.00505</td><td>0.08547</td><td>-0.49</td><td>-1</td><td>0.86</td></tr><tr><td>2–4 Hz</td><td>0.00569</td><td>0.00388</td><td>0.00775</td><td>0.00620</td><td>0.00045</td><td>0.02378</td><td>-0.29</td><td>-1</td><td>0.88</td></tr><tr><td>4–8 Hz</td><td>0.003218</td><td>0.00219</td><td>0.00413</td><td>0.00173</td><td>0.00035</td><td>0.00728</td><td>-0.03</td><td>-1</td><td>0.86</td></tr></tbody></table>
<caption>Table 3</caption>

Mean Absorption ($Q_i^{-1}$), Effective Scattering ($Q_{sc}^{*^{-1}}$) and Total Attenuation ($Q_i^{-1}$) in the Three Similar Tectonic Units in Taiwan

<table><thead><tr><th rowspan="2">Tectonic units</th><th colspan="3">Q<sub>sc</sub><sup>*</sup> <sup>-1</sup> × 10<sup>-3</sup></th><th colspan="3">Q<sub>i</sub><sup>-1</sup> × 10<sup>-3</sup></th><th colspan="3">Q<sub>i</sub><sup>-1</sup> × 10<sup>-3</sup></th></tr><tr><th>1–2 Hz</th><th>2–4 Hz</th><th>4–8 Hz</th><th>1–2 Hz</th><th>2–4 Hz</th><th>4–8 Hz</th><th>1–2 Hz</th><th>2–4 Hz</th><th>4–8 Hz</th></tr></thead><tbody><tr><td>HR + CR</td><td>11.688</td><td>3.660</td><td>1.663</td><td>8.744</td><td>5.445</td><td>3.521</td><td>20.431</td><td>9.105</td><td>5.184</td></tr><tr><td>CoR + LV</td><td>48.199</td><td>11.582</td><td>2.446</td><td>10.689</td><td>7.200</td><td>3.647</td><td>58.888</td><td>18.782</td><td>6.093</td></tr><tr><td>WF + CP + IP</td><td>42.453</td><td>10.390</td><td>2.124</td><td>10.019</td><td>5.852</td><td>2.586</td><td>52.472</td><td>16.242</td><td>4.710</td></tr></tbody></table>

Nishigami, 1997; Takahashi et al., 2009; Wegler, 2003). Crystalline massifs show rather low scattering and absorption, again in agreement with results from previous studies (Carcolé & Sato, 2010; Lacombe et al., 2003; Mayor et al., 2016; Sens-Schönfelder et al., 2009). Active faults are marked by strong scattering as previously observed in Japan or in the San Andreas Fault (Asano et al., 2004; Nishigami, 1997, 2000). The back-scattering signature of fault zones and volcanic areas appears to be unique to Taiwan.

# 6. Discussion

## 6.1. Comparison With Others Studies

Y.-J. Wang et al. (2010) have previously determined the attenuation structure of Taiwan, from $t^*$ measurements on P- and S-wave spectra. In their study, the authors hypothesize that seismic attenuation is (a) frequency independent, (b) mainly controlled by anelastic processes (effects of temperature and/or fluid saturation). In the first 10 km, they observe that central and northern CR are globally less attenuating than southern CR, CoR/LV, and CP. Northern WF shows lower attenuation than HR while the southern WF (mainly the Pingtung Plain) and Southern CR are more attenuating than the surrounding regions. Y.-J. Wang et al. (2010) proposed that the strong attenuation anomaly in southern CR may be ascribed to the lower crust exhumation associated with an orogen-parallel extrusion. They also systematically observed a sharp contrast of the total attenuation across the Chelungpu, Kaoping and Chaochou faults. Our results are rather different. First, Figures 7–10 and Tables 2 and 3 clearly show the frequency dependence of attenuation: (a) the mean total attenuation varies over about one order of magnitude from $3.467 \times 10^{-2}$ at 1.5 Hz to $4.95 \times 10^{-3}$ at 6 Hz, (b) the frequency dependence is stronger for scattering than for absorption. Second, scattering is generally dominant in CP, WF, LV, and CoR (see Table 3), so that $t^*$-measurements cannot be directly interpreted as absorption only. Third, even if our spatial coverage is not perfect, most active faults in the WF (from north to south) are characterized by both strong scattering and absorption at low frequency. In particular, the Chelungpu Fault, as well as the southern part of the Kaoping fault, are mainly marked by strong absorption and scattering (at low frequency) both in the hangingwall and the footwall of the fault. By contrast, the Chaochou Fault is rather characterized by low scattering and low absorption. Lastly, northern WF is systematically more attenuating than HR, and CR has overall both low scattering and low absorption compared to surrounding regions. Consequently, the absence of an attenuation anomaly in southern CR may not support the collisional tectonic model with a restricted subduction of the Eurasian crust discussed by Y.-J. Wang et al. (2010).

The seismic attenuation structure of Taiwan is clearly original compared to other regions of the world. In Figure 11, we have reported $1/Q_{sc}^*$ and $1/Q_i$ in the frequency band 1–2 Hz for two kinds of geological structures which correspond to both extremes in terms of seismic attenuation: volcanoes and old orogens. Volcanoes are characterized by strong scattering ($Q_{sc}^{*-1} \sim [0.01 - 1]$) and strong absorption ($Q_i^{-1} \sim [0.01 - 0.2]$) while old orogens are the least attenuating tectonic structures on Earth ($1/Q_{sc}$ and $1/Q_i$ smaller than 0.001). To make the comparison meaningful, we have distinguished between the main geological units of Taiwan as follows: crystalline massifs (HR and CR) are shown in blue, sedimentary materials (CP, IP, and WF) are in green, and areas affected by volcanism (mainly LV/CoR) are in yellow. Globally, the values of scattering and absorption in Taiwan are very close to the one of volcanoes. Sediments in CP and IP, highly tectonized areas in WF, and volcanic materials in LV have a similar level of scattering as volcanoes while absorption is close to the higher bound of $Q_i^{-1}$ for volcanoes. By contrast, Taiwanese crystalline massifs (CR and HR) are more attenuating by one order of magnitude (both in scattering and absorption) than old orogens. In Figures 11 and 12, we also compare Taiwanese attenuation data to other regions (Japan, North America and Korea).
<caption>Table 4</caption>

Results of Multiple Lapse Time Window Analysis: Total Scattering Attenuation (1/Q\*\_sc), Absorption (1/Q\_i), Scattering Anisotropy (g) and Albedo (B\*\_0) Are Given for the Three Frequency Bands

<table><thead><tr><th rowspan="2">Station</th><th colspan="3">1/Q<sub>i</sub> × 10<sup>-3</sup></th><th colspan="3">1/Q<sub>sc</sub> × 10<sup>-3</sup></th><th colspan="3">Scattering anisotropy (g)</th><th colspan="3">B<sub>0</sub><sup>*</sup></th></tr><tr><th>[1-2]</th><th>[2-4]</th><th>[4-8]</th><th>[1-2]</th><th>[2-4]</th><th>[4-8]</th><th>[1-2]</th><th>[2-4]</th><th>[4-8]</th><th>[1-2]</th><th>[2-4]</th><th>[4-8]</th></tr></thead><tbody><tr><td>ALS</td><td>8.62</td><td>4.45</td><td>2.47</td><td>9.83</td><td>3.05</td><td>1.20</td><td>-0.13</td><td>0.77</td><td>0.75</td><td>0.53</td><td>0.41</td><td>0.33</td></tr><tr><td>CHGB</td><td>7.83</td><td>5.09</td><td>3.80</td><td>8.23</td><td>3.47</td><td>1.43</td><td>-0.15</td><td>0.45</td><td>-0.92</td><td>0.51</td><td>0.41</td><td>0.27</td></tr><tr><td>CHK</td><td>11.02</td><td>7.75</td><td>3.76</td><td>33.14</td><td>14.44</td><td>2.76</td><td>-1.00</td><td>-1.00</td><td>0.75</td><td>0.75</td><td>0.65</td><td>0.42</td></tr><tr><td>CHN1</td><td>10.38</td><td>4.30</td><td>2.55</td><td>17.83</td><td>2.38</td><td>0.87</td><td>-0.91</td><td>-0.04</td><td>0.28</td><td>0.63</td><td>0.36</td><td>0.26</td></tr><tr><td>CHN2</td><td>10.64</td><td>6.77</td><td>-</td><td>53.67</td><td>23.78</td><td>-</td><td>-0.99</td><td>-0.99</td><td>-</td><td>0.83</td><td>0.78</td><td>-</td></tr><tr><td>CHN3</td><td>9.38</td><td>6.64</td><td>-</td><td>30.54</td><td>12.91</td><td>-</td><td>-0.99</td><td>-0.95</td><td>-</td><td>0.77</td><td>0.66</td><td>-</td></tr><tr><td>CHN4</td><td>11.59</td><td>4.49</td><td>2.45</td><td>25.56</td><td>2.82</td><td>1.17</td><td>-0.95</td><td>0.28</td><td>0.85</td><td>0.69</td><td>0.39</td><td>0.32</td></tr><tr><td>CHN5</td><td>10.13</td><td>5.81</td><td>2.38</td><td>20.23</td><td>7.26</td><td>0.85</td><td>-0.97</td><td>-0.95</td><td>0.26</td><td>0.67</td><td>0.56</td><td>0.26</td></tr><tr><td>CHN8</td><td>5.99</td><td>5.17</td><td>3.60</td><td>54.74</td><td>22.71</td><td>7.28</td><td>-0.99</td><td>-0.99</td><td>-0.99</td><td>0.90</td><td>0.81</td><td>0.67</td></tr><tr><td>CHY</td><td>10.38</td><td>6.72</td><td>-</td><td>56.49</td><td>22.82</td><td>-</td><td>-1.00</td><td>-0.99</td><td>-</td><td>0.84</td><td>0.77</td><td>-</td></tr><tr><td>DPDB</td><td>12.58</td><td>6.53</td><td>-</td><td>6.58</td><td>1.29</td><td>-</td><td>-0.98</td><td>-0.96</td><td>-</td><td>0.34</td><td>0.16</td><td>-</td></tr><tr><td>EAS</td><td>8.07</td><td>6.07</td><td>3.78</td><td>11.94</td><td>4.69</td><td>1.78</td><td>0.85</td><td>0.88</td><td>0.85</td><td>0.60</td><td>0.44</td><td>0.32</td></tr><tr><td>ECL</td><td>8.19</td><td>5.71</td><td>3.75</td><td>12.43</td><td>4.31</td><td>1.83</td><td>0.35</td><td>-0.20</td><td>0.45</td><td>0.60</td><td>0.43</td><td>0.33</td></tr><tr><td>EGFH</td><td>10.77</td><td>6.58</td><td>-</td><td>65.72</td><td>9.52</td><td>-</td><td>-0.96</td><td>-0.36</td><td>-</td><td>0.86</td><td>0.59</td><td>-</td></tr><tr><td>EGS</td><td>7.54</td><td>5.81</td><td>-</td><td>38.95</td><td>8.16</td><td>-</td><td>-0.88</td><td>-0.65</td><td>-</td><td>0.84</td><td>0.58</td><td>-</td></tr><tr><td>EHY</td><td>12.62</td><td>7.38</td><td>3.39</td><td>43.11</td><td>13.67</td><td>1.77</td><td>-0.96</td><td>-0.96</td><td>-0.13</td><td>0.77</td><td>0.65</td><td>0.34</td></tr><tr><td>ELD</td><td>8.78</td><td>5.26</td><td>3.18</td><td>8.87</td><td>3.12</td><td>1.69</td><td>-0.04</td><td>-0.15</td><td>-0.15</td><td>0.50</td><td>0.37</td><td>0.35</td></tr><tr><td>ENA</td><td>10.37</td><td>6.73</td><td>3.91</td><td>17.91</td><td>8.21</td><td>4.15</td><td>-0.98</td><td>-1.00</td><td>-1.00</td><td>0.63</td><td>0.55</td><td>0.52</td></tr><tr><td>ENT</td><td>7.85</td><td>5.14</td><td>3.44</td><td>6.83</td><td>2.49</td><td>0.94</td><td>0.05</td><td>0.37</td><td>-0.05</td><td>0.47</td><td>0.33</td><td>0.22</td></tr><tr><td>ESL</td><td>9.29</td><td>5.34</td><td>3.42</td><td>13.07</td><td>3.46</td><td>2.25</td><td>-0.82</td><td>-0.15</td><td>0.55</td><td>0.58</td><td>0.39</td><td>0.40</td></tr><tr><td>ETLH</td><td>9.17</td><td>5.19</td><td>4.10</td><td>9.40</td><td>2.93</td><td>2.64</td><td>-0.76</td><td>-0.25</td><td>-0.94</td><td>0.51</td><td>0.36</td><td>0.39</td></tr><tr><td>FULB</td><td>10.22</td><td>6.90</td><td>-</td><td>17.12</td><td>6.02</td><td>-</td><td>0.86</td><td>0.85</td><td>-</td><td>0.63</td><td>0.47</td><td>-</td></tr><tr><td>FUSS</td><td>8.94</td><td>5.33</td><td>3.31</td><td>15.41</td><td>3.86</td><td>1.56</td><td>0.84</td><td>0.65</td><td>0.85</td><td>0.63</td><td>0.42</td><td>0.32</td></tr><tr><td>HEN</td><td>8.26</td><td>6.79</td><td>-</td><td>22.83</td><td>12.74</td><td>-</td><td>-0.85</td><td>-0.99</td><td>-</td><td>0.73</td><td>0.65</td><td>-</td></tr><tr><td>HSN</td><td>10.95</td><td>-</td><td>-</td><td>85.47</td><td>-</td><td>-</td><td>-0.97</td><td>-</td><td>-</td><td>0.89</td><td>-</td><td>-</td></tr><tr><td>HWA</td><td>9.18</td><td>7.23</td><td>-</td><td>51.08</td><td>9.64</td><td>-</td><td>-0.99</td><td>-1.00</td><td>-</td><td>0.85</td><td>0.57</td><td>-</td></tr><tr><td>ILA</td><td>7.61</td><td>5.55</td><td>-</td><td>30.13</td><td>7.43</td><td>-</td><td>-0.99</td><td>-0.97</td><td>-</td><td>0.80</td><td>0.57</td><td>-</td></tr><tr><td>LAY</td><td>10.22</td><td>-</td><td>-</td><td>32.02</td><td>-</td><td>-</td><td>-0.99</td><td>-</td><td>-</td><td>0.76</td><td>-</td><td>-</td></tr><tr><td>LIOB</td><td>9.28</td><td>-</td><td>-</td><td>20.62</td><td>-</td><td>-</td><td>-0.99</td><td>-</td><td>-</td><td>0.69</td><td>-</td><td>-</td></tr><tr><td>MASB</td><td>5.61</td><td>5.13</td><td>-</td><td>5.05</td><td>2.55</td><td>-</td><td>-0.18</td><td>-0.25</td><td>-</td><td>0.47</td><td>0.33</td><td>-</td></tr><tr><td>NACB</td><td>9.81</td><td>6.57</td><td>3.72</td><td>17.46</td><td>5.98</td><td>1.78</td><td>-0.95</td><td>-0.99</td><td>-0.87</td><td>0.64</td><td>0.48</td><td>0.32</td></tr><tr><td>NANB</td><td>7.88</td><td>5.01</td><td>-</td><td>16.07</td><td>3.59</td><td>-</td><td>0.85</td><td>0.20</td><td>-</td><td>0.67</td><td>0.42</td><td>-</td></tr><tr><td>NDS</td><td>10.63</td><td>-</td><td>-</td><td>12.20</td><td>-</td><td>-</td><td>-0.85</td><td>-</td><td>-</td><td>0.53</td><td>-</td><td>-</td></tr><tr><td>NDT</td><td>7.75</td><td>4.92</td><td>4.13</td><td>5.82</td><td>1.43</td><td>0.66</td><td>-0.15</td><td>-0.10</td><td>-0.98</td><td>0.43</td><td>0.23</td><td>0.14</td></tr><tr><td>NHDH</td><td>11.67</td><td>-</td><td>-</td><td>45.21</td><td>-</td><td>-</td><td>-0.87</td><td>-</td><td>-</td><td>0.79</td><td>-</td><td>-</td></tr><tr><td>NMLH</td><td>10.00</td><td>-</td><td>-</td><td>19.50</td><td>-</td><td>-</td><td>-0.95</td><td>-</td><td>-</td><td>0.66</td><td>-</td><td>-</td></tr><tr><td>NNS</td><td>8.70</td><td>5.17</td><td>3.28</td><td>10.51</td><td>3.35</td><td>1.64</td><td>0.25</td><td>0.65</td><td>0.85</td><td>0.55</td><td>0.39</td><td>0.33</td></tr><tr><td>NNSB</td><td>8.49</td><td>5.22</td><td>3.33</td><td>7.92</td><td>2.83</td><td>0.82</td><td>-0.55</td><td>0.46</td><td>-0.16</td><td>0.48</td><td>0.35</td><td>0.20</td></tr><tr><td>NNSH</td><td>9.75</td><td>6.17</td><td>-</td><td>13.76</td><td>2.80</td><td>-</td><td>-0.95</td><td>-0.99</td><td>-</td><td>0.59</td><td>0.31</td><td>-</td></tr><tr><td>NSK</td><td>8.00</td><td>4.51</td><td>3.73</td><td>6.59</td><td>2.00</td><td>1.24</td><td>-0.25</td><td>0.05</td><td>-0.92</td><td>0.45</td><td>0.31</td><td>0.25</td></tr><tr><td>NST</td><td>9.93</td><td>5.27</td><td>2.19</td><td>14.49</td><td>1.86</td><td>0.35</td><td>-0.91</td><td>-0.76</td><td>0.15</td><td>0.59</td><td>0.26</td><td>0.14</td></tr></tbody></table>
<caption>**Table 4**</caption>

Continued

<table><thead><tr><th rowspan="2">Station</th><th colspan="3">1/Q<sub>i</sub> × 10<sup>-3</sup></th><th colspan="3">1/Q<sub>sc</sub> × 10<sup>-3</sup></th><th colspan="3">Scattering anisotropy (g)</th><th colspan="3">B<sub>0</sub><sup>*</sup></th></tr><tr><th>[1-2]</th><th>[2-4]</th><th>[4-8]</th><th>[1-2]</th><th>[2-4]</th><th>[4-8]</th><th>[1-2]</th><th>[2-4]</th><th>[4-8]</th><th>[1-2]</th><th>[2-4]</th><th>[4-8]</th></tr></thead><tbody><tr><td>NSY</td><td>10.74</td><td>6.07</td><td>—</td><td>33.81</td><td>7.53</td><td>—</td><td>-0.98</td><td>-0.98</td><td>—</td><td>0.76</td><td>0.55</td><td>—</td></tr><tr><td>NWF</td><td>8.79</td><td>5.38</td><td>—</td><td>12.86</td><td>4.97</td><td>—</td><td>-0.25</td><td>0.85</td><td>—</td><td>0.59</td><td>0.48</td><td>—</td></tr><tr><td>NWL</td><td>8.25</td><td>—</td><td>—</td><td>11.95</td><td>—</td><td>—</td><td>-0.05</td><td>—</td><td>—</td><td>0.59</td><td>—</td><td>—</td></tr><tr><td>SCL</td><td>6.34</td><td>4.78</td><td>—</td><td>42.00</td><td>16.53</td><td>—</td><td>-0.99</td><td>-0.94</td><td>—</td><td>0.87</td><td>0.78</td><td>—</td></tr><tr><td>SCZ</td><td>5.37</td><td>4.49</td><td>3.38</td><td>6.25</td><td>1.71</td><td>1.13</td><td>0.85</td><td>0.05</td><td>0.85</td><td>0.54</td><td>0.28</td><td>0.25</td></tr><tr><td>SEB</td><td>9.83</td><td>—</td><td>—</td><td>30.94</td><td>—</td><td>—</td><td>-1.00</td><td>—</td><td>—</td><td>0.76</td><td>—</td><td>—</td></tr><tr><td>SGL</td><td>9.43</td><td>—</td><td>—</td><td>47.65</td><td>—</td><td>—</td><td>-0.99</td><td>—</td><td>—</td><td>0.83</td><td>—</td><td>—</td></tr><tr><td>SGS</td><td>7.53</td><td>4.36</td><td>2.65</td><td>6.89</td><td>2.26</td><td>1.00</td><td>-0.14</td><td>0.07</td><td>0.45</td><td>0.48</td><td>0.34</td><td>0.27</td></tr><tr><td>SLG</td><td>7.68</td><td>5.05</td><td>—</td><td>12.35</td><td>3.34</td><td>—</td><td>0.25</td><td>0.65</td><td>—</td><td>0.62</td><td>0.40</td><td>—</td></tr><tr><td>SML</td><td>8.34</td><td>5.11</td><td>2.31</td><td>8.18</td><td>2.29</td><td>0.83</td><td>-0.35</td><td>-0.79</td><td>0.05</td><td>0.50</td><td>0.31</td><td>0.26</td></tr><tr><td>SNS</td><td>10.62</td><td>4.91</td><td>—</td><td>13.59</td><td>2.31</td><td>—</td><td>-0.89</td><td>0.02</td><td>—</td><td>0.56</td><td>0.32</td><td>—</td></tr><tr><td>SSD</td><td>6.59</td><td>6.80</td><td>4.07</td><td>9.27</td><td>7.54</td><td>2.55</td><td>0.26</td><td>-0.95</td><td>-0.96</td><td>0.58</td><td>0.53</td><td>0.39</td></tr><tr><td>SSLB</td><td>8.11</td><td>4.66</td><td>3.58</td><td>7.93</td><td>2.47</td><td>1.41</td><td>-0.25</td><td>0.05</td><td>-0.89</td><td>0.49</td><td>0.35</td><td>0.28</td></tr><tr><td>STY</td><td>8.52</td><td>4.83</td><td>2.61</td><td>8.96</td><td>2.82</td><td>0.62</td><td>0.01</td><td>0.65</td><td>-0.05</td><td>0.51</td><td>0.37</td><td>0.19</td></tr><tr><td>TAI</td><td>8.24</td><td>—</td><td>—</td><td>42.07</td><td>—</td><td>—</td><td>-0.99</td><td>—</td><td>—</td><td>0.84</td><td>—</td><td>—</td></tr><tr><td>TAII</td><td>8.25</td><td>—</td><td>—</td><td>45.21</td><td>—</td><td>—</td><td>-1.00</td><td>—</td><td>—</td><td>0.85</td><td>—</td><td>—</td></tr><tr><td>TAP</td><td>10.72</td><td>—</td><td>—</td><td>45.52</td><td>—</td><td>—</td><td>-0.95</td><td>—</td><td>—</td><td>0.81</td><td>—</td><td>—</td></tr><tr><td>TAP1</td><td>10.37</td><td>—</td><td>—</td><td>42.89</td><td>—</td><td>—</td><td>-0.93</td><td>—</td><td>—</td><td>0.81</td><td>—</td><td>—</td></tr><tr><td>TAW</td><td>6.71</td><td>5.60</td><td>—</td><td>12.84</td><td>5.20</td><td>—</td><td>0.85</td><td>0.83</td><td>—</td><td>0.66</td><td>0.48</td><td>—</td></tr><tr><td>TCU</td><td>10.16</td><td>5.69</td><td>2.41</td><td>47.29</td><td>5.28</td><td>1.52</td><td>-0.98</td><td>-0.96</td><td>0.85</td><td>0.82</td><td>0.48</td><td>0.39</td></tr><tr><td>TDCB</td><td>9.77</td><td>5.29</td><td>—</td><td>17.85</td><td>3.89</td><td>—</td><td>0.55</td><td>0.46</td><td>—</td><td>0.65</td><td>0.42</td><td>—</td></tr><tr><td>TPUB</td><td>7.94</td><td>4.78</td><td>2.85</td><td>9.26</td><td>3.35</td><td>1.10</td><td>-0.06</td><td>0.85</td><td>0.85</td><td>0.54</td><td>0.41</td><td>0.28</td></tr><tr><td>TTN</td><td>10.50</td><td>7.27</td><td>—</td><td>65.18</td><td>13.07</td><td>—</td><td>-1.00</td><td>-0.96</td><td>—</td><td>0.86</td><td>0.64</td><td>—</td></tr><tr><td>TWA</td><td>8.40</td><td>—</td><td>—</td><td>16.20</td><td>—</td><td>—</td><td>0.45</td><td>—</td><td>—</td><td>0.66</td><td>—</td><td>—</td></tr><tr><td>TWB1</td><td>8.22</td><td>4.94</td><td>—</td><td>13.87</td><td>3.23</td><td>—</td><td>-0.96</td><td>-0.48</td><td>—</td><td>0.63</td><td>0.40</td><td>—</td></tr><tr><td>TWC</td><td>9.27</td><td>6.26</td><td>—</td><td>19.80</td><td>8.95</td><td>—</td><td>-1.00</td><td>-1.00</td><td>—</td><td>0.68</td><td>0.59</td><td>—</td></tr><tr><td>TWD</td><td>8.55</td><td>6.49</td><td>3.25</td><td>20.22</td><td>9.39</td><td>1.62</td><td>-0.99</td><td>-1.00</td><td>-0.55</td><td>0.70</td><td>0.59</td><td>0.33</td></tr><tr><td>TWE</td><td>9.80</td><td>6.35</td><td>3.34</td><td>17.01</td><td>5.94</td><td>2.28</td><td>-0.99</td><td>-0.97</td><td>0.85</td><td>0.63</td><td>0.48</td><td>0.41</td></tr><tr><td>TWF1</td><td>12.45</td><td>7.47</td><td>3.53</td><td>56.96</td><td>16.79</td><td>2.13</td><td>-0.99</td><td>-0.95</td><td>-0.15</td><td>0.82</td><td>0.69</td><td>0.38</td></tr><tr><td>TWG</td><td>10.72</td><td>6.71</td><td>—</td><td>21.63</td><td>7.04</td><td>—</td><td>-0.96</td><td>-0.90</td><td>—</td><td>0.67</td><td>0.51</td><td>—</td></tr><tr><td>TWGB</td><td>8.19</td><td>—</td><td>—</td><td>12.41</td><td>—</td><td>—</td><td>0.85</td><td>—</td><td>—</td><td>0.60</td><td>—</td><td>—</td></tr><tr><td>TWK1</td><td>6.85</td><td>5.59</td><td>—</td><td>16.59</td><td>7.55</td><td>—</td><td>0.85</td><td>0.85</td><td>—</td><td>0.71</td><td>0.57</td><td>—</td></tr><tr><td>TWL</td><td>11.21</td><td>5.87</td><td>2.43</td><td>29.15</td><td>5.24</td><td>0.92</td><td>-0.98</td><td>-0.93</td><td>0.25</td><td>0.72</td><td>0.47</td><td>0.28</td></tr><tr><td>TWM1</td><td>12.36</td><td>—</td><td>—</td><td>58.19</td><td>—</td><td>—</td><td>-0.98</td><td>—</td><td>—</td><td>0.82</td><td>—</td><td>—</td></tr><tr><td>TWQ1</td><td>11.64</td><td>4.58</td><td>2.19</td><td>42.86</td><td>3.23</td><td>0.97</td><td>-0.95</td><td>-0.03</td><td>0.86</td><td>0.79</td><td>0.41</td><td>0.31</td></tr><tr><td>TWS1</td><td>9.81</td><td>7.03</td><td>—</td><td>35.67</td><td>0.45</td><td>—</td><td>-0.95</td><td>-0.96</td><td>—</td><td>0.78</td><td>0.06</td><td>—</td></tr><tr><td>TWT</td><td>8.92</td><td>6.44</td><td>3.89</td><td>9.06</td><td>5.58</td><td>2.62</td><td>0.00</td><td>-0.99</td><td>-1.00</td><td>0.50</td><td>0.46</td><td>0.40</td></tr><tr><td>TYC</td><td>11.43</td><td>4.45</td><td>2.50</td><td>24.78</td><td>2.51</td><td>0.88</td><td>-0.95</td><td>-0.05</td><td>-0.09</td><td>0.68</td><td>0.36</td><td>0.26</td></tr><tr><td>VWDT</td><td>8.99</td><td>6.04</td><td>3.73</td><td>10.71</td><td>3.06</td><td>0.95</td><td>0.22</td><td>-0.76</td><td>-0.91</td><td>0.54</td><td>0.34</td><td>0.20</td></tr><tr><td>WDLH</td><td>11.49</td><td>—</td><td>—</td><td>74.17</td><td>—</td><td>—</td><td>-0.98</td><td>—</td><td>—</td><td>0.87</td><td>—</td><td>—</td></tr><tr><td>WFSB</td><td>9.69</td><td>—</td><td>—</td><td>14.97</td><td>—</td><td>—</td><td>0.85</td><td>—</td><td>—</td><td>0.61</td><td>—</td><td>—</td></tr><tr><td>WGK</td><td>11.91</td><td>6.52</td><td>2.61</td><td>56.98</td><td>16.00</td><td>2.10</td><td>-0.99</td><td>-0.99</td><td>0.85</td><td>0.83</td><td>0.71</td><td>0.45</td></tr></tbody></table>
<table><caption>**Table 4**
*Continued*</caption><thead><tr><th rowspan="2">Station</th><th colspan="3">1/Q<sub>i</sub> × 10<sup>-3</sup></th><th colspan="3">1/Q<sub>sc</sub> × 10<sup>-3</sup></th><th colspan="3">Scattering anisotropy (g)</th><th colspan="3">B<sub>0</sub><sup>*</sup></th></tr><tr><th>[1-2]</th><th>[2-4]</th><th>[4-8]</th><th>[1-2]</th><th>[2-4]</th><th>[4-8]</th><th>[1-2]</th><th>[2-4]</th><th>[4-8]</th><th>[1-2]</th><th>[2-4]</th><th>[4-8]</th></tr></thead><tbody><tr><td>WHF</td><td>8.65</td><td>5.69</td><td>4.04</td><td>11.03</td><td>3.79</td><td>3.80</td><td>0.25</td><td>0.05</td><td>-0.99</td><td>0.56</td><td>0.40</td><td>0.48</td></tr><tr><td>WHP</td><td>7.55</td><td>3.88</td><td>-</td><td>8.74</td><td>1.58</td><td>-</td><td>0.85</td><td>0.75</td><td>-</td><td>0.54</td><td>0.29</td><td>-</td></tr><tr><td>WHY</td><td>8.28</td><td>-</td><td>-</td><td>13.20</td><td>-</td><td>-</td><td>0.34</td><td>-</td><td>-</td><td>0.61</td><td>-</td><td>-</td></tr><tr><td>WJS</td><td>11.92</td><td>6.35</td><td>-</td><td>50.53</td><td>2.48</td><td>-</td><td>-0.98</td><td>-0.97</td><td>-</td><td>0.81</td><td>0.28</td><td>-</td></tr><tr><td>WNT</td><td>10.82</td><td>4.19</td><td>2.49</td><td>18.26</td><td>2.06</td><td>1.22</td><td>-0.98</td><td>0.14</td><td>0.15</td><td>0.63</td><td>0.33</td><td>0.33</td></tr><tr><td>WTP</td><td>11.77</td><td>4.97</td><td>2.76</td><td>26.30</td><td>4.10</td><td>1.22</td><td>-0.95</td><td>0.65</td><td>0.55</td><td>0.69</td><td>0.45</td><td>0.31</td></tr><tr><td>YHNB</td><td>8.44</td><td>6.28</td><td>3.65</td><td>10.34</td><td>4.45</td><td>1.31</td><td>0.05</td><td>-1.00</td><td>-0.93</td><td>0.55</td><td>0.41</td><td>0.26</td></tr><tr><td>YULB</td><td>12.53</td><td>6.19</td><td>3.41</td><td>52.00</td><td>6.86</td><td>1.88</td><td>-0.96</td><td>-0.41</td><td>0.05</td><td>0.81</td><td>0.53</td><td>0.35</td></tr><tr><td>YUS</td><td>9.15</td><td>5.16</td><td>3.48</td><td>11.60</td><td>4.42</td><td>3.13</td><td>0.45</td><td>0.87</td><td>-1.00</td><td>0.56</td><td>0.46</td><td>0.47</td></tr></tbody></table>

As Taiwan and Japan are both in a subduction environment, one might expect similar levels of attenuation. Figures 11 and 12 shows that: (a) in the frequency band [1–8] Hz, absorption is systematically stronger in Taiwan than in Japan (at least by a factor 2 at low frequency $f < 2$ Hz), North America (in the whole frequency range) and Korea (except in the first frequency band); (b) scattering attenuation in Taiwan is one order of magnitude higher than in Japan, North America and Korea. Carcolé and Sato (2010) have also applied a

MLTWA approach to derive absorption and scattering maps, but they have used an isotropic scattering model (as well as Eulenfeld and Wegler (2017) for North America and Rachman et al. (2015) for Korea). As discussed in the previous section, an isotropic model may systematically underestimate attenuation (both absorption and scattering) compared to an anisotropic model with strong back-scattering. Consequently, distinct hypothesis on scattering anisotropy may partially explain some of the differences between Taiwan and Japan (or North America and Korea). Of course, this argument holds only in the case where the spatio-temporal distribution of energy is anomalous, an hypothesis for which we have no direct confirmation so far in the case of Japan. It is also noteworthy in Figure 11 that the spatial variations of attenuation in Taiwan are much larger than in Japan which covers a larger area -compared to Taiwan- but with a comparable diversity of geological structures.

## 6.2. On the Distribution of Fluids in the Taiwanese Crust

Huang et al. (2014) estimated the velocity structure ($V_P$ and $V_S$) of the taiwanese crust. They obtained a reliable $V_P/V_S$ model constrained at the surface (from 30–60–500 m depth) using both seismic data and logging data from the Engineering Geological Database for TSMIP (EGDT). From the surface to 10 km depth, the main plains and basins (Taipei Bassin, Western Plain and Pingtung Plain), the LV, and the highly tectonically active regions in WF and CP show low $V_P$ and $V_S$ and high $V_P/V_S$ ratios (Figure 13). These regions also correlate with areas of strong absorption. High $V_P/V_S$ ratio, low velocities and strong absorption are usually interpreted in terms of fluid content as many mechanisms of intrinsic attenuation are related to the presence of fluids within cracks and/or pores in crustal rocks (see Mavko et al. (2020) for a review). The Pingtung Plain, characterized by low velocities, high $V_P/V_S$ and strong attenuation, is associated with a rapid deposition environment and accumulates unconsolidated coastal and estuarine sediments over more than 5 km thickness. It is reasonable to assume that such

<figure><img src="image_13.png" /><figcaption>**Figure 11.** $Q_{sc}^{-1}$ and $Q_i^{-1}$ in the frequency band 1–2 Hz for Taiwan (circles) and other regions of the world (diamonds and boxes). Colored circles indicate similar tectonic units in Taiwan: crystalline massifs (HR and CR) are in blue, sedimentary and highly tectonized areas (Western Foothills, Coastal Plain, and Ilan Plain) are in green, and stations in Coastal Range and Longitudinal Valley (volcanic materials) are in yellow. Old orogens (mainly western Europe) are marked by blue diamonds (Lacombe et al., 2003; Sens-Schönfelder & Wegler, 2006; Sens-Schönfelder et al., 2009) and volcanoes by red diamonds (Del Pezzo et al., 2001; Mayeda et al., 1992; Ugalde et al., 2010; Wegler, 2003). Minimum and maximum values for absorption and scattering in Japan (Carcolé & Sato, 2010), Korea (Rachman et al., 2015) and North America (Eulenfeld & Wegler, 2017) are also shown.</figcaption></figure>
<figure><img src="image_14.png" /><figcaption>**Figure 12.** $Q_i^{-1}$ (top), $Q_{sc}^{* -1}$ (middle) and $Q_i^{-1} = Q_{sc}^{* -1} + Q_i^{-1}$ (middle) as a function of the frequency. Colored thin lines show attenuation parameters for seismic stations in similar tectonic units in Taiwan: crystalline massifs (HR and CR) in blue, sedimentary and highly tectonized areas (Western Foothills, Coastal Plain, and Ilan Plain) in green, and stations in Coastal Range and Longitudinal Valley (volcanic materials) in yellow. Average attenuation models for Taiwan (Table 3) are indicated by thick lines. Results for Japan (Carcolé & Sato, 2010), Korea (Rachman et al., 2015) and North America (Eulenfeld & Wegler, 2017) are also shown.</figcaption></figure>
<figure><img src="image_15.png" /><figcaption>**Figure 13.** Comparison between the scattering anisotropy *g* at 1.5 Hz (colored circles) and the velocity ratio *V*<sub>P</sub>/*V*<sub>S</sub> at the surface (background color on the right) (Huang et al., 2014). Models with strong backscattering (resp. strong forward scattering) are shown in red (resp. blue). Isotropic models are in white. The size of the circles correspond to the misfit ratio between isotropic and anisotropic models (χ<sub>iso</sub>/χ<sub>aniso</sub>). Pink lines are surface traces of active faults from Shyu et al. (2016) and the five major lithotectonic units are delineated with thick dark lines (see Figure 1 for more details).</figcaption></figure>

materials contain a significant amount of fluids. Similarly, high V<sub>P</sub>/V<sub>S</sub>, low velocities and strong attenuation in LV may also reflect the presence of fluids or melts related to the volcanic activity. These areas are also characterized by strong scattering (Figure 10) but more importantly by strong back-scattering (g < -0.5) at low frequency additional information on the nature of the crustal heterogeneities. Takemura and Yoshimoto (2014) previously demonstrate a strongly heterogeneous structure in the Kanto area (Japan) interpreted as the result of the dehydration of the subducting oceanic crust of the PSP. They assume that these strong small-scale velocity heterogeneities may be related to a random distribution of fluid which may induce low propagation velocities, high -ratio, strong attenuation and generate wide-angle scattering. Wu and Aki (1985) demonstrate that an elastic perturbation may be decomposed into a sum of “velocity-type” and “impedance-type” perturbations. A “velocity-type” heterogeneity corresponds to an object marked by a contrast of velocity without a contrast of impedance with the environments and vice-versa for an “impedance-type” heterogeneity. In the low frequency regime, “impedance-type” perturbations scatter seismic energy in the backward direction producing a negative *g*. By contrast, “velocity-type” perturbations scatter seismic energy in the forward direction resulting in a positive *g*. Fluid inclusions correspond to strong “impedance-type” perturbations. As a consequence, negative-values of *g* at low frequency, combined with other observations such as a low velocities, high V<sub>P</sub>/V<sub>S</sub> ratio and high absorption, could also mark the presence of fluids in LV, CP, WF and in the main sedimentary basins. It is also worth noting that *g* globally increases with frequency in these areas. This observation is compatible with the hypothesis of the presence of fluids since scattering theory predicts *g* that becomes positive at higher frequency even for “impedance-type” heterogeneities.

7. Conclusions

Using a MLTWA approach and a multiple anisotropic scattering model, we have produced the first crustal atten-
uation maps in Taiwan that distinguishes between scattering and anelastic attenuation. Despite a very limited
geographical extent, the Taiwanese crust shows a strong spatial variability of seismic attenuation. Across the
island, scattering varies over more than one order of magnitude and absorption fluctuations are of the order of
a few tens of percent. The spatial distribution and the amplitude of the two attenuation mechanisms vary very
strongly with frequency, but scattering is overall dominant at low frequencies (f < 2 Hz) while anelasticity domi-
nates at higher frequencies. Our results confirm that attenuation maps previously obtained from t* measurements
cannot be directly interpreted as absorption only. The frequency dependence in scattering may be interpreted in
terms of the scale and nature of the heterogeneities (volcanic/magmatic or clastic materials, melts, etc...), and the
frequency dependence of absorption probably depends on the composition and distribution of melts or fluids. In
this study, we demonstrate that scattering anisotropy - quantified by the mean cosine of the scattering angle g- is
needed to properly describe the energy distribution in the S-wave coda of local earthquakes at stations located in
sedimentary, tectonized and volcanic areas. g also strongly varies with frequency. At low frequency, we observe
that strong back-scattering (g < -0.5) is generally correlated with strong scattering, strong absorption, low veloc-
ities and high Vp/Vs ratio. As back-scattering is observed only for impedance fluctuations, it suggests that the
parameter g could be used as a marker of the presence of fluids. Hence, our maps possibly highlight the spatial
distribution of fluids in the Taiwanese crust in particular in the vicinity of active faults.

Our study also demonstrates that the widely-used multiple isotropic scattering model may in some cases be too
simplistic. In particular, neglecting the anisotropy of scattering leads to an underestimation of the attenuation
(both in absorption and in scattering). Even if scattering anisotropy significantly improves the agreement between
the multiple-scattering model and observations in Taiwan, this anisotropic model is probably still too simplistic
given the complexity of the Taiwanese crust. First, the vertical stratification of the velocity and the depth variation
of the Moho as a function of the station location have not been considered in this work. Margerin et al. (1998); Margerin (2017) have demonstrated that the partial confinement of seismic energy in a crustal waveguide gives rise to energy leakage which mimics the effect of absorption and potentially biases the estimation of anelasticity. Second, our MLTWA approach only gives access to the lateral variations of attenuation without anelasticity in the depth structure. In the spirit of Ogiso (2019), future works should explore in more details the sensitivity of the energy ratio in the coda to the vertical stratification and lateral variations of both scattering loss and absorption.

## Data Availability Statement

Data used in this study are available on the CWB datacenter (https://gdmsn.cwb.gov.tw and Citation: https://doi.org/10.7914/SN/T5) and on the BATS datacenter (https://bats.earth.sinica.edu.tw/ and Citation: https://doi.org/10.7914/SN/TW). Softwares used in this study are available at https://doi.org/10.5281/zenodo.7433822.

## Acknowledgments

This research was supported by the International Joint Project for Scientific Corporation between CNRS, France and the Ministry of Science and Technology (MOST), Taiwan under the Grant CNRS-INSU PRC 2146 (LM and MC) and 107-2911-1-002-546- and 108-2911-1-002-514- (SHH). This work was performed using HCP resources from CALMIP supercomputing center (Project P16046). Valuable comments on the paper were provided by the Associate Editor and two anonymous reviewers.

## References

Asano, Y., Obara, K., Nakajima, J., & Hasegawa, A. (2004). Inhomogeneous crustal structure beneath northern Miyagi prefecture, northeastern Japan, imaged by coda envelope inversion: Implication for fluid distribution. *Geophysical Research Letters*, 31(24), L24615. https://doi.org/10.1029/2004GL021261

Brune, J. N. (1970). Tectonic stress and the spectra of seismic shear waves from earthquakes. *Journal of Geophysical Research*, 75(26), 4997–5009. https://doi.org/10.1029/jb075i026p04997

Calvet, M., & Margerin, L. (2013). Lapse-time dependence of coda Q: Anisotropic multiple-scattering models and application to the pyrenees. *Bulletin of the Seismological Society of America*, 103(3), 1993–2010. https://doi.org/10.1785/0120120239

Carcolé, E., & Sato, H. (2010). Spatial distribution of scattering loss and intrinsic absorption of short-period S waves in the lithosphere of Japan on the basis of the multiple lapse time window analysis of hi-net data. *Geophysical Journal International*, 180(1), 268–290. https://doi.org/10.1111/j.1365-246x.2009.04394.x

Chen, K. J. (1998). S-wave attenuation structure in the Taiwan area and its correlation to seismicity. *Terrestrial, Atmospheric and Oceanic Sciences*, 9(1), 97–118. https://doi.org/10.3319/tao.1998.9.1.97(t)

Chen, K.-J., Yeh, Y. H., & Shyu, C. T. (1996). Q<sub>p</sub> structure in the Taiwan area and its correlation to seismicity. *Terrestrial, Atmospheric and Oceanic Sciences*, 7(4), 409–429. https://doi.org/10.3319/tao.1996.7.4.409(acct)

Chen, W.-S., Chung, S.-L., Chou, H.-Y., Zugeerbai, Z., Shao, W.-Y., & Lee, Y.-H. (2017). A reinterpretation of the metamorphic Yuli belt: Evidence for a middle-late Miocene accretionary prism in eastern Taiwan. *Tectonics*, 36(2), 188–206. https://doi.org/10.1002/2016tc004383

Chen, Y.-L., Hung, S.-H., Jiang, J.-S., & Chiao, L.-Y. (2016). Systematic correlations of the earthquake frequency-magnitude distribution with the deformation and mechanical regimes in the Taiwan orogen. *Geophysical Research Letters*, 43(10), 5017–5025. https://doi.org/10.1002/2016gl069020

Cheng, W.-B. (2013). Three-dimensional seismic attenuation structure beneath the Taiwan region and its tectonic implication. *Journal of Asian Earth Sciences*, 65, 86–99. https://doi.org/10.1016/j.jseaes.2012.09.022

Ching, K.-E., Hsieh, M.-L., Johnson, K. M., Chen, K.-H., Rau, R.-J., & Yang, M. (2011). Modern vertical deformation rates and mountain building in Taiwan from precise leveling and continuous GPS observations, 2000–2008. *Journal of Geophysical Research*, 116(B8), B08406. https://doi.org/10.1029/2011JB008242

Chung, J.-K., Chen, Y.-L., & Shin, T.-C. (2009). Spatial distribution of coda Q estimated from local earthquakes in Taiwan area. *Earth Planets and Space*, 61(9), 1077–1088. https://doi.org/10.1186/bf03352958

Dadson, S. J., Hovius, N., Chen, H., Dade, W. B., Hsieh, M.-L., Willett, S. D., et al. (2003). Links between erosion, runoff variability and seismicity in the Taiwan orogen. *Nature*, 426(6967), 648–651. https://doi.org/10.1038/nature02150

Del Pezzo, E., Bianco, F., & Saccorotti, G. (2001). Separation of intrinsic and scattering q for volcanic tremor: An application to Etna and Masaya volcanoes. *Geophysical Research Letters*, 28(16), 3083–3086. https://doi.org/10.1029/2001gl013372

Eulenfeld, T., & Wegler, U. (2017). Crustal intrinsic and scattering attenuation of high-frequency shear waves in the contiguous United States. *Journal of Geophysical Research: Solid Earth*, 122(6), 4676–4690. https://doi.org/10.1002/2017jb014038

Fehler, M., Hoshiwa, M., Sato, H., & Obara, K. (1992). Separation of scattering and intrinsic attenuation for the Kanto-Tokai region, Japan, using measurements of S-wave energy versus hypocentral distance. *Geophysical Journal International*, 108(3), 787–800. https://doi.org/10.1111/j.1365-246x.1992.tb03470.x

Gaebler, P. J., Eulenfeld, T., & Wegler, U. (2015). Seismic scattering and absorption parameters in the W-Bohemia/Vogtland region from elastic and acoustic radiative transfer theory. *Geophysical Supplements to the Monthly Notices of the Royal Astronomical Society*, 203(3), 1471–1481. https://doi.org/10.1093/gji/ggv393

Galassi, M., Davies, J., Theiler, J., Gough, B., Jungman, G., Alken, P., et al. (2002). GNU scientific library. Citeseer.

Gusev, A. A., & Abubakirov, I. R. (1996). Simulated envelopes of non-isotropically scattered body waves as compared to observed ones: Another manifestation of fractal heterogeneity. *Geophysical Journal International*, 127(1), 49–60. https://doi.org/10.1111/j.1365-246x.1996.tb01534.x

Heller, G., Margerin, L., Sèbe, O., Mayor, J., & Calvet, M. (2022). Revisiting multiple-scattering principles in a crustal waveguide: Equipartition, depolarization and coda normalization. *Pure and Applied Geophysics*, 179(6-7), 1–35. https://doi.org/10.1007/s00024-022-03063-3

Henyey, L. G., & Greenstein, J. L. (1941). Diffuse radiation in the galaxy. *The Astrophysical Journal*, 93, 70–83. https://doi.org/10.1086/144246

Ho, C. (1988). *An introduction to the geology of Taiwan*. Ministry of Economic Affairs Taiwan, Central Geological Survey.

Hoshiwa, M. (1993). Separation of scattering attenuation and intrinsic absorption in Japan using the multiple lapse time window analysis of full seismogram envelope. *Journal of Geophysical Research*, 98(B9), 15809–15824. https://doi.org/10.1029/93jb00347

Hoshiwa, M. (1995). Estimation of nonisotropic scattering in western Japan using coda wave envelopes: Application of a multiple nonisotropic scattering model. *Journal of Geophysical Research*, 100(B1), 645–657. https://doi.org/10.1029/94jb02064

Hsu, T. L. (1976). The Lichi melange in the coastal range framework. *Bulletin of the Geological Survey of Taiwan*, 25, 87–96.
Huang, H.-H., Wu, Y.-M., Song, X., Chang, C.-H., Lee, S.-J., Chang, T.-M., & Hsieh, H.-H. (2014). Joint $V_p$ and $V_s$ tomography of Taiwan: Implications for subduction-collision orogeny. *Earth and Planetary Science Letters*, **392**, 177–191. https://doi.org/10.1016/j.epsl.2014.02.026

Ko, Y.-T., Kuo, B.-Y., & Hung, S.-H. (2012). Robust determination of earthquake source parameters and mantle attenuation. *Journal of Geophysical Research*, **117**(B4), B04304. https://doi.org/10.1029/2011JB008759

Ko, Y.-T., Kuo, B.-Y., Wang, K.-L., Lin, S.-C., & Hung, S.-H. (2012). The southwestern edge of the Ryukyu subduction zone: A high Q mantle wedge. *Earth and Planetary Science Letters*, **335**, 145–153. https://doi.org/10.1016/j.epsl.2012.04.041

Lacombe, C., Campillo, M., Paul, A., & Margerin, L. (2003). Separation of intrinsic absorption and scattering attenuation from $l_g$ coda decay in central France using acoustic radiative transfer theory. *Geophysical Journal International*, **154**(2), 417–425. https://doi.org/10.1046/j.1365-246x.2003.01976.x

Lin, A., Watts, A. B., & Hesselbo, S. (2003). Cenozoic stratigraphy and subsidence history of the South China Sea margin in the Taiwan region. *Basin Research*, **15**(4), 453–478. https://doi.org/10.1046/j.1365-2117.2003.00215.x

Lin, Y.-P., Zhao, L., & Hung, S.-H. (2011). Assessment of tomography models of Taiwan using first-arrival times from the TAIGER active-source experiment. *Bulletin of the Seismological Society of America*, **101**(2), 866–880. https://doi.org/10.1785/0120100244

Lux, I., & Koblinger (2018). *Monte Carlo particle transport methods*. CRC Press.

Margerin, L. (2017). Breakdown of equipartition in diffuse fields caused by energy leakage. *The European Physical Journal: Special Topics*, **226**(7), 1353–1370. https://doi.org/10.1140/epjst/e2016-60165-6

Margerin, L., Campillo, M., & Tiggelen, B. (1998). Radiative transfer and diffusion of waves in a layered medium: New insight into coda Q. *Geophysical Journal International*, **134**(2), 596–612. https://doi.org/10.1111/j.1365-246x.1998.tb07142.x

Margerin, L., & Nolet, G. (2003). Multiple scattering of high-frequency seismic waves in the deep Earth: PKP precursor analysis and inversion for mantle granularity. *Journal of Geophysical Research*, **108**(B11), 2514, https://doi.org/10.1029/2003jb002455

Margerin, L., Planès, T., Mayor, J., & Calvet, M. (2016). Sensitivity kernels for coda-wave interferometry and scattering tomography: Theory and numerical evaluation in two-dimensional anisotropically scattering media. *Geophysical Journal International*, **204**(1), 650–666. https://doi.org/10.1093/gji/ggv470

Mavko, G., Mukerji, T., & Dvorkin, J. (2020). *The rock physics handbook*. Cambridge University Press.

Mayeda, K., Koyanagi, S., Hoshiwa, M., Aki, K., & Zeng, Y. (1992). A comparative study of scattering, intrinsic, and coda Q-1 for Hawaii, Long Valley, and central California between 1.5 and 15.0 Hz. *Journal of Geophysical Research*, **97**(B5), 6643–6659. https://doi.org/10.1029/91jb03094

Mayor, J., Calvet, M., Margerin, L., Vanderhaeghe, O., & Traversa, P. (2016). Crustal structure of the Alps as seen by attenuation tomography. *Earth and Planetary Science Letters*, **439**, 71–80. https://doi.org/10.1016/j.epsl.2016.01.025

Mayor, J., Margerin, L., & Calvet, M. (2014). Sensitivity of coda waves to spatial variations of absorption and scattering: Radiative transfer theory and 2-D examples. *Geophysical Journal International*, **197**(2), 1117–1137. https://doi.org/10.1093/gji/ggu046

Nishigami, K. (1997). Spatial distribution of coda scatterers in the crust around two active volcanoes and one active fault system in central Japan: Inversion analysis of coda envelope. *Physics of the Earth and Planetary Interiors*, **104**(1–3), 75–89. https://doi.org/10.1016/s0031-9201(97)00058-7

Nishigami, K. (2000). Deep crustal heterogeneity along and around the San Andreas fault system in central California and its relation to the segmentation. *Journal of Geophysical Research*, **105**(B4), 7983–7998. https://doi.org/10.1029/1999jb900381

Ogiso, M. (2019). A method for mapping intrinsic attenuation factors and scattering coefficients of S waves in 3-D space and its application in southwestern Japan. *Geophysical Journal International*, **216**(2), 948–957. https://doi.org/10.1093/gji/ggy468

Rachman, A. N., Chung, T. W., Yoshimoto, K., & Son, B. (2015). Separation of intrinsic and scattering attenuation using single event source in South Korea. *Bulletin of the Seismological Society of America*, **105**(2A), 858–872. https://doi.org/10.1785/0120140259

Roth, E. G., Wiens, D. A., Dorman, L. M., Hildebrand, J., & Webb, S. C. (1999). Seismic attenuation tomography of the Tonga-Fiji region using phase pair methods. *Journal of Geophysical Research*, **104**(B3), 4795–4809. https://doi.org/10.1029/1998jb900052

Sato, H. (1982). Amplitude attenuation of impulsive waves in random media based on travel time corrected mean wave formalism. *Journal of the Acoustical Society of America*, **71**(3), 559–564. https://doi.org/10.1121/1.387525

Sato, H. (2019). Isotropic scattering coefficient of the solid Earth. *Geophysical Journal International*, **218**(3), 2079–2088. https://doi.org/10.1093/gji/ggz266

Sato, H., Fehler, M. C., & Maeda, T. (2012). *Seismic wave propagation and scattering in the heterogeneous earth*. Springer.

Sens-Schönfelder, C., Margerin, L., & Campillo, M. (2009). Laterally heterogeneous scattering explains Lg blockage in the Pyrenees. *Journal of Geophysical Research*, **114**(B7), B07309. https://doi.org/10.1029/2008JB006107

Sens-Schönfelder, C., & Wegler, U. (2006). Radiative transfer theory for estimation of the seismic moment. *Geophysical Journal International*, **167**(3), 1363–1372. https://doi.org/10.1111/j.1365-246x.2006.03139.x

Shin, T.-C., Chang, C.-H., Pu, H.-C., Hsiao-Wei, L., & Leu, P.-L. (2013). The geophysical database management system in Taiwan. *TAO: Terrestrial, Atmospheric and Oceanic Sciences*, **24**(1), 11. https://doi.org/10.3319/tao.2012.09.20.01(t)

Shyu, J. B. H., Chuang, Y.-R., Chen, Y.-L., Lee, Y.-R., & Cheng, C.-T. (2016). A new on-land seismogenic structure source database from the Taiwan earthquake model (tem) project for seismic hazard analysis of Taiwan. *Terrestrial, Atmospheric and Oceanic Sciences*, **27**(3), 311. https://doi.org/10.3319/tao.2015.11.27.02(tem)

Takahashi, T., Sato, H., Nishimura, T., & Obara, K. (2009). Tomographic inversion of the peak delay times to reveal random velocity fluctuations in the lithosphere: Method and application to northeastern Japan. *Geophysical Journal International*, **178**(3), 1437–1455. https://doi.org/10.1111/j.1365-246x.2009.04227.x

Takemura, S., & Yoshimoto, K. (2014). Strong seismic wave scattering in the low-velocity anomaly associated with subduction of oceanic plate. *Geophysical Journal International*, **197**(2), 1016–1032. https://doi.org/10.1093/gji/ggu031

Takeuchi, N. (2016). Differential Monte Carlo method for computing seismogram envelopes and their partial derivatives. *Journal of Geophysical Research: Solid Earth*, **121**(5), 3428–3444. https://doi.org/10.1002/2015jb012661

Teng, L. S. (1990). Geotectonic evolution of late Cenozoic arc-continent collision in Taiwan. *Tectonophysics*, **183**(1–4), 57–76. https://doi.org/10.1016/0040-1951(90)90188-e

Teng, L. S. (1996). Extensional collapse of the northern Taiwan mountain belt. *Geology*, **24**(10), 949–952. https://doi.org/10.1130/0091-7613(1996)024<0949:ecotnt>2.3.co;2

Teng, L. S., & Lin, A. T. (2004). Cenozoic tectonics of the China continental margin: Insights from Taiwan. *Geological Society, London, Special Publications*, **226**(1), 313–332. https://doi.org/10.1144/gsl.sp.2004.226.01.17

Ugalde, A., Carcolé, E., & Vargas, C. A. (2010). S-wave attenuation characteristics in the Galeras volcanic complex (south western Colombia). *Physics of the Earth and Planetary Interiors*, **181**(3–4), 73–81. https://doi.org/10.1016/j.pepi.2010.04.009
Wang, K.-L., Chung, S.-L., Chen, C.-H., Shinjo, R., Yang, T. F., & Chen, C.-H. (1999). Post-collisional magmatism around northern Taiwan and its relation with opening of the Okinawa trough. *Tectonophysics*, **308**(3), 363-376. https://doi.org/10.1016/s0040-1951(99)00111-0

Wang, Y.-J., Ma, K.-F., Mouthereau, F., & Eberhart-Phillips, D. (2010). Three-dimensional *Q*<sub>p</sub>-and *Q*<sub>t</sub>-tomography beneath taiwan orogenic belt: Implications for tectonic and thermal structure. *Geophysical Journal International*, **180**(2), 891-910. https://doi.org/10.1111/j.1365-246x.2009.04459.x

Wegler, U. (2003). Analysis of multiple scattering at Vesuvius volcano, Italy, using data of the TomoVes active seismic experiment. *Journal of Volcanology and Geothermal Research*, **128**(1-3), 45-63. https://doi.org/10.1016/s0377-0273(03)00246-4

Wu, R.-S. (1982). Attenuation of short period seismic waves due to scattering. *Geophysical Research Letters*, **9**(1), 9-12. https://doi.org/10.1029/gl009i001p00009

Wu, R.-S., & Aki, K. (1985). Scattering characteristics of elastic waves by an elastic heterogeneity. *Geophysics*, **50**(4), 582-595. https://doi.org/10.1190/1.1441934

Yu, S.-B., Chen, H.-Y., & Kuo, L.-C. (1997). Velocity field of GPS stations in the Taiwan area. *Tectonophysics*, **274**(1-3), 41-59. https://doi.org/10.1016/s0040-1951(96)00297-1

Zeng, Y. (2017). Modeling of high-frequency seismic-wave scattering and propagation using radiative transfer theory. *Bulletin of the Seismological Society of America*, **107**(6), 2948-2962. https://doi.org/10.1785/0120160241

Zhang, F., & Papageorgiou, A. S. (2010). Attenuation characteristics of Taiwan: Estimation of coda Q, S-wave Q, scattering Q, intrinsic Q, and scattering coefficient. *Seismological Research Letters*, **81**(5), 769-777. https://doi.org/10.1785/gssrl.81.5.769