Slip Model for the 25 November 2016 Mw 6.6 Aketao Earthquake, Western China, Revealed by Sentinel-1 and ALOS-2 Observations

Abstract

On 25 November 2016 (UTC 14:24:30), an Mw 6.6 dextral strike-slip earthquake ruptured Aketao county in the northwestern portion of the Kongur Shan extensional system, western China. We extracted surface deformation maps and investigated the distribution of the coseismic slip of the 2016 Aketao earthquake by exploiting the Interferometric Synthetic Aperture Radar data imaged by the Sentinel-1 satellites of the European Space Agency and the ALOS-2 satellite of the Japanese Aerospace Exploration Agency. Assuming the crust of the earth is an elastic half-space homogeneous medium, the best fitting slip model suggests a dip angle of 78° for the seismogenic fault. The rupture of the 2016 Aketao earthquake may have consisted of two sub-events that occurred in rapid succession within a few seconds, resulting in two large discrete asperities with maximum slip of ~0.85 m, which were separated by a ~6 km-wide small slip gap. The maximum slip for the sub-event near the epicenter was mainly concentrated at a depth of ~10 km and that of the other at a depth of ~5 km. The estimated total seismic moment from the optimal slip model is 1.58 × 10¹⁹ N·m, corresponding to an event with a moment magnitude of 6.74. More than 65% of the aftershocks occurred in the areas of increased Coulomb failure stress, in which the stress was estimated to have been increased by at least 0.1 bar. Matching the potential barrier on the fault with the depth distribution of aftershocks implies that friction on the causative fault was heterogeneous, which may play a primary role in controlling the active behavior of the Muji fault.

1. Introduction

On 25 November 2016 (UTC 14:24:30), an Mw 6.6 earthquake jolted Aketao county in northeastern Pamir, western China (Figure 1). This earthquake was the largest to occur beneath mainland China in 2016, and it was followed by ~386 aftershocks of M 0.2–5.0 within five days of the main shock, which were relocated using the double-difference algorithm [1,2]. The largest aftershock, with a magnitude of M 5.0, occurred on 26 November 2016, located to the southeast of the main shock [1]. The solutions from U.S. Geological Survey (USGS) [3] and global Centroid Moment Tensor (GCMT) [4] both suggest that this rupture occurred along a NW-trending fault with a primarily dextral mechanism, which is also attested to by the preliminary finite fault inversion of seismic data [5]. The epicenter of the quake was initially determined to be 74.04°E, 39.27°N, with a ~10 km focal depth. Combined with the regional tectonic background, focal mechanism solutions from the seismic data [3,4,5] suggest that the causative fault possesses a steep dip angle.

Field investigations conducted after the shock indicated that the 2016 Mw 6.6 Aketao earthquake occurred on the Muji fault [1] in the northwestern section of the Kongur Shan extensional system [6]. Due to harsh natural conditions such as extreme terrain and snowy mountains, few seismic and geodetic instruments (e.g., GPS) have been set up surrounding the seismogenic zone. Fortunately, several Interferometric Synthetic Aperture Radar (InSAR) satellites, including Sentinel-1A, Sentinel-1B and ALOS-2, were in operation at the time of this earthquake to collect the measurements due to the rupture of fault. The high accuracy and good spatial coverage of coseismic surface displacements measured by SAR satellites provide us with a valuable opportunity to determine the rupture properties of the 2016 Aketao event and to study the mode of collision interactions among Pamir, India and Eurasia plates [7,8,9].

In this study, we used the SAR observations from the latest satellites, including Sentinel-1A, Sentinel-1B and ALOS-2, to map the coseismic deformation caused by the 2016 Aketao earthquake. Based on the assumption that the medium of the earth crust is an elastic half-space homogeneous model [10], the coseismic slip model for the 2016 Aketao earthquake was derived by using the calculated Line-of-Sight (LOS) deformation. Finally, the static Coulomb failure stress (CFS) triggered by this shock was evaluated to discuss the implications and potential seismic hazards in the study area.

2. Tectonic Setting

