# Bulletin of the Seismological Society of America

- Copyright © 2006 Bulletin of the Seismological Society of America

## Abstract

The attenuation of shear waves propagating in the crust of northwestern Turkey has been investigated in the frequency range 1–10 Hz. A standard spectral inversion scheme is applied to a data set of 245 aftershocks (*M*_{L} <4.5) of the 1999 Izmit earthquake. The obtained attenuation-with-distance curves have been described in terms of the *t** cumulative attenuation parameter and its dependence on frequency and distance investigated. At 1 Hz, *Q*^{−1}, evaluated by normalizing *t** to the travel time, is generally larger than 0.025 for source-to-station distances smaller than 40 km, indicating the presence of a highly attenuating upper crust in the area. Over longer distances, *Q*^{−1} decreases, suggesting a decrease in the attenuation with depth. By contrast, the normalized *t** computed for earthquakes recorded at stations having almost the same distance from the sources do not show a strong dependence on the backazimuth. These results suggest that the decrease of *Q*^{−1} with depth is more significant than its lateral variations. Regarding its frequency dependence, *Q*^{−1} almost linearly decreases with frequency.

Finally, the near-surface-attenuation parameter *k* is evaluated at 12 stations and the results discussed in terms of site, event, and propagation contributions. The event contribution is not negligible and shows a significant positive correlation with magnitude. The site term is smaller than 0.020 sec for rock or topographic sites, while it assumes values of 0.036 sec and 0.042 sec for two stations installed over thick soft sedimentary layers.

## Introduction

The description and quantification of seismic-wave attenuation plays a fundamental role in many seismological analyses. Attenuation estimates allow the investigation of the structure and physical state of the Earth's interior, and provides corrections that permit a more accurate description of processes at the source. Moreover, studies on hazard assessment require a quantitative evaluation of seismic-wave attenuation in order to predict ground motion from earthquakes.

We present the results of a study on crustal attenuation performed in northwestern Turkey. Previous studies (Bindi *et al.*, 2006, and the reference cited therein) were devoted to characterizing the attenuation in this region by estimating the average quality factor, *Q*(*f*), from the *S*-wave window or using the coda portion of the seismograms. In particular, Akinci, Del Pezzo, and lbanez (1995); Akinci, lbanez, *et al.* (1995); and Akinci *et al.* (2004) investigated the attenuation in the western portion of the North Anatolian fault (naf), while Gündüz *et al.* (1998), Kaslilar-Özcan *et al.* (2002), and Horasan and Boztepe-Güney (2004) analyzed the attenuation properties of the Marmara region.

We analyze the aftershocks (*M*_{L} <4.5) of the Izmit earthquake (*M*_{w} 7.4) that occurred on 17 August 1999. The attenuation-with-distance curves, obtained by applying the generalized inversion technique (Castro *et al.*, 1990), are used to infer the frequency-dependent attenuation characteristics of *S* waves in the frequency range 1–10 Hz. We estimate the *t** parameter (Kanamori, 1967) along each ray path and investigate its dependence on frequency, distance, and backazimuth. The spatial distribution of attenuation is also investigated, carrying out a two-dimensional tomographic inversion for a simplified model in which depth is not accounted for. Finally, the near-surface attenuation is evaluated by computing the *k* (kappa) parameter (Anderson and Hough, 1984). Following Purvance and Anderson (2003), the parameterization of *k* in terms of source, site, and propagation effects is investigated and the results compared with the area's geology, which is characterized by the presence of several alluvial basins (e.g., eeri, 2000) surrounded by complex topography.

## Data

We analyze 245 aftershocks following the 1999 Izmit earthquake, recorded at the Sapanca-Bolu (sabo) and German Task Force (gtf) networks (Grosser *et al.*, 1998; Baumbach *et al.*, 2003). The seismological equipment used in this study consists of Mark L4-3D 1-Hz geophones, 24- bit digitizers with a sampling rate of 100 samples per second, and gps timing. The analyzed earthquakes have magnitudes, *M*_{L}, from 1 to 4.5, the hypocentral distances ranging between 10 and 142 km, with most events having a source depth of 5–15 km. We analyze the spectra of *S*-wave windows 5 sec wide. Each window is cosine tapered (5%) and Fourier transformed. Instrumental corrections are applied and the spectral amplitudes are smoothed using the Konno–Ohmachi window, *b* = 20 (Konno and Ohmachi, 1998). With respect to the pre-event noise, the selected *S* windows have signal-to- noise ratios (s/n) greater than 3 over the frequency band between 1 and 10 Hz. Finally, the two horizontal (east–west and north–south) and the vertical (*Z*) spectra are vectorially summed.

