Vertical high-velocity structures and seismic activity in western Shandong Rise, China: Case study inspired by double-difference seismic tomography

1 Introduction

Shandong, located in the eastern part of mainland China, is divided into eastern and western regions by the Tancheng–Lujiang fault zone oriented in the NNE direction (Figure 1a). The eastern part is primarily composed of the Jiaobei uplift, Jiaolai basin, and Jiaonan uplift (Figure 1b). The western region of Shandong, referred to as the Luxi region is further divided into the Luxi uplift zone and Jiyang depression zone by the Qihe-Guangrao fault [1,2]. The Jiyang depression is a significant oil-bearing basin in China, making the western Shandong area a focus of attention among geologists due to its unique geological features and abundant mineral resources. The Luxi uplift zone is one of the most important tectonic units in the eastern part of the North China Craton (NCC), bounded by the Tancheng-Lujiang fault to the east, the Liaohe-Lankao fault to the west, the Chengning uplift to the north, and the Fengpei fault to the south (Figure 1b). Surrounding the Luxi uplift zone are basins or plains, and the internal main fault structures in the basins on both the eastern and western sides are oriented in the NE-NNE direction. In contrast, the western Shandong area is a topographic uplift zone, with its main faults oriented in the NW–NWW direction, controlling the distribution of Cenozoic and Neogene sedimentary landforms, topography, and stratigraphic sequences [3,4]. We collected the focal mechanisms of 49 M ≥ 6.5 earthquakes in and around the study area, which are mainly of the strike-slip type, as shown in Figure 1a [5].

Figure 1

Earthquakes and tectonics on a topographic map. (a) Seismicity in the study region (black pane) and surrounding regions. The focal mechanisms were obtained from the previous studies [5]. (b) Simplified geological map of Western Shandong, modified after [4,6]. The black lines indicate the faults. F1: Taishan-Tongyedian fault; F2: Xintai-Duozhuang fault; F3: Mengshan fault; F4: Wensi fault; F5: Ganlin fault and F6: Nishan fault. The thin black dotted line represents depression. B1: Laiwu Basin; B2: Feicheng Basin; B3: Wenkou Basin; B4: Mengyin Basin; B5: Ningyang Basin; B6: Sishui Basin; and B7: Pingyi Basin; B8: Jiyang depression.

The Luxi area is situated in the eastern part of the NCC, and belongs to the first-level tectonic unit within the NCC. It has a typical North China-type crystalline basement structure. The crystalline basement comprises metamorphic rocks of the Archean Taishan Complex, while the sedimentary cover includes Paleozoic, Mesozoic, and Cenozoic carbonate rocks or clastic rocks. In particular, the Paleozoic strata are predominantly marine carbonate rocks, while the Mesozoic and Cenozoic basins are characterized by terrestrial clastic rocks, volcaniclastic rocks, and lacustrine mudstones [7,8]. The western Shandong region has undergone a highly complex tectonic evolution, particularly since the Late Mesozoic. Tectonic activity in the study area has been dynamic, coupled with regional geotectonic changes and structural deformations in surrounding basin groups. The region has primarily experienced three tectonic evolution stages during the Indosinian, Yanshanian, and Himalayan periods [9,10]. The Indosinian period is the initial stage of the formation of the geotectonic pattern of the Chinese continent. In the Middle and Late Triassic, the Yangzi plate was subducted along the NNE to the North China plate, which suffered from uplift and denudation under the strong near-SN extrusion and produced a series of reverse faults spreading in the near-EW direction [9,11]. Since the Early-Middle Jurassic, as the collision and extrusion activities between the Yangzi Plate and the North China Plate began to gradually weaken, North China entered the Yanshan Movement stage. In the Late Jurassic-Early Cretaceous, the Izanagi Plate began to subduct the Eurasian Plate at a high speed and low angle, resulting in mantle upwelling in the Early Cretaceous, and the study area changed from a weakly extruded to a back-arc extensional environment [12]. During the Himalayan orogeny, under the subduction of the Pacific Plate, the remote effect of the collision of the Indian Plate on the Eurasian Plate, and the conversion of the Tan-Lu fault from left-rotation strike-slip to right-rotation strike-slip and other stresses, the lithosphere of eastern China was thinned, and the Luxi Rise was rapidly uplifted at this time, and at the end of the Eocene, the Rise and the Jiyang Depression separated from each other [13,14,15]. The Luxi region development and evolution of Mesozoic–Cenozoic extensional structures share a similar tectonic background with the NCC. The development is mainly controlled by three stress interactions: interplate interactions along the periphery, manifested by the post-arc extension caused by the subduction of the Western Pacific Plate beneath the Eurasian Plate, forming a trench-arc-basin system [16,17]; local stress fields generated by boundary faults, particularly the strike-slip activity along the Tancheng–Lujiang Fault Zone in the eastern boundary, inducing block movements and generating an extensional-strike-slip deformation field in its hinterland [18,19]; and the gravitational tension on the overlying lithosphere generated by upwelling and convection of deep mantle material [20]. Since the Mesozoic, the crust in the western Shandong region has undergone compression, extension, shearing, and other tectonic processes, resulting in the structural combination of the Luxi uplift and the Jiyang depression, forming the overall tectonic pattern of the western Shandong region. Therefore, the Luxi region has become one of the optimal study areas for researching continental internal dynamics.