The 25 November 2016 Mw 6.6 Aketao earthquake occurred on the Holocene active Muji fault (Figure 1), northeast of Pamir, which strikes west-northwest and dips south-southwest [7]. Pamir, situated north of the western Himalayan syntaxis [11], is surrounded by the Tajikistan depression, Tarim Basin, Tibet Plateau and Tien Shan fold system (Figure 1). Pamir was pieced together from a series of tectonic blocks during the Paleozoic and Mesozoic periods [12,13]. The Cenozoic strong shortening and thickening of Pamir is considered to be a response to the continuous collision and compression between the Indian and Eurasian plates [12,13,14], with a maximum crustal thickness of approximately 70 km [10,15], which approaches the average crustal thickness of Tibet Plateau estimated from seismic wave data [16,17]. The kinematic venation of Pamir is dominated by ~N-S shortening and orogen-parallel extrusion flow of material [8,9]. The reported geodetic measurement, using GPS, of the north-northwestward contraction rate across central Pamir was found to be nearly 25 mm/year [18].

In the past forty years, Pamir and its surroundings have experienced several moderate earthquakes, such as the 2008 Wuqia Mw 6.7 earthquake [19] and the 2016 Kyrgyzstan M 6.4 earthquake (Figure 1). The Muji fault belongs to one segment of the active Kongur Shan extensional system [7]. Satellite remotely sensed images and stream terrace studies show that the Muji fault is a dextral strike-slip structure [7,20]. Based on the cosmogenic nuclide dating method, Chevalier et al. [21] revealed that the right-lateral slip rate along the Muji fault is 4.5–11 mm/year, and the lithology of the footwall is characterized by mixed Paleozoic metamorphic [7]. As seen from Figure 1, in the Tien Shan and Pamir regions, the GPS-inferred crustal movement velocity exhibits a substantial decreasing trend from the south to north, a well-documented rapid crustal shortening phenomenon, characterized by features of the widespread active faulting and folding, as well as the recent crustal arching and high regional seismicity [22]. Hence, the tectonic activity and lithospheric deformation in Pamir is very intensive, as supported by the fact that Pamir is an earthquake-prone zone (Figure 1). The seismic moment deficit result, estimated from GPS observations from 1998 to 2013, shows a significant positive value (0.69 × 10²⁰ N·m) along the seismogenic zone of the Mw 6.6 Aketao earthquake, indicating that the Muji fault has the capability to produce an Mw 7.2 earthquake [23].

3. Coseismic Deformation Mapped by InSAR Observations

3.1. InSAR Observations and Processing

As a part of the Copernicus mission, Sentinel-1, operated by the European Space Agency (ESA), is a constellation of two satellites (Sentinel-1A and Sentinel-1B) that collects C-band SAR images by using the Terrain Observation with Progressive Scans (TOPS) mode. The constellation covers the entire earth every six days regardless of weather conditions [24]. Sentinel-1A was launched on 3 April 2014 and has been widely used to capture the surface deformation associated with earthquakes [25,26], landslides [27], and anthropogenic activities [28]. Sentinel-1B was launched on 25 April 2016 and can be used in a variety of geoscience fields. Both Sentinel-1A and Sentinel-1B observations, which are provided as free and open-access data through the Scientific Data Hub, were used in this study.

As the successor of the Advanced Land Observing Satellite (ALOS), ALOS-2 launched on 24 May 2014 carrying a Phased Array type L-band Synthetic Aperture Radar-2 (PALSAR-2) instrument to collect the L-band SAR images. Many aspects, including resolution and time lag for data provision, have been significantly improved [29]. Image products can be delivered from the ALOS-2/ALOS User Interface Gateway of the Japan Aerospace Exploration Agency (JAXA). Compared with C-band SAR sensors, the L-band ALOS-2 can preserve high coherence even in regions with mountainous topography covered by dense vegetation or glaciers as a result of the specific sensor that is equipped with a long-wavelength radar (0.229 m) [30]. In this study, we utilized three interferograms (

Table 1. Interferograms used in this study.

Satellite	Track	Reference Date	Repeat Date	Perp.B * (m)	Mean Coherence	Σ (mm)	λ (km)
Sentinel-1A	T027A	13 November 2016	7 December 2016	−100	0.51	8.19	5.3
Sentinel-1B	T107D	1 November 2016	19 December 2016	−24	0.60	9.90	9
ALOS-2	T162A	20 July 2016	7 December 2016	92	0.77	7.00	10

Perp.B * represents the perpendicular baseline between orbits; σ is the standard deviation estimated from the pixels far from the deformation region; λ is the e-folding correlation length scale of the 1D-covariance function.

) to map the responses of coseismic displacement due to slip at depth, among which the descending pass interferogram is from the Sentinel-1B satellite (T107D) and the two ascending pass interferograms are from the Sentinel-1A satellite (T027A) and ALOS-2 satellite (T162A).