The spectra are collected from 23 stations (Table 1). Figure 1 shows the station locations and the earthquake epicenters. The earthquakes are grouped into three different sets: set *s*1 consists of the earthquakes that occurred between Izmit Bay and Sapanca Lake, including the Gölcük and Sapanca segments of the naf (eeri, 2000); set *s*2 mainly includes the earthquakes that occured south of the Adapazari basin, including the Sakarya segment of the naf; and set *s*3 is the earthquakes that occurred close to the Düzce basin, including part of the Karadere segment of the naf. A sketch of some of the alluvial basins is shown in Figure 1 (eeri, 2000).

## Parametric Description: Attenuation along Single Paths

Following Bindi *et al.* (2006), the attenuation-with- distance curves *A*(*f*,*r*) of *S* waves are obtained by applying the generalized spectral inversion, git (e.g., Castro *et al.*, 1990) to the velocity spectra. The curves *A*(*f*,*r*) can be parameterized in term of attenuation along a single path as (Sanders, 1993) 1 where *G*(*r*) is the geometrical spreading, *t _{ij}* is the travel time relevant to event

*i*recorded at station

*j*, and

*Q*is the apparent quality factor along the ray path. The cumulative attenuation along the path

_{ij}*ij*can be described by introducing the

*t** operator (Kanamori, 1967; Cormier, 1982; among others): 2 In this work, we processed

*t** as a frequency-dependent quantity, and therefore equation (1) can be rewritten as 3 We evaluated the attenuation along all ray paths, assuming

*G*(

*r*) = 10

*r*

^{−1}for

*r*less than 60 km and

*G(r)*= 10

*r*

^{−1.2}for

*r*≥ 60 km, following Bindi

*et al.*(2006). Figures 2, 3, and 4 show the results for three selected stations. Following Haberland and Rietbrock (2001), the

*t** is normalized to the travel time, to represent the average attenuation

*Q*

^{−1}along each ray path, being dependent upon the backazimuth, frequency (top panels), and event location (middle and bottom panels).

The frequency dependency of the attenuation parameter is clear: *Q*^{−1} decreases by about a factor of 10 for frequencies increasing from 1 to 10 Hz. In the polar diagrams of Figures 2–4, the dependence of *Q*^{−1} on distance is also evident: at each frequency, the highest values correspond to the shortest distances (black filled circle, 10 ≤ *r* ≤ 38 km), while the lowest values correspond to the longest paths (empty triangles, 60 ≤ *r* ≤ 80 km). At 1 Hz, *Q*^{−1} for 10 ≤ *r* ≤ 38 km is generally between 0.04 and 0.08 (i.e., 12.5 < *Q* < 25), while for 60 ≤ *r* ≤ 80 km, *Q*^{−1} is generally between 0.001 and 0.04 (25 < *Q* < 1000). The attenuation along each ray path does not show a strong dependence on the backazimuth, and a comparison between the maps showing *Q*^{−1} as a function of the event location (middle and bottom panels of Figs. 2–4) shows that the highest attenuation values correspond to the closest earthquakes. The decrease of *Q*^{−1} with distance is particularly evident by comparing the projection of *Q*^{−1} on a vertical plane for different stations (e.g., station 007 and 026, bottom panels). Comparing *Q*^{−1} evaluated at several stations surrounding the two main clusters of events (e.g., sets *s*2 and *s*3 of Fig. 1), the decrease in *Q*^{−1} with source-to-station distance appears to dominate the variations related to different directions of propagation. Figure 5, top panel, shows *Q*^{−1} at 1 Hz relevant to earthquakes belonging to set *s*3 and recorded at stations 039, 026, and 033, which have different backazimuths, but almost the same distance from the selected earthquakes. The values of *Q*^{−1} range from 0.02 to 0.05 sec, regardless of the azimuth. The bottom panel of Figure 5 shows *Q*^{−1} at 1 Hz against distance for data sets *s*1, *s*2, and *s*3.