In recent years, significant progress has been made in geophysical exploration in the eastern part of China [21,22,23]. Abundant exploration results have provided a basis for the study of lithospheric structure and deep structures in the NCC, also supporting the evidence for lithospheric thinning in the Mesozoic in the North China region. Zhang inferred the thickness of the lithosphere in the Bohai Bay Basin and surrounding areas using geothermal field data and electromagnetic deep sounding, revealing the correspondence between the basin and the crust-mantle regional uplift [24]. Gong et al. calculated the current average geothermal gradient of the Jiyang depression based on temperature data from over 700 wells, obtaining a value of 35.5°C/km [25]. Jiang et al. used large-scale Bouguer gravity data to invert the crustal structure of the Bohai Bay Basin and adjacent areas [26]. The middle and upper crust are controlled by NNE-oriented structures, and many faults extend down to the middle crust. Shear wave splitting is one of the effective means to study the anisotropic characteristics of crustal media [27]. The fast-wave polarization in the crust is always a few degrees above ten with the SKS fast-wave in the Shandong area. This difference may show the decoupling of crustal movement and mantle deformation in the Shandong region [28]． Seismic tomography can obtain the velocity structure characteristics of the study area, which is important for understanding the geodynamic background and the environment of seismicity [29]. The S-wave velocity structure at a resolution of 2° × 2° reveals key evidence that the lithosphere of the overlying plate in the NCC region has been subjected to thermochemical erosion/transformation, leading to lithospheric destruction [30]. Tang et al. analyzed the crustal structure of the Jiyang depression based on a V _p velocity structure diagram along the Yiyuan–Laoling–Dacheng profile [31]. They identified a three-layered structure: the upper crust with V _p ranging from 5.7 to 6.3 km/s and V _s from 3.0 to 3.7 km/s, dominated by brittle fractures; the middle crust with V _p ranging from 6.45 to 6.55 km/s and V _s from 3.8 to 3.9 km/s, representing a low-velocity zone, suggesting the possible presence of magma chambers, magma sheets, or fluids, and serving as the boundary layer where the upper crustal fractures detach [31,32]. Based on a wide-angle reflection seismic profile, Wang et al. identified near-vertical, high-velocity anomalies in the crust of western Shandong [33]. Wang et al. inverted the 3D velocity structure of the crust and upper mantle in NCC and analyzed the differences in velocity structure between different tectonic blocks [34]. Li et al. utilized natural seismic P-wave arrival time data in the circum-Bohai Sea region. Employing a grid with dimensions of 0.5° × 0.6° in latitude and longitude, they inverted the three-dimensional P-wave velocity structure of the crust and upper mantle in the region [35]. Su et al. inverted the P-wave velocity structure of the Shandong region using a 0.5° × 0.5° grid. They identified large-scale low-velocity anomalies at depths below 20 km in the Luxi fault block, which is related to the large-scale detachment tectonics that occurred in the Eocene [36].