The Switzerland GAMMA software [31] was employed to process all of the SAR data used in this study, including the Sentinel-1 TOPS mode and ALOS-2 PALSAR-2 Fine model SAR data. All interferograms were generated from the Single Look Complex (SLC) products. The number of multi-looks in range and azimuth, which can effectively suppress noise to some extent, were specified as 10:2 and 8:8 for Sentinel-1 and ALOS-2, respectively [26]. The strategy for coping with ALOS-2 Fine Mode interferometry is the same as for the conventional differential InSAR data. For the Sentinel-1 TOPS SAR system, its antenna rotates from backward to forward during the observation of each burst. In this context, an accuracy better than 0.001 of a pixel is required for the azimuth co-registration to avoid the phase jumps between subsequent bursts [32]. To ensure such high co-registration accuracy, a method that combines the iterative matching procedure that considers the effects of the terrain height with the procedure of spectral diversity that considers the interferometric phase of the burst overlap regions was used [33]. Once a high-quality co-registration result between the TOPS SLC data had been achieved, the differential interferograms could be created by the removal of topography effects from the interferograms with the help of a 90 m resolution digital elevation model implemented by the Shuttle Radar Topography Mission [34]. To reduce the phase noise and achieve a high signal-to-noise ratio result, an adaptive filtering method was employed for the interferograms [35]. Then, all the interferograms were unwrapped using the Minimum Cost-Flow method [36]. We also estimated and removed a quadratic ramp offset, which may be ascribed to the orbital errors, from the ALOS-2 unwrapping data using the far-field displacement. Simultaneously, an empirical approach involving the estimate of a best linear fit of the topography-correlated atmospheric delays (APS) was used to mitigate the effects of atmospheric contributions [37]. Finally, all interferograms were geolocated to the geographic coordinate at 90 m posting and the relevant unwrapped phase data were converted to line-of-sight (LOS) displacement (Figure 2).

3.2. Coseismic Deformation Mapped by InSAR

As shown in Figure 2, the entire deformation field induced by Aketao earthquake is covered by InSAR images used in this study, with strong coherence (Table 1) except in regions of snow cover. The continuous deformation stripes are clearly visible in both ascending and descending interferograms. From this figure, we can see that coseismic deformation due to the 2016 Mw 6.6 event is distributed along the Muji seismogenic fault, and more interestingly, there exist two major interferometric fringe rings, with one region near the epicenter (labeled as region A in Figure 2) and the other concentrated on the southeastern Muji fault (labeled as region B in Figure 2). The observed maximum LOS displacements are 10.6 cm, 19.2 cm and 10.5 cm for the T027A, T107D and T162A differential interferograms, respectively. The slight discrepancy between the two ascending interferograms can be attributed to the different look angles between Sentinel-1A (~33.7°) and ALOS-2 (~36.3°). The obvious deformation in region A is much broader than that in region B in the sense of spatial extension. However, interestingly, region B is characterized by more densely distributed fringes. These characteristics suggest that the earthquake may rupture a deeper depth along the dip of the fault in the region near the epicenter compared with that in region B, but allow for a shallow rupture for the B region. From focal mechanism solutions, we know that the horizontal coseismic deformation of the 2016 Aketao Mw 6.6 earthquake was essentially parallel to the fault, and the vertical deformation should also be substantial. However, InSAR measurement is less sensitive to northward motion due the sensors’ viewing geometry [38]. Therefore, it is difficult to derive a full and reasonable 3-dimensional (3D) coseismic deformation map of the Aketao earthquake.

By combining measurements from two independent look directions, the range change along the LOS direction (v_i) acquired by InSAR technique can be projected onto any desired set of two orthogonal basis vectors, with the assumption that the third component is zero [39]. Furthermore, coseismic displacements associated with a strike-slip event are essentially characterized by fault-parallel and vertical motions. Therefore, the fault-parallel (v_f) and vertical (v_z) directions of motion for a shearing structure can be reasonably expressed as follows [40]:

$$\left(\begin{array}{l}{v}_f\\ v_z\end{array}\right)={\left(\begin{array}{l}{e}_1\sin \alpha +n_1\cos \alpha u_1\\ e_2\sin \alpha +n_2\cos \alpha u_2\end{array}\right)}^{-1}\left(\begin{array}{l}{v}_1\\ v_2\end{array}\right)$$

(1)

where (e_i n_i u_i) are the satellite look vectors and α represents the average strike angle of the fault. In this study, with the assumption of 105° average fault azimuth, we decomposed the LOS deformation maps, derived from T027A and T107D SAR measurements, into the fault-parallel and vertical directions, respectively. The estimated maximum fault-parallel displacement is ~16.3 cm, with the maximum uplift and subsidence of 7.2 and 14.0 cm, respectively (Figure 3). From the coseismic displacement distribution, we can easily distinguish these two regions that share a similar deformation pattern (Figure 3b), suggesting that the rupture of the 2016 Aketao earthquake may have consisted of two sub-events that occurred in rapid succession within a few seconds.

4. Finite Slip Model and Static Stress Changes

4.1. Methods

Focal mechanism solutions derived by GCMT [4] and USGS [5] show that the seismogenic fault of the 2016 Aketao event has a steep dip angle (close to a vertical fault). Considering research results from geology [7,41], field investigations [1], and our analyses of the coseismic deformation field (Figure 2 and Figure 3), we adopted a simple 1-segment fault geometry model in which the fault dips to the southwest. In addition, based on the surface trace determined from field investigations [1], the strike angle along the trend of the fault was allowed to change. Finally, based on the above analysis, we built a simple 1-segment, strike-variable fault model composed of 880 discrete 2 km × 2 km rectangular sub-patches. We used this fault geometry model, which has an average strike of 105°, 88 km length along the strike and 40 km width along the down-slip, to invert the slip model of the 2016 Aketao event. In our study, neither the dip angle of the fault plane nor the rake angle of the slip were fixed during the inversion, the optimal dip angle was determined according to the relationship between the weighted residual sum of squares and the variable dip angles. The optimal model can be achieved when the mathematical model, summarized as follows, reaches the minimum [42]:

$$\mathsf{\varphi }(\mathrm{m})={\Vert {\displaystyle {\sum }_{\mathrm{d}}^{-1/2}(\mathrm{d}-\mathrm{G}\mathrm{m})}\Vert }^2+{\mathsf{\kappa }}^2{\Vert {\nabla }^2\mathrm{m}\Vert }^2$$

(2)

where d is the observation vector, m is the slip vector that we want to invert, G is the surface deformation response to a unit slip on fault, and κ² is a smoothing factor to balance the variables between weighted misfit data ${\Vert {\displaystyle {\sum }_{\mathrm{d}}^{-1/2}(\mathrm{d}-\mathrm{G}\mathrm{m})}\Vert }^2$ and slip roughness ${\Vert {\nabla }^2\mathrm{m}\Vert }^2$. ${\nabla }^2$ is a Laplacian smoothing operator. We smoothed the slip on fault plane with an improved Laplacian constraint [43] for this ill-conditioned problem, which is more reasonable than the classical method of Jónsson et al. [44]. In this improved method, for rectangular dislocation elements (RDEs) with four adjacent patches, the operator is the same as in Jónsson et al. [44]:

$${\nabla }^{2}m=\frac{m_{i,j-1}-2m_{i,j}+m_{i,j+1}}{{(\Delta l_r)}^2}+\frac{m_{i-1,j}-2m_{i,j}+m_{i+1,j}}{{(\Delta l_c)}^2}$$

(3)

where m_i,j represents slip on the patch located in the ith row and jth column of the fault. Δl_c and Δl_r are the distances between adjacent patches along the strike-parallel and down-dip-parallel directions, respectively. For the patches in the top and bottom rows of the fault model, the operator is modified as follows [43]:

$${\nabla }^{2}m=\frac{m_{i,j-1}-2m_{i,j}+m_{i,j+1}}{{(\Delta l_r)}^2}+\frac{m_{ii,j}-m_{i,j}}{{(\Delta l_c)}^2}$$