The behavior of *Q*^{−1} with distance is similar for the three sets, providing a further indication that lateral contrasts of attenuation play a minor role with respect to source-to- station distance (and therefore the vertical variation of *Q*^{−1}) in determining the average attenuation along each single ray path. The dispersion in the attenuation values at distances in the range 60–80 km is worth attention. This feature could be ascribed to the secondary arrivals whose energy contributes to the direct *S* waves, leading to an apparent attenuation lower than the actual value for of the direct waves, consistent with the results of Bindi *et al.* (2006).

## Spatial Distribution of Attenuation

The analysis of parameter *t** has indicated that seismic attenuation in northwestern Turkey is frequency dependent and varies mainly with depth. Since most of the earthquakes are shallow (hypocentral depths less than 15 km) and clustered into two main groups (sets *s*2 and *s*3 in Fig. 1), and with most of the hypocentral distances within 60 km, a three- dimensional tomographic inversion could provide an accurate reconstruction only for the upper crust in the Adapazari and Düzce regions. Therefore, in order to infer some characteristics of the attenuation in northwestern Turkey within the range 1–10 Hz, we invert the cumulative attenuation along each ray path by using a simplified approach where the *t** values are spread throughout the straight line connecting each epicenter to the recording station. A fully three- dimensional inversion, however, will be the goal of future studies after obtaining a larger data set. The inversion performed in this study aims simply to check whether the spatial distribution of attenuation correlates with the main geologic features of the region. From equations (1) and (2), the anelastic attenuation can be expressed in terms of quality factor and velocity as 4 The source-to-receiver path *r _{ij}* is approximated as a straight line and the traveled medium is approximated as a plane. Discretizing the plane into rectangular pixels with index

*k*having constant

*Q*, equation (4) becomes 5 where

_{k}v_{k}*γ*= 1/

_{k}*Q*, and

_{k}v_{k}*l*is the length of the ray

_{ijk}*ij*in each pixel

*k*, normalized such that

*r*= ∑

_{ij}*, where*

_{k}l_{ijk}*r*is the hypocentral distance. For all stations and sources, and for a given frequency

_{ij}*f*, equation (5) describes a tomographic problem where the unknowns

*γ*must be reconstructed from their projections along some straight lines. The projections are the attenuation curves corrected for the geometrical spreading, and the coefficients of the projection matrix are equal to −

_{k}*πf*times the length

*l*. Several techniques can be exploited to find approximate solutions of (5). We used an iterative method known as the row action maximum likelihood algorithm, ramla (Brown and De Pierro, 1996; Bindi and Caponnetto, 2001), which has the advantage of furnishing nonnegative solutions without any additional constraints, under some quite general conditions. After some tests, we set the number of ramla iterations to 10, and the iterative scheme is initialized with a uniform solution equal to

_{ijk}*γ*= 0.01.

_{k}Figure 6, top-left panel, shows the station locations (triangles), epicenters (black filled circles), and ray-path coverage relevant to the data set considered for the tomographic inversion. The top-right panel shows the parameterization and number of rays crossing each pixel *k*. Each pixel is 0.15 and 0.075 degrees wide in longitude and latitude, respectively. The pixels in region *s*2 (see Fig. 1) have the highest sampling, while some pixels at the edges are not sampled at all (cross symbol). The tomographic reconstruction will not be shown for unsampled pixels.

The middle and bottom panels show the results of the resolution estimate in terms of the point-spread function (e.g., Cong and Mitchell, 1998). Synthetic attenuation values are computed by solving equation (5) in the forward direction by applying the matrix for the actual source and station geometry to a point-input function that assumes a value of 0.02 in a single pixel and zero elsewhere. In particular, we show the results relevant to an input function located in region *s*2 (middle panels) and in region *s*3 (bottom panels). The middle- and bottom-left panels show the rays that sampled the selected pixel (star). The right panels show the results of the tomographic inversion. In both cases, the reconstruction of the point-input function is satisfactory. At the edge of the considered region, a few artifacts with amplitudes between 5% and 50% of the input-point function appear.