The deployment of numerous digitized seismic stations and the acquisition of high-precision seismic wave time arrival data have positioned seismic stratigraphic imaging as a crucial tool for investigating the structure of the Earth’s crust and mantle. With strong tectonic activity and active magmatic activity, the Luxi area is an important study area for analyzing the destruction of the NCC. Most of the previous research results have focused on large-scale tomography of the NCC, and the resolution of the velocity structures is only about 0.5° due to constraints such as station distribution and data sample limitations. A comprehensive understanding of the crustal velocity structure in the Luxi area is still insufficient, and there is a need to improve the imaging resolution. The systematic geophysical interpretation of the deep structure, crustal morphology, and media properties in the Luxi area needs further clarification and discussion. In this work, we have made great efforts to collect a large number of high-quality seismic phase data that covers the Luxi area densely and uniformly. We applied the advanced double-difference tomography method to obtain a more detailed model of the three-dimensional P- and S-wave velocity structure of the Luxi region. Focusing on the characteristics of the crustal structure of typical tectonic structures such as uplifts and faults in the Luxi region, the deep tectonic environment, and controlling factors of seismicity, this provides new information for understanding the deep structure and dynamic evolution of the Luxi region.

2 Data and method

In this article, we collected the seismic phase data of earthquakes above magnitude 1.0 recorded by 128 seismometers in the research region (115°–120°E, 34°–39°N) and covered the period from January 2008 to October 2023 provided by China Digital Seismic Network (Figure 2a). After excluding non-natural seismic events such as explosions and collapses, a total of 2,833 natural earthquake events were identified. To ensure the reliability of the inversion results, a preliminary screening was applied to the obtained phase arrival times. This involved selecting events with travel time residuals (TRES) TRES ≤ 0.5 s and requiring a minimum of 8 seismic stations to have observed each seismic event. Furthermore, based on the travel time distance curves, we eliminated the data that deviated significantly from the theoretical travel time. Figure 2b shows the travel-time distance curves of P-wave and S-wave. Then, we use ph2dt executable program in the HypoDD software to combine seismic events into seismic pairs and extract relative travel-time data [37]. It is required that the interval between event pairs should be within the range of 0.1–10 km and the maximum distance between events and stations should be no more than 500 km. Finally, we obtained a total of 26,351 absolute P- and 26,349 absolute S-wave arrivals and 99,627 P- and 99,625 S-differential arrival times from 2,429 local earthquakes.

Figure 2

Seismic events and initial velocity model in this study. (a) Spatial distribution of the data set and five profile location. Red triangles indicate seismic stations. (b) Travel time distance curves of P and S wave. (c) The one-dimensional initial velocity model.

The tomoDD is an efficient method for characterizing the local and regional velocity structures and earthquake hypocenters [38,39]. As the name suggests, the tomoDD integrates the double-difference seismic location technique into tomographic inversion. This approach involves jointly inverting seismic source location parameters and three-dimensional velocity structures by utilizing both relative travel time data and absolute travel time data. In accordance with ray theory, the residual between the observed arrival time at seismic station k from seismic source i and the theoretical arrival time is expressed as:

where T k i is the observed time of arrival of the body wave from the earthquake source i to the seismic observation station k. Δ m = ( Δ x , Δ y , Δ z ) signifies variations or disturbances in the coordinates of the seismic source location. ∆ τ represents variations or disturbances in the timing of the seismic event. δ u denotes perturbations in slowness. d l is an element of ray path length. Assuming that these two events are near each other, the paths from the events to a common station are almost identical and the velocity structure is known. For a pair of events i and j recorded by stations k, the double difference d r k ij can be simplified as:

Equation (1) is employed to establish the large-scale velocity structure of the study area and the absolute seismic source location using absolute arrival time data. Equation (2) is utilized to characterize the fine-scale velocity structure near the seismic source zone and the relative location based on arrival time data. In the double-difference relocation algorithm, it is assumed that neighboring events arriving at the same station share similar paths. However, when the distance between neighboring events exceeds the scale of velocity changes, the differences in paths caused by velocity heterogeneity become pronounced with significant variations in the seismic source locations. The tomoDD incorporates absolute arrival time data and accounts for spatial variations in the medium’s velocity structure, reducing errors caused by assuming constant velocity between station–event pairs in double-difference relocation. This approach allows for greater accuracy in seismic event localization, without imposing distance constraints on earthquake pairs [17,38].

3 Inversion details and regularization parameter