(4)

where ii is to be i − 1 or i + 1. For the boundary patches along the columns, the operator is expressed as follows [43]:

$${\nabla }^{2}m=\frac{m_{i,jj}-m_{i,j}}{{(\Delta l_r)}^2}+\frac{m_{i-1,j}-2m_{i,j}+m_{i+1,j}}{{(\Delta l_c)}^2}$$

(5)

where jj is to be j − 1 or j + 1.

We related the surface displacements due to slip at depth by using the elastic half-space homogeneous dislocation model [10] and assuming a Poisson ratio of 0.25. The covariance structure estimated by the 1D-covariance function [45] was adopted to represent the uncertainties of the SAR measurements (Table 1). Upon the weighted ratio of 1:1:1 for the three SAR datasets, the amount of slip on each patch was jointly inverted using the bounded variable least squares method [46]. Meanwhile, the slip distribution uncertainties were determined by using the standard deviation of the solutions inverted from the 100 perturbed datasets that were generated using a Monte Carlo method [47]. In the interest of calculation efficiency, the LOS measurements were downsampled using the quadtree technique [45] while retaining as much information as possible (Figure 4).

4.2. Finite Slip Model

We found that the fault model with a dip angle of 78° fits the LOS displacements from all three tracks very well (Figure 5b) when a 0.33 smoothing factor was chosen (Figure 5a). The best-fit inverted slip model from InSAR measurements suggests that the rupture process of this event was dominated by a dextral strike-slip motion, with a slight dip-slip component at depth (Figure 6). The maximum strike-slip and dip-slip components are ~0.85 m and ~0.26 m, respectively. Our preferred model features a heterogeneous slip pattern with a significant slip area extending 55 km along the strike and 28 km along the down-dip, and the peak slip magnitude can approach 0.85 m over the larger size asperity located near the epicenter (Figure 6). The asperity II near the southeastern Muji fault shares a similar slip maximum (0.85 m). These two geodetic inverted asperities coincide well on locations with that of substantial surface deformation regions analyzed above (Section 3.2). The estimated geodetic moment magnitude, released through the rock fracture, is 1.58 × 10¹⁹ N·m, corresponding to an event of Mw 6.74, assuming a shear modulus of 30 GPa. This estimate, however, may include a minor contribution from the post-seismic deformation caused by afterslip and/or poro-elastic relaxation, as postearthquake SAR measurements were observed with delays of several weeks (Table 1). The slip uncertainties for the main slipping patches are generally less than 4 cm and the average slip uncertainty of the optimal slip distribution is less than 2.0 cm, indicating that the inverted result is reliable.

The model results for the T027A fit the observations best (Figure 4c and Figure 7). Figure 4 and Figure 8 show the residuals from a comparison of InSAR observations with those of predictions by our best fitting model. The relatively poor fitting result for the track T107D may result, to some extent, from the longer time delays of repeated SAR measurements (Table 1), which contain parts of the post-seismic deformations. Furthermore, the APS, caused by the ionospheric effects, increases for the longer-wavelength L-band ALOS-2 SAR observations. These components of APS may exist in the SAR measurements, especially for the track T162A, and are difficult to remove with the two-pass differential interferogram method [37], and thus resulting in a relatively poor fitting. Taken as a whole, the residuals in the simulated LOS displacements are generally smaller than 4 cm, and most smaller than 2 cm, with standard deviations of 1.02, 1.87 and 1.56 cm for the track T027A, T107D and T162A LOS displacements, respectively (Figure 8). Therefore, it is reasonable to conclude that the general coseismic displacement patterns reflected by the Sentinel-1A, Sentinel-1B and ALOS-2 satellites could be accurately reproduced by the optimal slip model.

4.3. Static Coulomb Stress Changes