Figure 7 shows the results of the tomographic inversion for frequencies 1 and 10 Hz. Before discussing the results obtained with real data, we recall that the values of *γ* in Figure 7 do not represent the attenuation averaged over a certain fixed depth interval, since the sampling of actual rays is not uniform with depth. Therefore, if a shallow high *Q*^{−1} zone exists in a region crossed by only long (and hence deep) rays, this anomaly cannot be imaged in the restored model. On the contrary, if the crust in a certain area is mainly sampled by rays propagating mostly at shallow depth, the resulting values of *γ* will be determined by the *t** values corresponding to these rays. Consequently, if the *Q*^{−1} factor for a certain area varies (decreases) with depth but the rays sample mostly the uppermost layers, the value of *γ* inside the pixels will be dominated by the high *Q*^{−1} of the uppermost layers. The lower *Q*^{−1} values (sampled by a few deep rays) will be mapped either at the edge of the investigated area (where the solution is less constrained) or in zones sampled exclusively by these rays.

The values of *γ* depicted in Figure 7 confirm that a high attenuation is present in the area, in particular in association with the alluvival basins and in the Izmit Bay region. At 1 Hz, values of *Q _{k}* down to 125/

*v*are obtained inside the Adapazari basin, and down to 100/

_{k}*v*and 85/

_{k}*v*inside the Düzce basin and Izmit Bay, respectively. These results also outline the existence of lateral variations of attenuation in the area. Assuming the seismic velocity to be frequency independent,

_{k}*Q*increases with frequency in all basins and reaches the value of about 330/

_{k}*v*at 10 Hz. It is worth noting that high attenuation affects the region surrounding the Izmit Bay at all frequencies. High attenuation for the regions at the borders, as well as the high contrast affecting the maps in Figure 7, might have been emphasized by the simplified model adopted, as previously explained.

_{k}## Near-Surface Attenuation

In this section, the spectral amplitudes at frequencies higher than 10 Hz are analyzed to investigate the effect of local geology on seismic-wave attenuation by computing the *k* parameter (Anderson and Hough, 1984). The *k* parameter is calculated as the high-frequency asymptote of the acceleration spectrum *A*(*f*): 6 where *A*_{0} depends upon the source, epicentral distance, and other factors. Since both source (Papageorgiou and Aki, 1983) and propagation effects (Hanks, 1982) can cause deviations from the flat high-frequency acceleration spectra predicted by the Brune model (Brune, 1970), the interpretation of *k* is usually based on empirical models.

Recently, Purvance and Anderson (2003) showed that 7 is an effective parametrization of *k*. In equation (7), **( r)** describes the distance dependence of

*k*while

*k*

^{site}and

*k*

^{event}are the parameters of the model that capture the dependence of

*k*on the specific recording site and earthquake, respectively. Several empirical observations showed that

*k*

^{site}is generally greater for sites on soft sediments than for sites on rock. It is interpreted as a measure of the attenuation due to the propagation through the shallow subsurface geology;

**(**shows the tendency to increase with increasing distance, and it can be considered as a regional effect due to lateral propagation; and

*r*)*k*

^{event}is mainly a source effect, as shown by Purvance and Anderson (2003), who analyzed strong- motion records obtained by the Guerrero Array (Mexico). They observed that

*k*

^{event}varies systematically with focal mechanism, being lower for normal faulting than for thrust faulting.

Since site amplification, within the frequency range considered for evaluating *k*, could affect the slope of the spectrum decay (Parolai and Bindi, 2004), we discarded stations that exhibited peaks of amplification in the high- frequency range, following the results of Parolai *et al.* (2004). Examples of stations excluded due to amplification peaks at frequencies higher than 10 Hz are stations 004 and 008 in Parolai *et al.* (2004, fig. 5, p. 1102). Eventually, the analysis was performed on a subset of 12 stations. We evaluated the high-frequency decay for 312 recordings from 60 earthquakes triggered by at least 3 of the 12 selected stations.

From equation (7), a linear system is obtained by considering all the *k _{ij}* values evaluated for the

*j*th earthquake recorded at the