The double-difference tomography algorithm relies fundamentally on linear inversion, emphasizing the critical necessity of a robust initial velocity model [40]. In order to minimize the influence of the initial one-dimensional velocity model on the inversion results, a comprehensive analysis was performed based on the initial model and results used in previous seismic stratigraphic imaging studies in the region [33,34,41]. Finally, we employed the wide-angle reflection seismic profile across the western Shandong region, specifically the Liaocheng-Lianyungang, as the initial P-wave velocity model [33]. Within the depth range of 0–30 km, we set the P-wave to S-wave velocity ratio to 1.73. The initial P and S velocity model is shown in Figure 2c. The tomoDD is based on parameterizing the model using a three-dimensional regular grid. We transformed the initial one-dimensional velocity model into a series of 3D grid points. The grid nodes of the velocity model are spaced at intervals of 0.2° in both latitude and longitude directions (shown in Figure 3 in the blue cross grid), and inversion nodes in the vertical direction are located at 2, 6, 10, 15, 20, and 25 km depth. To minimize occurrences of seismic events at depths less than 0 km, the depth of the uppermost layer is set to −100 km, with a P-wave velocity assigned to 0.3 km/s. In the current version of tomoDD, we use the pseudo-bending ray-tracing algorithm to find the rays and calculate the travel-times between events and stations [42]. The model is represented as a regular set of 3-D nodes, and the velocity values are interpolated by using the trilinear interpolation method.

Figure 3

Seismic ray paths of P-wave (a) and S-wave (b) in the studied. Red triangles indicate seismic stations, green solid dots indicate earthquake locations, and blue crosses indicate inversion grid points.

Two types of data, namely absolute arrival times and catalog differential arrival times, are utilized in the inversion process. To amalgamate these data types into unified data, we implement a hierarchical weighting scheme akin to hypoDD [43]. Initially, greater weight is assigned to the absolute data, facilitating the determination of large-scale velocity structures. Subsequently, we systematically diminish the weight of the absolute data while concurrently increasing the weight assigned to the differential data. This approach enables the extraction of fine-scale velocity structures within the source region. The tomoDD algorithm alternates between joint inversion of hypocenter parameters and velocity; however, the convergence rate of the velocity structure is faster than that of the hypocenter. To balance the inversion coordination between the two, it is customary to add a separate positioning inversion after each joint inversion. In this article, a total of 10 iterations are carried out in the whole inversion process, with two simultaneous seismic relocation and velocity tomography and one seismic relocation only performing alternately. The tomoDD inverts iteratively for hypocenter locations and velocity structure using Least Square QR factorization (LSQR) [44]. It takes the second-order norm of the travel-time residual as the objective function and iterates many times until a stable solution is obtained.

In the tomoDD, the LSQR algorithm is employed, utilizing the two-norm of the total travel-time residuals as the objective function for an iterative solution. The trade-off of damping and smoothing weight parameters has a greater impact on the stability of the inversion results [45,46,47]. Therefore, in order to make the inversion more stable, we tested the damping factor (from 10 to 1,500) and the smoothing weight factor (from 0.001 to 10,000), and the equilibrium curve of the data variance and the model variance is drawn. The smoothing factor is utilized to restrict variations in slowness. Consequently, a balance analysis is conducted by weighing the two-norm of slowness variation against the two-norm of total travel-time residuals (refer to Figure 4a). The damping parameter concurrently imposes constraints on changes in both earthquake locations and slowness. Hence, a balance analysis for the damping parameter is carried out by considering the two-norms of earthquake locations, slowness variation, and total travel-time residuals (refer to Figure 4b). The optimal regularization parameter corresponds to the point on the L-curve with the maximum curvature, referred to as the L-corner. We solve for the numerical gradient of the discrete points and then find the point with the largest change in the maximum gradient, from which we determine the location of the L-corner. In this article, we have ultimately selected a smoothing factor of 40 and a damping factor of 150 as control parameters for the inversion process.

Figure 4

Trade-off curves of (a) smoothing and (b) damping weight parameters. The horizontal axis represents the normalized norm of slowness variation, and the vertical axis represents the normalized norm of the overall travel-time residual.

4 Results and discussion