There is adequate evidence to support the hypothesis that coseismic slip on a fault may affect the state of neighboring faults and seismicity rates via the transfer of stress [48,49,50]. The stress state of a fault, or of those faults around it, is a crucial indicator of the likelihood that further earthquakes will occur [50]. To evaluate the effect of this event on the surrounding seismic activity, the CFS [48,49,50], a quantitative description of stress state that can be used to characterize the potential failure in rocks, was estimated, which could be helpful in further exploring the mechanisms of relevant earthquakes and may also provide new insight into earthquake predictions [50]. The CFS on optimally oriented faults at a depth of 10 km triggered by this earthquake was calculated with the optimal slip model we inverted above. The effective coefficient of friction was assumed to be 0.3, and the shear modulus was assumed to be 30 GPa. The regional tectonic principal stress, which is indispensable to the determination of the optimal failure directions, was set to 100 bars with the orientation of N25°E, which is estimated from the regional broadband long-period waveform records [51]. Figure 9 shows that the Aketao earthquake increased CFS by approximately 0.5 bar along the strike of the Muji fault, except for regions surrounding the northeastern and southeastern ends of the fault, where little elastic strain energy was released during the brittle failure of the earthquake (Figure 6). Our statistical results show that more than 65% of the aftershocks occurred in the regions where stress was estimated to have increased by at least 0.1 bar. As shown in Figure 6, the zone on the west side of the epicenter experienced many aftershocks, but the slip was quite minor. This phenomenon was observed in several cases, for example in the 2014 Kangding earthquake [52]. However, the calculated static CFS was increased by at least 0.1 bar, therefore these aftershocks may have been triggered by the static CFS. Interestingly, aftershock distributions along the Muji fault seem to have been arrested before they propagated into the southeastern end stress shadows, where no aftershocks occurred. For the surrounding substantially enhanced CFS region, where there is an absence of historical earthquakes, the disturbed stress, triggered by the 2016 Aketao earthquake, may bring them close to failure with further tectonic loading and cause an exogenous deformation event.

Geological studies indicate that the branched normal fault of the Kongur Shan extensional system has the highest metamorphic grade (650–750 °C, ~0.8 GPa), consistent with an area of the highest exhumation rates [9], indicating that vigorous geodynamic processes of crustal motion have occurred in this region. This region, as calculated above, shows positive stress increases (Figure 9). Furthermore, the GPS-inferred seismic moment deficit result shows that this region is capable of producing an Mw 7.2 earthquake [23]. Therefore, this region is an earthquake-prone area that might experience several earthquakes in the future and is thus worth paying close attention to.

5. Discussion

5.1. Implications for the Potential Barrier

Our geodetic inversion suggests that there exist two separate asperities along the fault strike, with the larger one located near the epicenter (Figure 6), this characteristic is well constrained by using both the ascending and descending SAR data and the full coverage of data for the entire rupture area, which is especially important for refinement of the shallow slip [53]. The spatial pattern of our slip model agrees, to some extent, with the result from preliminary finite fault inversion of seismic data [5]. Our slip model clearly shows that 2016 Aketao event had propagated into a shallow part of the Muji fault and ruptured the earth surface in several places (Figure 6), which matches the observations from field investigations conducted after the earthquake [1]. Notably, there is an area with a small slip between the two asperities (Figure 6), which is indicative of the complex source rupture process of this event. The spatial distribution of slip from geodetic data joins the source time function derived from seismic data [5], providing evidence for the plausible rupture pattern of this earthquake: that it consists of two sub-events that occurred in rapid succession within a few seconds, resulting in two areas that share a similar deformation pattern in their spatial responses (Figure 3b). This type of multi-event rupture pattern was observed in several cases in China, such as the 1990 Ms 6.9 Gonghe earthquake [54] and the 2010 Ms 7.1 Yushu earthquake [55]. Field investigations [1] show that the Muji fault may represent a negative flower structure in the segment between these two asperities, which is characterized by a small slip (Figure 6). Although it is relatively small in size (approximately 6 km-wide along the strike), this segment resisted being completely ruptured when the dynamic ruptures propagated into this area, indicating that this section of fault may manifest as a strong structural barrier, with a great friction force on the fault, playing a primary role in impeding the fault rupture.

5.2. Seismic Energy Source