*i*th station (e.g., Purvance and Anderson, 2003). A nonparametric description is used for

**(Anderson, 1991) by discretizing the distance range from 0 to 110 km into 22 intervals 5 km wide (**

*(r*_{k})*k*= 0, … , 22). The system is solved constraining

**to be a smooth function of**

*(r*_{k})*r*, and setting equal to zero the site term for two

*a priori*assumed reference stations (039 and 002) installed on rock sites. These stations exhibited almost flat high-frequency spectra and a site response less than 2 in the range 1–20 Hz (Parolai

*et al.*, 2004, fig. 7, p. 1105). Moreover, these stations allow

*k*to be constrained over the whole analyzed distance range. The system is solved performing 50 iterations of the Paige and Saunders (1982) least squares algorithm (lsqr). The bootstrap technique (Efron, 1979) is used for an assessment of the model parameter uncertainties, with the standard error of , and

**computed from 200 replications. The two horizontal components are not composed but used together in the inversion. In Figure 8, the high-frequency spectral fitting for an earthquake recorded at six stations is shown as an example.**

*(r*_{k})

Figure 9 shows the results of the regression (7). The *k*^{event} term ranges from −0.014 sec to 0.019 sec. An analysis we performed (but not reported here) shows that positive and negative values do not correlate with earthquake location. The high scatter in the distribution of *k*^{event} against magnitude (Fig. 9, top-left panel) implies a small linear correlation of *R* = 0.30, but a null hypothesis of no correlation between *k*^{event} and magnitude can be rejected at a 0.98 level of confidence by performing a Student's *t*-test.

Table 1 and Figure 9 (bottom-left panel) show the values of *k*^{site}. The highest values are found for stations installed in areas underlaid by deep sedimentary cover, such as station 013 (*k*^{site} = 0.036 sec) and station 029 (*k*^{site} = 0.042 sec). Rock sites show *k*^{site} values less than 0.010 sec. The dependence of *k* on distance shown in Figure 9 is consistent with the decrease in attenuation with depth observed in the previous sections. For distances smaller than 30 km, ** (r_{k})** increases with distance at a rate of nearly 0.0007 sec/km while, for longer ray paths, its dependency on distance weakens. Finally, for distances larger than 80 km,

**slightly decreases with distance.**

*(r*_{k})## Discussions and Conclusions

The attenuation of seismic waves propagating in the crust of northwestern Turkey has been investigated by applying the generalized inversion scheme to a data set of aftershocks of the 1999 Izmit earthquake. The estimated attenuation-distance curves have been described in terms of the *t** attenuation parameter and its dependency on frequence and distance has been investigated.