The resolution of inversion results is influenced by factors such as the spatial distribution of earthquakes, station coverage, seismic ray density, and intersections among rays. To assess the reliability and spatial resolution capability of the velocity structure inversion, we conducted a checkerboard test [48,49]. In this test, velocity model nodes were arranged in a spatially alternating “checkerboard” pattern, and random errors were introduced during the computation of theoretical travel times to examine the stability under this model configuration. Positive and negative velocity anomalies of ±5% are assigned to the grid nodes adjacent to each other for P and S velocity. Figures 5 and 6 show the results of the checkerboard tests for P- and S-wave tomography with a grid spacing of 0.2° × 0.2°. Both the checkerboard pattern and velocity anomaly amplitudes at different depths are well reconstructed. In contrast, the northern part of the study area is not well restored to the checkerboard test. The northern segment of the study area is predominantly characterized by a faulted basin featuring a substantial cover layer. Compounded by a limited number of permanent stations (as depicted in Figure 2a), the recorded seismic phases from existing stations are notably unclear, resulting in poor data quality and the consequent absence of seismic phase travel time data. This absence, in turn, diminishes the resolution for the area. The distribution of seismic wave rays in Figure 3 also shows a clear lack of ray coverage in the northern region. In order to quantitatively assess the reliability of the inversion results, we chose the Distribution of Derivative Weight Sum (DWS) parameters to evaluate. The parameter DWS is an indicator of the average relative ray density surrounding a given node. Evidently, the relationship graph between model resolution values and DWS values illustrates that higher DWS values correspond to a denser distribution of rays, thereby leading to enhanced resolution [39,50]. In this study, we adopted parameter selections from prior research on double-difference tomography, and in conjunction with chessboard test results, we ultimately opted to retain regions where the DWS values exceeded 100. This meticulous criterion ensures the reliability of our imaging results [51,52,53,54].

Figure 5

Results of a checkerboard resolution test for P-wave tomography. The velocity perturbation scale is shown at the right. The amplitude of velocity perturbations is 5 percent in the input model. The layer depth is shown at the lower-left corner of each map. The white line indicates the provincial boundary of Shandong Province. The black contours represent DWS > 100 distribution.

Figure 6

The same as Figure 5 but for the S-wave checkerboard test.

Figure 7 shows the RMS histogram of travel-time residuals before and after earthquake relocation. Notably, seismic phase travel time residuals exhibit a considerable increase for stations located in the northern part of the study area relative to those in the southern part. This discrepancy is primarily attributed to the thicker cover layer in the northern basin. Importantly, there is a discernible and significant reduction in residuals for each station after relocation. In addition, Figure 7b and c give the distribution of the traveled residuals for the P- and S-waves, respectively. Before the inversion, the travel-time residuals are mainly in the range of −2 to 2 s. After the inversion, most of the seismic wave travel-time residuals are concentrated in a narrow range of −0.3 to 0.3 s. The root mean square (rms) P-wave travel-time residual is reduced to 0.06 s from 0.54 s. Similarly, the rms of the S-wave travel-time is reduced from 0.71 to 0.10 s.

Figure 7

(a) The travel-time errors of seismic stations before and after inversion; (b) The histogram of P-wave and (c) S-wave travel-time residuals before (gray) and after (colored) inversion.

4.1 Crustal structural features

Figures 8 and 9 show the results of high-resolution 3D stratigraphic imaging of P- and S-wave velocities within the middle and upper crust of western Shandong. Only the distribution of velocity structures with DWS >100 is plotted in the figures. These maps also depict the distribution of earthquakes that occurred within a layer of 10 km thickness bound to each slice. From the overall analysis, the crustal velocity structure in the study area shows obvious unevenness in the horizontal directions. The distribution of surface basins and uplifted areas often correlates with variations in seismic velocities within the crust. The results of the horizontal stratigraphic imaging at a depth of 2 km represent a shallow velocity structure that agrees well with known geologic formations. The delineation between the uplift zone and the basin is relatively obvious. For instance, in Figure 1b, the Laiwu Basin, Feicheng Basin, and Ningyang Basin all exhibit discernible low-velocity body anomalies, a characteristic attributed to the substantial thickness of sedimentary layers. In order to more intuitively recognize the relationship between seismic wave velocity structure and geomorphology, we obtained the vertical velocity structure distribution of the CC’ profile along the longitudinal direction vertically from south to north through the Luxi region, as shown in Figure 10CC’. We can clearly see that the faulted basins correspond to low-velocity zones within the crust, while there is a clear intrusion of high-velocity anomalies beneath Taishan and Pingyi. The Liaocheng–Lianyungang deep seismic sounding (DSS) profile reveals the presence of a nearly upright, high-speed anomaly in the crustal velocity structure in the Luxi region [33,55]. We found velocity structure features similar to the DSS using natural seismic tomography. We extracted the vertical cross-sections of velocity by following the location of the DSS profile from Liaocheng to Lianyungang (as shown in Figure 10AA′). The crustal structure section (Figure 10AA′) shows that there are near-upright high velocity bodies in the crust, which extend from the lower crust into the upper crust, and the top reaches almost to the base of the sedimentary basin. The high-velocity bodies are about 40 km wide between Liaocheng and Sishui, and about 60 km wide between Sishui and Linyi. In addition, our results reveal the existence of vertical high-velocity bodies beneath Linyi with a width of about 25 km. The other profiles are also characterized by high-velocity body intrusion. These high-velocity bodies may be taken as the seismological evidence for the uprising of lower crust material into the upper crust and reaching the bottom of the crystalline basement.