Receiver function images from the seismic data revealed that in Pamir there exists a southerly dip, 10–15 km thick, with a dip angle that increases with the depth of the low-velocity zone that extends from 50 km depth near the base of the crust to more than 150 km in depth [11]. In mid-crustal (at a depth of ~30 km) Pamir, the existence of a weak zone is expected, which allows for the ductile shortening of crustal material [11]. Geological studies of the geometry and kinematics of major deformation zones show that the prominent ~N-S shortening, which created Pamir, is accompanied by lateral orogen-parallel flow of material into the depressions west of Pamir Plateau [8], and this may be closely related to the seismotectonics of Pamir. Thus, we speculate that the accumulated strain energy of this event may be loaded by the coupling interactions between lateral extrusion flow of the mantle material and locking of fault in the mid-upper crust [23], that is, the elastic strain energy is a result of the flowing of deep material beneath the locking depth.

5.3. Implications of Frictional Properties on the Fault Plane

As shown in Figure 10, most aftershock sequences of the 2016 Aketao earthquake occurred at depths of 3–20 km (Figure 10a), with nearly one-quarter of these sequences located at ~10 km depth (Figure 10b). Previous studies indicate that the frictional properties of fault planes are usually inhomogeneous, allowing the alternative and interwork of aseismic slip in the frictional stability zone and seismic sliding in the velocity-weakening zone from the top edge of the fault down to the deeper crust [56]. As shown in Figure 6 and Figure 10a, there is little or no slip on the shallow upper (depths at 0–3 km) and deep lower part (depth > 20 km) of the fault plane, which agrees broadly with the transition depth of the crustal faults that the upper transition depth is typically located at 3–4 km and the lower transition depth generally occurring at 15–20 km [56]. Therefore, we speculated that the friction in regions (depths of 0–3 km and depth > 20 km) surrounding the 2016 Aketao earthquake seismogenic zone (depths of 3–20 km) manifested as velocity-strengthening properties. These findings support the fact that in most places, the fault movements do not rupture the earth surface due to the velocity-strengthening frictional properties of the shallow part of the fault. In addition, the potential barrier discussed in Section 5.1 may indicate that friction on the fault planes of the seismogenic zone is also heterogeneous. However, the crustal faults have a complex dependence of friction on rate and state, which can also be affected by material properties and temperatures.

In general, earthquakes can cause an increase in stress around regions that have slipped coseismically. As a response to this type of stress disturbance, afterslip may slip aseismically on faults considering the rate-strengthening characteristics of rock friction, and this process can be accompanied by viscoelastic relaxation in rheologically weak regions within the lower crust and asthenosphere [56] as well as poroelastic rebound caused by the transition of Poisson’s ratio from undrained to drained conditions after the shock [57]. All three post-seismic mechanisms, including afterslip, viscoelastic relaxation and poroelastic rebound, can cause post-seismic deformation. Thus, post-seismic studies, by using the geodetic-measured surface displacements with high resolutions in both time and space following this earthquake, might provide us with a powerful lever for further exploration of the viscoelastic structure of the Kongur Shan extensional system and may shed light on the properties of rate- and state-dependent friction on the Muji fault.

In this study, slip inversions were performed by adopting Okada’s analytical solutions that proposed for quantitatively describing the movement of the ground surface due to the dislocation of a fault buried in the uniform elastic crustal model [10]. Hearn et al. [58] showed that the homogeneous dislocation model is adequate to some extent to calculate the static CFS in the upper crust for the exploration of stress-triggered seismicity. However, depth-dependent elasticity should be incorporated for the further study of post-seismic viscoelastic relaxation of strike-slip earthquakes, especially the large strike-slip earthquakes [58].

6. Conclusions

We extracted three entirely coseismic deformation maps from the TOPS Mode of Sentinel-1A and Sentinel-1B SAR satellites data and the Fine Mode of the ALOS-2 SAR satellite data. These coseismic deformation maps were used to model the coseismic slip distribution and determine the dip angle of the fault geometry of the 2016 Aketao earthquake in western China. The best fitting slip model suggests a relatively steep dip angle of 78°. The preferred model features a predominantly dextral strike-slip motion with a slight dip-slip component. The maximum slip is ~0.85 m, with the estimation of total seismic moment as 1.58 × 10¹⁹ N·m, corresponding to an event of Mw 6.74. Our slip model combined with the field investigations suggest that the segment of the fault between the inverted two separated asperities may have acted as a structural barrier. More than 65% of the aftershocks were triggered by the 2016 Aketao main shock, with the stress estimated to have been increased by at least 0.1 bar.