At 1 Hz, we found values of *Q*^{−1} =*t**/*T*, where *T* is the travel time, generally greater than 0.025 (i.e., *Q* < 40) for distance up to 40 km, indicating that the seismic waves undergo a strong attenuation in the upper crust. The decrease of *Q*^{−1} with distance suggests a diminishing of the attenuation with depth. Moreover, the attenuation-versus-distance curves estimated from the inversion show a bump for distances between 40 and 60 km, suggesting that secondary arrivals significantly contribute to the spectral amplitudes in this distance range. This issue should be accounted for in the description of attenuation that is adopted in hazard- oriented studies. Attenuation is strongly frequency dependent in the range between 1 to 10 Hz, and it decreases almost linearly with frequency. In the case of using a simplified two-dimensional tomographic inversion, high attenuations are concentrated in the main alluvial basins, such as the Adapazari (region *s*2) and Düzce (region *s*3) basins, where rays spent a large portion of their travel times, and in the region surrounding the Izmit Bay. Coinciding with the sedimentary basins, strong vertical variations in attenuation must be expected.

These results are in agreement with the conclusions of previous studies that showed that several regions of the analyzed area exhibit low and strongly frequency-dependent *Q* values. Examples are *Q _{s}* = 46.59

*f*

^{0.67}for the Bursa region in the distance range 5–60 km, and in the frequency range 0.5–25 Hz (Akyol

*et al.*, 2002);

*Q*= 50.7

_{c}*f*

^{1.01}in western Anatolia considering a lapse time of 30 sec (Akinci

*et al.*, 1994); and

*Q*= 41

_{c}*f*

^{1.08}and

*Q*= 50

_{s}*f*

^{1.09}for the Marmara region in the distance range 20–110 km and the frequency range 1.5–24 Hz (Gündüz

*et al.*, 1998). Recently, Horasan and Boztepe-Güney (2004) found regional differences in the attenuation evaluated for five regions in the Sea of Marmara. Their estimated

*Q*values range from 13

_{s}*f*

^{1.22}to 94

*f*

^{0.83}.

The *Q*^{−1} values, estimated by normalizing the cumulative attenuation parameter *t** to the travel time, show a stronger decrease with distance than a dependence on backazimuth. These results also suggest that, in the analyzed area, the depth variations in attenuation could be larger than the lateral ones. The almost linear decrease of *t** with frequency is also in agreement with the results of Bindi *et al.* (2006).

The near-surface attenuation parameter *k* has been evaluated, and its dependence on source, site, and distance was investigated. Most of the stations show values of *k*^{site} smaller than 0.01 sec. Only for stations 029 and 013 is *k*^{site} greater than 0.030 sec. Station 029 was installed close to the Adapazari basin, near the fault of the 1967 Mudurnu earthquake (*M*_{s} 7.1). Goto and Sawada (2004) proposed a three-dimensional velocity model for the Adapazari basin, consisting of three soil layers and two rock mediums. According to this model (Goto and Sawada, 2004, fig. 4, p. 6), the sediments below station 029 should be between 200 and 400 m thick, with a shear velocity of the order of 200 m/sec. These values are also in agreement with the peak values in the h/v ratios shown in Figure 9. Assuming that *k*^{site} = 0.042 sec at station 029, as determined by the attenuation inside the soft sedimentary cover, a quality factor *Q* ∼25 is found. Station 013 is installed over the thick sedimentary fill of the Yalova area (southern side of the Izmit Bay). The *k*^{site} = 0.036 sec for this station is higher than the values found for rock sites. The estimated *k*^{site} is consistent, although lower, with the average value of 0.056 sec found by Durukal and Catalyurekli (2004) in northwestern Turkey for the sites belonging to class D of the National Earthquake Hazards Reduction Program (nehrp) soil classification (Boore and Joyner, 1997). The same authors found an average value of 0.041 sec for sites corresponding to nehrp class C.

The high attenuation found in northwestern Turkey can have a strong implication on hazard assessments for the area. Recently, Erdik *et al.* (2004) observed that the amplitudes of the recorded ground motions in the Marmara region induced by the Izmit earthquake were lower than those predicted using standard attenuation relations. They suggested that this discrepancy could be related to a source effect. Since the high attenuation found might partially explain the observed discrepancies, the calibration of scaling laws for this region is a matter that deserves future attention.

Following Purvance and Anderson (2003), we introduce a term related to the seismic source for describing the observed *k* values. We found that *k*^{event} correlates with magnitude, in agreement with the results of Purvance and Anderson. In the same region we analyzed, Durukal and Catalyurekli (2004) found a correlation between magnitude and *k* values, although they calculated the latter by averaging, for each event, the values relevant to different stations without isolating the site and propagation contributions. Their average *k* values showed an increase with magnitude at a rate of 0.0038. The best least squares fit of the results in Figure 9 is *k*^{event} = 0.0051 *M*_{L} − 0.014.

Finally, the agreement of the overall tomographic features with the surface geology should stimulate the development of future studies devoted to achieving a more detailed image of the attenuation properties in the region that properly take into account the depth-dependent structure of the region. Such a goal motivates efforts to enlarge the data set, hence making it suitable to perform a three-dimensional tomographic inversion.

## Acknowledgments

The authors express their gratitude to the Hannover Rückversicherung AG for their significant financial support of the field mission. Helpful suggestions were provided by F. Scherbaum. We are thankful to R. Milkereit for drawing some figures using the Generic Mapping Tools software (Wessel and Smith, 2000). K. Fleming kindly improved our English. Comments from two anonymous reviewers also improved the manuscript. Part of this work was conducted during the visits of D. Bindi at the GeoForschungsZentrum Potsdam that were partially funded by the gfz-Potsdam.

- Manuscript received 4 March 2005.