Figure 8

Map views of V _p structures. The layer depth is shown at the lower-left corner of each map. White dots represent small earthquakes since 2008, and solid gray dots indicate earthquakes of magnitude 5.5 or greater since 780 B. C. White circles denote the earthquakes that occurred within a layer of 10 km thickness bound to each slice.

Figure 9

The same as Figure 8 but for V _s structures.

Figure 10

Vertical cross-sections of P-wave velocity (V _p) and S-wave velocity (V _s) tomography along the five profiles shown in Figure 2a. The gray spheres represent earthquakes of magnitude M ≥ 3.5 within 25 km on either side of the profile. White dots represent earthquakes of magnitude M1.0–3.4. The color scale in the lower right corner indicates elevation and seismic wave velocity perturbation, respectively.

Field survey performed in the Luxi region indicates that extensional faults are well developed and characterized by complicated structural styles. The main structure style is the compound extensional fault system composed of high-angle faults and low-angle décollement faults [4,6]. The décollement depths of different faults in the western part of China are different, with a general trend of deepening from south to north [6,56]. There exists a low-velocity layer about 3–4 km thick in the upper and middle crust, which is highly plastic and is an important décollement detachment surface in the Luxi area [56,57]. Previous studies have shown that a large-scale low-velocity anomaly exists at 15–20 km depth [3,56]. Our imaging results indicate the presence of a significant low-velocity anomaly at depths of 13–17 km, which is in general agreement with the location of the high-conductivity layer obtained by magnetotelluric sounding [56,58,59,60]. In addition, the AA′–DD′ profile clearly observed the presence of high-velocity bodies in the 7–13 km depth range. Previous research findings indicate the presence of a distinct lithospheric layer-cake structure in the Luxi region. At the bottom of the upper crust, specifically at depths ranging from 7–11 km, there is identified a brittle–ductile transition zone [59]. The high-velocity anomalous zone corresponds essentially to the identified brittle–ductile transition zone. Within this transition zone, the upper segment is characterized by a brittle deformation structural layer, while the lower segment exhibits properties indicative of a tough deformation structural layer. The presence of this transition zone assumes significance in facilitating the transfer of lithospheric stress from the deeper regions to the shallower domains.

5 Conclusions

We used advanced double-difference tomography to determine the 3-D velocity structure of the P- and S-waves under the Luxi region. The results of seismic stratigraphic imaging in this article reveal the prevalence of near-vertical high-velocity anomalies in the Luxi region, reflecting the possibility of intrusion of mantle material into the crust. There is a transition zone of high and low velocities at a depth of 7–13 km, which is also a brittle and ductile transition zone, where most of the historical earthquakes above M3.5 occurred. Meanwhile, our imaging results also clearly indicate the presence of a significant low-velocity anomaly at depths of 13–17 km, which corresponds to the presence of a high-conductivity layer at that depth obtained from magnetotelluric sounding. The formation of décollement structures in the western Luxi region is related not only to the activity of the Tan-Lu fault zone, but also to the thinning of the lithosphere due to the subduction of the Pacific slab and the small-scale convection caused by the upwelling of the mantle. Combined with the relationship between the distribution of historical earthquakes and the velocity structure, most of the moderate and strong earthquakes occurred in the high- and low-velocity transition zone, which corresponds to the brittle–ductile transition zone. Exploration of the distribution and physical properties of such transition zones will yield important insights into the medium- and long-term prediction of earthquakes.