Coseismic displacement fields and the slip mechanism of the 2021 Mw 6.7 Hovsgol earthquake in Mongolia constrained by Sentinel-1 and ALOS-2 InSAR

ABSTRACT

On 11 January 2021, an Mw 6.7 earthquake, which involved a complex rupture mechanism, occurred in Lake Hovsgol, Mongolia. This study measured the earthquake’s 3-D coseismic surface displacement fields using interferometric synthetic aperture radar (InSAR) pairs of three radar look directions: Sentinel-1B from the descending path and Advanced Land Observing Satellite-2 from the ascending and descending paths. The three DInSARs showed that the maximum coseismic surface displacement appeared east of the Northern Hovsgol Fault (NHF), where the displacement components were 18 cm, 5 cm, and 33 cm in the east–west, north–south, and up–down directions, respectively (precision: 1.7 cm, 20.4 cm, and 5.2 cm). However, the 3-D displacements indicated that although the earthquake induced a combination of normal faulting with strike slip motions, the displacements in the north–south direction had very large uncertainty owing to the similar geometry of the InSARs in the descending path. Later, we performed inversions of the DInSAR-measured coseismic surface displacement fields in the line-of-sight direction, which assumes a fault plane’s uniform and distributed slips in analyzing the slip mechanism. We conducted slip distribution estimations on the ruptured fault plane in accordance with the optimal fault geometry and source parameters determined by the uniform slip model. Investigations revealed that, with a correlation of 95.3%, the simulated displacements from the best-fitted distributed slip model were consistent with the observed displacements from DInSARs. Besides, our slip distribution model showed two distinctive slip patches, which include differences in their magnitudes and directions on the fault plane. We also observed ruptured faults experiencing a predominant right-lateral strike slip with a significant dip slip, according to the slip distribution, which caused two distinct slips due to the dramatic bending of the fault strike. Then, by analyzing the Coulomb stress change, our findings proposed that the seismic risk potential of active faults in the Hovsgol Basin increased after the earthquake. Overall, the great potential of multi-track DInSAR observations in the identification of complex slip mechanisms was demonstrated.

KEYWORDS:

• Hovsgol earthquake
• InSAR
• surface displacement
• slip distribution
• Sentinel-1
• ALOS-2

1. Introduction

Earthquakes are one of the most devastating natural disasters as a result of the various causes, which alter the stress levels in the Earth’s crusts, such as the movement and collision of tectonic plates, volcanic activities, and human activities. Of the various earthquake types, tectonic earthquakes are classified as intraplate and interplate earthquakes (Scholz, Aviles, and Wesnousky 1986; Skordas et al. 1991). Intraplate earthquakes specifically occur inside a tectonic plate, which release more stress across faults and show greater earthquake intensity than interplate earthquakes. However, it is difficult to predict the epicenter, depth, and magnitude of intraplate earthquakes because of their long recurrence time, thereby causing difficulties in understanding the fault mechanisms and rupture processes of intraplate earthquakes in regions with low convergence rates (Liu and Stein 2016; Sasajima and Ito 2016; Qingyun et al. 2021).

Earthquakes that are sufficiently large and sufficiently shallow can cause surface displacement (Hirose et al. 1999; Jin et al. 2022; Jin and Fialko 2021). The source parameters of an earthquake and the rupture fault geometry can be determined by inverting the coseismic surface displacement fields using a model. This can in turn assist in understanding the fault rupture and slip mechanism. Global Positioning System (GPS) and satellite-based synthetic aperture radar (SAR) have been widely used to measure the coseismic surface displacement fields. GPS provides accurate information, with high temporal resolution, on 3-dimensional (3-D) displacement vector components in the east–west (E–W), north–south (N–S), and up–down (U–D) directions, which account for its wide use in the measurement of earthquake-induced surface deformations (Houlié, Dreger, and Kim 2014; Cheloni et al. 2016; Jiang et al. 2018; Alif et al. 2021; Crowell 2021; Hu et al. 2022). However, regional data with sparse GPS networks are insufficient in estimating the earthquake faulting mechanism. Interferometric SAR (InSAR) can differentiate phases between two or more SAR observations in the same area (Massonnet et al. 1993). Surface displacements, with precision at the sub-centimeter level, can be derived by removing topographic phases from the InSAR signals. This technique is called differential InSAR (DInSAR). Owing to its advantages, studies have adopted DInSAR as an extremely effective tool observing coseismic surface displacement (Massonnet et al. 1993; Zebker et al. 1994; Lakhote et al. 2021; Carboni et al. 2022), particularly in regions with sparse or no GPS networks. Generally, DInSAR-measured surface displacements are present in the radar look (line-of-sight, LOS) direction. Hence, at least three InSAR pairs obtained in different radar look directions should be combined in order to calculate the 3-D displacement fields and perform better-constrained inversions for the source parameters and fault geometry because ground deformation caused by an earthquake involves both horizontal and vertical displacements (Wright, Parsons, and Lu 2004; Fialko et al. 2005; Zhang, Shan, and Feng 2016; Liu et al. 2022; Fialko, Simons, and Agnew 2001).

On 11 January 2021, an Mw 6.7 intraplate earthquake occurred at Lake Hovsgol in Mongolia, which was reported to be the largest earthquake in the instrumental record history of Hovsgol graben (Liu et al. 2021, 2022). A surface crack formed by this earthquake was found on the west bank of Lake Hovsgol (Battogtokh et al. 2021), whose strike was oblique compared to that of the Hovsgol Fault of which the rupture is estimated to have caused the earthquake, which implies that a complex rupture mechanism possibly caused by the January 2021 earthquake. Several studies attempted to explain the rupture mechanism and slip distribution of a fault ruptured by the January 2021 earthquake using DInSAR (Liu et al. 2021, 2022; He, Wang, and Zhao 2022; Timoshkina et al. 2022). Liu et al. (2021) and Liu et al. (2022) determined the optimal source parameters of the earthquake using the joint inversion of teleseismic records and only one DInSAR observation. Interestingly, they both reported that the earthquake was caused by normal faulting with a right-lateral strike slip at the plane of the Hovsgol Fault. However, the pattern and magnitude of slip presented in these studies were different. Liu et al. (2021) showed two major slip patches at different depths: a predominant dip slip at deep depths (10–13 km depth, near the hypocenter) with a maximum slip of 3.0 m and a significant dextral slip at shallow depths (4–8 km depth). Liu et al. (2022) also presented a complex slip which was characterized by comparable dip slip and right-lateral slip of the ruptured fault plane; however, the maximum slip was estimated 1.2 m at a 7-km depth. The differences in the inversion results, especially for the slip distribution, could be possibly caused by insufficient DInSAR observations of the coseismic surface displacement. Notably, the coseismic displacements caused by the earthquake should have been measured in the E–W, N–S, and U–D directions since the faults in Hovsgol Basin had both normal faulting and strike slip motions. Moreover, even if teleseismic data were combined, the convergence instability of the inversion’s numerical solution would have occurred by using single DInSAR observations. DInSAR-measured coseismic surface displacements can be used more effectively than teleseismic data in resolving the fault slip mechanisms of near-field earthquakes at depths of 0–15 km (Zhang et al. 2013). Therefore, the source parameters and rupture mechanism of earthquakes occurring at shallow depths, such as the January 2021 earthquake, could also be reasonably analyzed using a multi-track DInSAR dataset. In the most recent study (He, Wang, and Zhao 2022), the slip distribution of the earthquake was analyzed according to the inversion solution estimated using the coseismic DInSAR images obtained from one descending and two ascending satellite paths. The optimal fault geometry and source parameters inverted from the multiple DInSAR observations here were comparable to those of previous studies by using the joint inversion of teleseismic data and a single DInSAR observation (Liu et al. 2021, 2022). However, unlike the previous studies, three distinctive major slip patches were estimated at different depths: a purely right-lateral slip at a 0–5 km depth and the predominant dip slips at a 5–10 km depth with a maximum slip magnitude of 1.5 m. Remarkably, this recent study provided a new perspective on understanding the earthquake’s complex slip mechanism. Still, the region of coeseismic displacement was analyzed after combining the two DInSAR observations in the ascending path; therefore, the 3-D coseismic displacement fields could not be estimated. In addition, their inversion results did not successfully simulate the DInSAR-observed displacements in the east bank of the Hovsgol Lake. Thus, DInSAR observations that are capable of estimating 3-D displacement fields in the main region of coseismic displacement and conducting inversion modeling in order to simulate the DInSAR observations accurately are required for more reliable analyses of the slip mechanism.

Based on the elucidative facts presented above, this study observed the coseismic 3-D surface displacement fields of the January 2021 Mw 6.7 Hovsgol earthquake using three DInSAR pairs which were obtained from different radar look directions. Then, we performed inversions of the DInSAR-measured displacement fields by assuming uniform and distributed slips. The optimal fault geometry and source parameters were derived, after which the earthquake’s rupture mechanism and slip distribution were analyzed. Furthermore, the potential seismic risk of the major faults in the Baikal Rift Zone was evaluated by calculating the changes in Coulomb stress after the earthquake.

2. Study area and data

2.1 Study area

The Baikal Rift Zone (Figure 1), which is one of the most active detachment rifts in Asia, comprises the Hovsgol Basin as one of its main basins (Bushingol, Darkhad, and Hovsgol) (Petit and Déverchère 2006). Specifically, the Bushingol, Darkhad, and Hovsgol basins are located where the Baikal Rift Zone expands and the India-Eurasia plate collides. Consequently, faults at the edges of these three basins exhibit a dominant normal faulting motion (Jolivet et al. 2013) with slight motions in the N–S direction (Tapponnier and Molnar 1979). However, while the Mondy Fault (MF), Tunka Fault (TF), and Bulnay Fault are located near Lake Hovsgol and show Sinistral Strike Slip motions, which suggests that the three basins expand in the E–W direction with Baikal Rift Zone expansion (Jolivet et al. 2013), the Hovsgol Fault shows both N–S motion and E–W expansion, according to the GPS observations of 1994–2002 (Calais et al. 2003). These findings indicate that the evolution of the Baikal Rift Zone is not yet clearly understood. To this end, Molnar and Tapponnier (1975) and Petit and Déverchère (2006) argued that the Baikal Rift Zone was expanding as a result of the impact of the India-Eurasia plate collision, whereas Logatchev (1993) claimed that it was not affected by plate collision and was instead expanding because of mantle convection.

Figure 1. (a) A topographic map of the Baikal Rift Zone. Black lines indicate the traces of active faults reported by the Global Earthquake Model (GEM) Global Active Fault, and white rectangles represent the imaging coverages of the Sentinel-1 and ALOS-2 InSAR pairs. (b) Enlarged topographic map around the epicenter of the January 2021 earthquake (red dotted box in (a)). The dark blue color represents Lake Hovsgol, the red stars represent the epicenter of the Mw 6.7 Hovsgol earthquake reported by different seismological agencies (Global Centroid Moment Tensor, GCMT; United States Geological Survey, USGS; Geoscience Australia, GeoAu; German Research Centre for Geosciences, GFZ; and Institut de physique du Globe de Paris, IPGP), and the yellow dots indicate the aftershock epicenters reported by USGS.

An Mw 6.7 earthquake (mainshock) occurred on 11 January 2021, at Lake Hovsgol in the Hovsgol Basin, Mongolia (Figure 1). Lake Hovsgol (1645 m above sea level), which is located at the boundary between Mongolia and Russia, is 136 km long and 20–40 km wide, which makes it the second-largest lake in Asia (Kouraev et al. 2016). It is typically frozen from November to June, and the mean thickness of the lake ice during the winter is ~1.5 m (Zhang et al. 2021). The lake ice exhibited an extensive displacement anomaly of approximately 10 m in the E–W direction and 5 m in the N–S direction due to the earthquake, which was supposed to be strongly related to the slip direction of the earthquake (He and Wen 2022). The Hovsgol Basin experiences cold and dry weather with marginal precipitations (Etzelmüller et al. 2006).

The source parameters for the mainshock of the January 2021 Hovsgol earthquake as reported by several organizations and studies are listed in Table 1. The focal depth and magnitude of the strike slip component remain inconsistent, although similar epicenter locations and consistent normal faulting have been reported, which suggests a complex faulting mechanism of the earthquake (Liu et al. 2022). After the mainshock, more than 52 aftershocks (> Mw 4.0) were recorded, and their epicenters are depicted as yellow dots in Figure 1(b).

Table 1. Source parameters of the mainshock of the January 2021 Hovsgol earthquake as reported by different organizations and studies.

	Strike (°)	Dip (°)	Rake (°)	Depth (km)	Mw
GCMT	354/236	43/66	−143/−53	13.9	6.8
USGS	16/219	32/60	−110/−78	11.5	6.74
IPGP	358/237	46/62	−139/−52	13.0	6.84
GeoAu	53/172	43/65	322/−126	18.0	6.7
GFZ	4/226	47/51	−121/−60	18	6.7
Liu et al. (2021)	345	42	−140	11	6.75
Liu et al. (2022)	340.4	53.9	−146.4	8.9	6.75
He, Wang, and Zhao (2022)	340	53	−116	-	6.71
This study	343	59	−143	-	6.78

2.2 DInSAR dataset

In this study, the Sentinel-1B C-band (5.4 GHz center frequency) and Advanced Land Observing Satellite-2 (ALOS-2) Phased Array type L-band Synthetic Aperture Radar (PALSAR-2) L-band (1.2 GHz center frequency) InSAR pairs were used to measuring the coseismic 3-D displacement fields of the earthquake’s mainshock (Table 2). Sentinel-1B is the second satellite unit of the Sentinel-1 constellation operated by the European Space Agency (ESA) comprising Sentinel-1A and −1B, launched on 3 April 2014 and 25 April 2016, respectively. Both satellites have a 12-day revisit cycle and share similar orbital planes with a phasing difference of 180°, which allows for observations at 6-day intervals (Torres et al. 2012). Sentinel-1 is operated in four acquisition modes: Stripmap, interferometric wide (IW) swath, extra-wide swath, and wave. Considering that the Sentinel-1 SAR data are distributed free of charge via the Copernicus Open Access Hub (https://scihub.copernicus.eu/) by the ESA, it has become a valuable tool in monitoring surface deformation by geological hazards (Bui et al. 2021; Festa et al. 2022; Cigna and Tapete 2021; Lyu et al. 2020; Raspini et al. 2018). In this study, we downloaded an InSAR pair of Sentinel-1B at the descending path of the IW mode (250 km swath width and 20 m nominal resolution) and VV polarization (Table 2) through the Copernicus Open Access Hub. We then concatenated Sentinel-1B SAR images of two consecutive frames acquired on the same dates to cover the entire coseismic displacement area of the earthquake (Figure 1). Conversely, Sentinel-1A was not used because it has not observed the Hovsgol Basin since March 2017 according to ESA’s Sentinel-1 observation plans.

Table 2. Sentinel-1B and ALOS-2 InSAR pairs used to measure the coseismic surface displacements of the mainshock of the January 2021 earthquake.

	Sentinel-1B	ALOS-2
Dates of InSAR pair (yyyy/mm/dd)	2021/01/07 2021/01/19	2020/06/21 2021/06/20	2020/11/02 2021/11/01
Orbit direction	Descending	Ascending	Descending
Acquisition mode	Interferometric Wide	Fine Beam	Wide Beam
Polarization	VV	HH	HH
Perpendicular baseline (m)	−3.40	−7.43	84.86
Temporal Baseline (day)	12	364	364
Incidence angle (°)	39.4	40.6	29.2
Heading (°)	195.3	−9.8	194.8

ALOS-2 PALSAR-2, launched on 24 May 2014, as the successor of ALOS PALSAR, is operated by the Japan Aerospace Exploration Agency (JAXA). ALOS-2 PALSAR-2 operates in the Spotlight, Stripmap, and ScanSAR modes with a 14-day revisit cycle (Rosenqvist et al. 2014). The Stripmap mode is subdivided into ultrafine, highly sensitive, and fine modes according to the ground resolution. We obtained an ALOS-2 InSAR pair of the Stripmap fine mode (70 km swath width and 10 m nominal resolution) from the ascending path, and that of the ScanSAR mode (350 km swath width and 100 m nominal resolution) from the descending path, under our collaboration with the second Earth Observation Research Announcement (EO-RA2) of JAXA. All PALSAR-2 images were acquired with HH polarization. Notably, the ALOS-2 ScanSAR image comprised five sub-swaths, of which the first sub-swath that captured the Hovsgol Basin was used for InSAR processing. The long temporal baseline of ALOS-2 InSAR pairs (364-day) can cause decorrelation and challenge the measurement of the surface displacement using the InSAR technique. Nevertheless, we expected that high interferometric coherence could be maintained owing to dry weather conditions in the study area and the long wavelength of the L-band.

All SAR images were delivered in a single-look complex (SLC) format in promoting the application of the InSAR technique. We calculated the 3-D displacement fields by combining the DInSAR results of each pair since the three InSAR pairs exhibited different radar look directions and incidence angles. A Copernicus GLO-30 digital elevation model (DEM) with a spatial resolution of 30 m was used for the co-registration of SAR images, removal of topographic phases from InSAR signals, and terrain corrections of the acquired images.

3. Methodology

3.1 Measurement of the coseismic surface displacements by DInSAR

Using the commercial software GAMMA, developed by GAMMA Remote Sensing in Switzerland, we applied DInSAR technique in deriving the coseismic 3-D displacement fields from the Sentinel-1B and ALOS-2 PALSAR-2 data. A multi-look process of SLC images was performed using a window size of 10 (range) × 2 (azimuth) for Sentinel-1B and 3 × 7 and 3 × 20 for ALOS-2 Fine and ScanSAR images, respectively, to reduce speckle noise in the SAR images. Each InSAR pair was accurately co-registered using the GLO-30 DEM, and an interferogram was generated. To generate the DInSAR images, topographic phases had been removed from the interferograms using the GLO-30 DEM. Further, an adaptive spectral filter (Goldstein and Werner 1998) was applied to the DInSAR images to reduce phase noise. Subsequently, phase unwrapping was performed, in which pixels with a < 0.2 coherence value were masked and not used. Longas presumed to be due tto orbital and atmospheric-phase delays. These were simulated by 2-D quadratic polynomial fitting of the unwrapped DInSAR phases in the region without coseismic surface displacement and were eliminated from the DInSAR signals.

The 3-D displacement fields can theoretically be obtained by exploiting multiple DInSAR observations from at least three different radar look directions since the DInSAR-measured surface displacement is usually in the LOS direction (Hu et al. 2014). Therefore, we exploited Sentinel-1B and ALOS-2 DInSAR observations to produce the coseismic 3-D displacement fields of the earthquake. The coseismic 3-D displacements were determined as follows (Hu et al. 2014; Fuhrmann and Garthwaite 2019):

(1) $\left(\begin{array}{c}{d}_E\\ d_N\\ d_U\end{array}\right)=\Gamma \cdot \left(\begin{array}{c}{d}_{S\mathrm{\_}dsc}\\ d_{A\mathrm{\_}dsc}\\ d_{A\mathrm{\_}asc}\end{array}\right)$(1)

with

$\Gamma =\left(\begin{array}{c}{\Gamma }_{S_{dsc}}\\ {\Gamma }_{A_{dsc}}\\ {\Gamma }_{A_{asc}}\end{array}\right)={\left(\begin{array}{c}-sin{\theta }_{S_{dsc}}cos{\alpha }_{S_{dsc}}sin{\theta }_{S_{dsc}}sin{\alpha }_{S_{dsc}}cos{\theta }_{S_{dsc}}\\ -sin{\theta }_{A_{dsc}}cos{\alpha }_{A_{dsc}}sin{\theta }_{A_{dsc}}sin{\alpha }_{A_{dsc}}cos{\theta }_{A_{dsc}}\\ -sin{\theta }_{A_{asc}}cos{\alpha }_{A_{asc}}sin{\theta }_{A_{asc}}sin{\alpha }_{A_{asc}}cos{\theta }_{A_{asc}}\end{array}\right)}^{-1}$

where $d_E$, $d_N$, and $d_U$ are the displacement vector components in the E–W, N–S, and U–D directions, respectively, $\theta$ and $\alpha$ are the radar incidence angle and heading of the satellite, respectively, and $d_{S\mathrm{\_}dsc}$, $d_{A\mathrm{\_}dsc}$, and $d_{A\mathrm{\_}asc}$ indicate the LOS displacements measured at the Sentinel-1B at descending, ALOS-2 at descending, and ALOS-2 at ascending paths, respectively.

We determined the errors of the estimated 3-D displacements by InSAR pairs’ incidence and heading angles and the precision of the displacements in the different LOS directions (He et al. 2019; Hu et al. 2014). The precision of the LOS displacement was estimated by calculating the mean and standard deviation values in the region with no coseismic surface displacement (He et al. 2019; Hu et al. 2014). The uncertainties in the E–W, N–S, and U–D displacement components (${\sigma }_E,{\sigma }_N.\mathrm{a}\mathrm{n}\mathrm{d}{\sigma }_U$) were determined as follows (He et al. 2019; Hu et al. 2014):

(2) $\left(\begin{array}{c}{\sigma }_E^2\\ {\sigma }_N^2\\ {\sigma }_U^2\end{array}\right)=\left(\begin{array}{c}{\Gamma }_{S_{dsc}}{ }^{.}\mathrm{\Sigma }{ }^{.}{{\Gamma }_{S_{dsc}}}^T\\ {\Gamma }_{A_{dsc}}{ }^{.}\mathrm{\Sigma }{ }^{.}{{\Gamma }_{A_{dsc}}}^T\\ {\Gamma }_{A_{asc}}{ }^{.}\mathrm{\Sigma }{ }^{.}{{\Gamma }_{A_{asc}}}^T\end{array}\right),{\sum }^{=}\left(\begin{array}{ccc}{\sigma }_{S_{dsc}}^2& {\sigma }_{S_{dsc}{A}_{dsc}}& {\sigma }_{S_{dsc}{A}_{asc}}\\ {\sigma }_{S_{dsc}{A}_{dsc}}& {\sigma }_{A_{dsc}}^2& {\sigma }_{A_{dsc}{A}_{asc}}\\ {\sigma }_{S_{dsc}{A}_{asc}}& {\sigma }_{A_{dsc}{A}_{asc}}& {\sigma }_{A_{asc}}^2\end{array}\right)$(2)

where ${\sigma }_{S_{dsc}}$, ${\sigma }_{A_{dsc}}$, and ${\sigma }_{A_{asc}}$ are the standard deviation of the LOS displacements in the non-deforming area measured at the Sentinel-1B at descending, ALOS-2 at descending, and ALOS-2 at ascending paths, respectively. In matrix $\mathrm{\Sigma }$, the off-diagonal components are the covariance between different LOS measurements, which can be regarded as zero because they are generally mutually independent (Hu et al. 2014).

3.2 Inversions of the coseismic DInSAR displacement fields

To estimate the fault geometry and source parameters of the earthquake, we performed inversions of the DInSAR-measured LOS displacement fields using open-source software called the geodetic Bayesian Inversion Software (GBIS, available for download at http://comet.nerc.ac.uk/gbis), which was provided by the Center for the Observation and Modeling of Earthquakes, Volcanoes, and Tectonics (Bagnardi and Hooper 2018). First, the inversion was conducted by simulating the DInSAR-measured displacements, assuming a rectangular dislocation with a uniform slip in the elastic half space (Okada 1985). Then, we determined the optimal fault geometry and source parameters through inversion when the simulated displacements showed the best-fit to the DInSAR-measured displacements.

When a complex slip mechanism is expected, like this earthquake, it is difficult to understand the homogeneous slip mechanism provided by the uniform slip model, although the uniform slip model can provide reasonable fault geometry and source parameter solutions, thereby requiring slip distribution modeling. Therefore, we extended the fault plane in the strike and dip directions, according to the optimal fault geometry, which was obtained from the uniform slip modeling, in order to estimate the slip distribution. Then, we discretized it into rectangle patches of the same size. We used the steepest descent method (SDM) and Laplacian smoothing derived from a finite difference approximation of the Laplace operator to solve the slip distribution on the fault plane and regularize the solution (Wang, Diao, and Hoechner 2013; Wang et al. 2012). A smoothing factor was applied to the slip distribution inversion to reduce the instability and uncertainty of the inversion and ensure the continuity of the slip between patches. An appropriate smoothing factor was determined by analyzing the trade-off curve between misfit and model roughness for a series of smoothing factors.

3.3 Evaluation of potential seismic risk

Since earthquakes occur most probably in fault zones where the static Coulomb failure stress change (∆CFS) value exceeds 0.01 MPa (King, Stein, and Lin 1994), and aftershocks can occur in faults having high ∆CFS value, thereby releasing the stresses accumulated in the faults, ∆CFS has been widely used in indicating the potential seismic risks of faults (Li et al. 2021; Parsons et al. 2006; Utkucu et al. 2013; Wan and Shen 2010). We evaluated the potential seismic risk of major faults around the Hovsgol Basin after the earthquake by analyzing the ∆CFS. ∆CFS was calculated by using the optimal source parameters of the mainshock of the 2021 Hovsgol earthquake. The ruptured fault predicted by the inversion was set as a source fault, while the other faults were considered as receiver faults. The fault strikes required for calculating ∆CFS were measured by dividing these strikes into several segments at points where the fault traces were bent since the fault traces were not linear. The dip and rake of the MF estimated by Delouis et al. (2002) and those of the other receiver faults registered in the Global Earthquake Model (GEM) Global Active Fault database were used in calculating ∆CFS by using Coulomb 3.3 software (Toda et al. 2005). For this calculation, we assumed an elastic half-space with a shear modulus of 32 GPa and a Poisson’s ratio of 0.25 (Okada 1992). The friction coefficient was set to 0.4 for all depths, which was considered appropriate for most faults (King, Stein, and Lin 1994; Harris 2000). To examine whether the calculated ∆CFS was reliable, its spatial distribution was compared to the epicenter > Mw 4 aftershocks, reported by the USGS, for depth ranges of 1–10 km, 10–15 km, and 15–20 km.

4. Results and discussion

4.1 Coseismic surface displacement measured by multi-track DInSAR

The DInSAR images generated from the Sentinel-1 and ALOS-2 InSAR pairs are shown in Figure 2. Although ALOS-2 InSAR pairs have a long temporal baseline (364-day), high coherence was maintained because of the dry weather conditions and long wavelengths of the L-band SAR. The maximum LOS displacement was observed northwest of Lake Hovsgol and east of Northern Hovsgol Fault (NHF) at 22 cm, 33 cm, and 28 cm from the DInSAR of Sentinel-1B descending, ALOS-2 ascending, and ALOS-2 descending geometries, respectively (Figure 2). Furthermore, the DInSAR fringes in the descending paths showed that the east bank of the Lake Hovsgol moved slightly (~8 cm) in the LOS direction, that is, in the opposite direction toward the western region of the NHF. This indicates that the ground moved horizontally as a result of the earthquake. We also observed that the fringes northwest of Lake Hovsgol, with a high fringe density, differed in the three DInSAR images, which suggests that the coseismic rupture mechanism was a complex process, and multi-track DInSAR observations were essential in analyzing the coseismic deformation and slip mechanism of the earthquake. In all DInSAR images, a low fringe rate over a wide area appeared west of the NHF, whereas high fringe rates indicating intense surface displacement in different directions were observed east of the NHF. Both normal faulting motions and strike slips likely caused this disparity.

Figure 2. Differential interferogram for the mainshock of the January 2021 earthquake generated from the (a) Sentinel-1B descending interferogram, (b) ALOS-2 ascending interferogram, and (c) ALOS-2 descending interferogram. The black lines represent the traces of the active faults reported by GEM.

Table 3 shows the precision of the DInSAR-measured displacement determined as the standard deviation of the displacements in the region assumed to be free from coseismic deformations. We observed that the errors of the displacements in the LOS direction of the three DInSAR pairs were small: 0.63 cm for Sentinel-1 descending, 1.11 cm for ALOS-2 descending, and 1.20 cm for ALOS-2 ascending geometries. However, ALOS-2 DInSARs showed displacement measurement errors about twice as large as Sentinel-1 DInSAR, which was possibly caused by the decorrelation as a result of the long temporal baselines and atmospheric-phase perturbations. Nevertheless, these errors were small and could be ignored when analyzing the coseismic surface displacements in the LOS directions.

Table 3. Precisions of DInSAR-measured surface displacements in the LOS direction.

Satellite	Orbit Direction	Statistics of LOS displacement error
		Mean (cm)	Standard deviation (cm)
Sentinel-1	Descending	0.63	0.63
ALOS-2	Descending	0.38	1.11
ALOS-2	Ascending	−0.40	1.20

Figure 3 shows coseismic 3-D displacement components (N–S, E–W, and U–D), which were calculated for the overlapped area in all three DInSARs and the region east of Lake Hovsgol was not covered. In the region east of the NHF, we observed horizontal displacements of approximately 18 cm eastward (Figure 3b) and 5 cm southward (Figure 3a), with a subsidence of up to 33 cm (Figure 3c). Contrastingly, the region west of the NHF, where a low fringe rate appeared over a wide area, shifted horizontally to the west by ~14 cm and was uplifted by ~3 cm. Such a surface displacement pattern strongly suggested that a mixed mechanism of normal faulting and strike slip faulting caused the earthquake.

Figure 3. Surface displacements in (a) N–S, (b) E–W, and (c) U–D directions calculated from the Sentinel-1B and ALOS-2 DInSAR observations. The black lines represent the traces of the active faults reported by GEM.

The N–S displacement component was smaller than the E–W and U–D components (Figure 3). Previous studies estimated that a significant right-lateral strike slip occurred at the shallow depths (0–5 km) near locations of maximum coseismic surface displacements through distributed slip modeling by inverting the DInSAR-measured coseismic surface displacement (Liu et al. 2021; He, Wang, and Zhao 2022). The strike slip at shallow depths could cause large surface displacements in the N–S direction, considering the strike of the NHF. The look angles of the Sentinel-1 and ALOS-2 SAR images at the descending path differed by ~ 10°, but they had almost the same heading angle. This could cause a large uncertainty when calculating the N–S surface displacements, even though an InSAR pair from the ascending path was used for the calculations (He et al. 2019). With the method for evaluating errors in multi-track DInSAR-based 3-D surface displacement described in Hu et al. (2014) and He et al. (2019), although, the E–W and U–D displacement components calculated in this study had uncertainties of 1.7 cm and 5.2 cm, respectively, the N–S displacement component had a very large uncertainty of 20.4 cm, which makes it possible that a larger surface displacement occurred than the calculated displacement as shown in Figure 3(a).

The surface displacements of the January 2021 earthquake observed from the three DInSAR observations and the calculated 3-D displacement fields could include some uncertainties due to several factors. First, all DInSAR epochs comprise the mainshock and several aftershocks of >Mw 5. Therefore, the surface displacements measured by DInSARs could include the ground motions caused by aftershocks. Second, Sentinel-1 and ALOS-2 DInSAR were affected by orbital, tropospheric, and ionospheric phase delays (Gomba, González, and Zan 2017; Liang et al. 2019; Xu et al. 2021; Zhao et al. 2021). We analyzed a DInSAR image generated from the Sentinel-1 SAR images obtained on 25 April 2021 and 7 May 2021, which captured an Mw 5.6 aftershock on 3 May 2021 (the largest aftershock), to examine whether an Mw 5 earthquake in the study area caused surface displacements (Figure 4). We observed no DInSAR signals that could be interpreted as the coseismic surface displacements of the aftershock. As C-band DInSAR is more sensitive to smaller surface displacements than L-band DInSAR, the surface displacements caused by the aftershocks were hardly reflected in both Sentinel-1 and ALOS-2 DInSAR observations. Furthermore, the effect of atmospheric-phase delay was rarely observed in the Sentinel-1 interferogram for the mainshock (Figure 2a). Conversely, ionospheric and orbital contributions to the InSAR signals were observed in the ALOS-2 interferograms, but they were effectively removed through phase trend modeling through the 2-D quadratic polynomial fitting of the unwrapped DInSAR phases (Figure 2(b,c). Therefore, we could assume that the effects of atmospheric-phase distortion on the calculation of 3-D coseismic surface displacement fields were negligible.

Figure 4. DInSAR image generated using the Sentinel-1 InSAR pair obtained on 25 April 2021 and 7 May 2021, which captured an Mw 5.6 aftershock on 3 May 2021. The black lines represent the active faults reported by GEM and the red star indicates the aftershock epicenter reported by the USGS.

4.2 Fault geometry and source parameters inverted from the coseismic DInSAR data

A complicated DInSAR fringe pattern was observed around the NHF (Figure 2), which suggests the possibility of a complex slip mechanism. We performed inversions for the DInSAR-measured coseismic displacement fields, assuming that a single rectangular fault plane existed with a uniform slip. The optimal parameters for the best-fit fault plane are listed in Table 4, which include the final maximum of posterior probability solutions at a 95% confidence interval. The optimal parameters present the ruptured plane experiencing normal faulting motion with a significant right-lateral strike slip, as reported in Liu et al. (2022) and He, Wang, and Zhao (2022). The fault geometry determined assuming the uniform rectangular dislocation was generally consistent with different seismological agencies and previous studies (Table 1). The displacement fields are measured using the multi-track DInSAR data (observation) and those simulated using the inversion result (simulation), and the residuals between the measured and simulated displacements are shown in Figure 5. The simulated displacement fields agreed well with the observed displacement fields, which show small root mean square error (RMSE) and normalized RMSE (NRMSE) values (Table 5). The simulated displacements also fitted the observations in the region east of Lake Hovsgol reasonably well. However, in the region northwest of Lake Hovsgol, where extremely complex and dense DInSAR fringes appeared, the simulated displacements largely deviated from the observed displacements, with residuals up to 12 cm. Moreover, the moment magnitude calculated by the inversion result was Mw 6.67, which was smaller than that reported by different organizations and studies (Table 1). These findings show that the uniform dislocation model could not fully explain the slip mechanism of the earthquake.

Figure 5. DInSAR-observed and simulated surface displacements, assuming a uniform slip. Line-of-sight (LOS) displacements were measured from the (a) Sentinel-1B descending path, (d) ALOS-2 ascending path, and (g) ALOS-2 descending path. Simulated LOS displacements for the (b) Sentinel-1, (e) ALOS-2 ascending path, and (h) ALOS-2 descending path are also presented, including the residuals between the measured and simulated displacements for the (c) Sentinel-1, (f) ALOS-2 ascending path, and (i) ALOS-2 descending path. Black lines represent the traces of the active faults reported by GEM.

Table 4. Optimal fault geometry and source parameters derived from the inversions performed by assuming a uniform slip.

	Optimal	Mean	Median	2.5%	97.5%
Length (m)	33,050.88	33,724.4	33822.1	31,164.49	36,635.84
Width (m)	13,128.05	12,791.7	12781.3	12,100.79	13,550.75
Depth (m)	14,147.12	13,926.1	14233.71	13,374.67	14,534.59
Strike (°)	343.4	343.421	343.43	343	343.84
Dip (°)	59.04	59.4911	59.48	57.67	61.37
Strike-slip (m)	−0.69	−0.70	−0.70	−0.76	−0.65
Dip-slip (m)	−0.41	−0.41	−0.41	−0.43	−0.40
Fault X center (m)	−6262.24	−6412.04	−6380.17	−4959.26	−6856.09
Fault Y center (m)	7697.42	7293.5	7248.4	5846.99	8558.79

Table 5. Statistics of the residuals between the DInSAR-measured and modeled surface displacements inverted by assuming the uniform and distributed slips.

Satellite	Orbit direction	Uniform slip			Distributed slip
		Mean error (cm)	RMSE (cm)	NRMSE (%)	Mean error (cm)	RMSE (cm)	NRMSE (%)
Sentinel-1	Descending	−0.0003	1.20	4.23	−0.0043	1.24	4.37
ALOS-2	Descending	−0.0019	1.80	4.63	−0.0026	1.58	4.06
ALOS-2	Ascending	−0.0209	2.34	5.08	−0.0123	2.19	4.75

4.3 Slip distribution inverted from the coseismic DInSAR data

We subsequently extended the ruptured fault plane determined from the uniform dislocation model (Table 4) to 50 km and 30 km in the strike and dip directions, respectively, then discretized them into 2-km × 2-km patches. Afterward, the slip on all patches was estimated by applying a series of smoothing factors. Figure 6 shows the trade-off curve between misfit and model roughness, from which we determined an appropriate smoothing factor as 0.02. The DInSAR-observed coseismic surface displacements, the simulated displacements, and the residual between the DInSAR and simulations are shown in Figure 7. The correlation between the observation and simulation was 95.3%, which indicates that the simulated displacements fitted well with the observed displacements. In particular, the simulation based on the slip distribution was much improved from the uniform slip modeling in the region northwest of Lake Hovsgol, showing a residual of less than 5 cm.

Figure 6. A trade-off curve between misfit and model roughness. The optimal smoothing factor that balances model misfit and roughness was determined as 0.02.

Figure 7. DInSAR-observed and simulated surface displacements, assuming the distributed slip. Line-of-sight (LOS) displacements were measured from the (a) Sentinel-1B descending path, (d) ALOS-2 ascending path, and (g) ALOS-2 descending path. Simulated LOS displacements for the (b) Sentinel-1, (e) ALOS-2 ascending path, and (h) ALOS-2 descending path are also presented, including the residuals between the measured and simulated displacements for the (c) Sentinel-1, (f) ALOS-2 ascending path, and (i) ALOS-2 descending path. Black lines represent the traces of the active faults reported by GEM.

Figure 8 shows the estimated slip distribution. We observed that two major slips occurred at depths of 4 km to 11 km, with the maximum slip being 1.63 m at a rake of−144.3° and an 8 km depth. The magnitude of the earthquake calculated through the distributed slip modeling was Mw 6.78, which is slightly larger than that reported by the USGS and other studies (Liu et al. 2021, 2022; He, Wang, and Zhao 2022) (Mw 6.7–6.75, Table 1), but comparable with the GCMT catalog (Mw 6.8). The rake determined in this study was larger than that determined from He, Wang, and Zhao (2022), who inverted the coseismic surface displacements observed in the ascending and descending DInSARs. The rake of the ruptured fault plane reported by the GCMT is known to be similar to that inverted by accurate DInSAR observations (Weston, Ferreira, and Funning 2011). Therefore, since the rake determined in this study was similar to the GCMT catalog (−143°), it is considered reliable.

Figure 8. Slip distribution of the fault plane from the preferred model. (a) Total, (b) strike, and (c) dip slips. The red box in (a) indicates the location where the strike of the NHF is dramatically bent. Black lines represent the traces of the active faults reported by GEM.

Our slip distribution model shows that the slip magnitudes and directions of the segmented patches differed on the fault plane, contrasting the slip distribution estimated by Liu et al. (2022), which suggested that the occurrence of slip by the earthquake occurred at similar directions on the whole fault plane. The Hovsgol Basin has a complex stress structure (Liu et al. 2021), which could be verified from the coseismic surface displacement patterns (Figures 2 and 3). Therefore, it is somewhat unrealistic that the slips caused by the earthquake were in similar directions on the whole fault plane. Liu et al. (2021) suggested a slip distribution in which the strike slip predominated in the shallow depths (0–10 km) of the northern NHF and a significant dip slip occurred at the deep depths (10–15 km) of the southern NHF. However, the correlation between the DInSAR observations and the simulations was 80%, which is lower than our results. Our slip distribution pattern is quite similar to that of Liu et al. (2021). However, it showed dominant strike slips (~1.3 m) at the 0–25 km depths with significant dip slips (~1.0 m) below 4 km depths of both the northern and southern NHF. This result suggests that the right-lateral strike slip occurred more significantly at a wider area than the normal faulting motion, which makes it possible for the surface displacement in the N–S direction to be greater than that in the E–W direction. The surface displacements in the N–S direction shown in Figure 3(a) present that the upper part of the NHF (east of the fault trace) moved southward, indicating dextral movement, and its magnitudes were smaller than those in the E–W and U–D directions. However, the displacement in the N–S direction calculated in this study had a large uncertainty (20.4 cm). The slip distribution showed that strike slip occurred predominantly at shallow depths, which suggests that the actual surface displacement in the N– direction was larger than that calculated by the multi-track DInSAR observations.

The slip distribution shows that two distinctive strong strike slips and dip slips occurred at the plane of the northern and southern NHF. The surface corresponding to the depths where two strong dip slips occurred was consistent with a large eastward surface displacement (Figure 3b). The locations of the two major total slips are distinguished based on where the strike of the NHF is abruptly bent. The strike of the NHF ranges from 340° to 10° (Liu et al. 2022). Taking into consideration that irregular fault shapes (strike bending, fault-line breaks, etc.) cause the slip termination and restart (King and Nábělek 1985; Elliott et al. 2015; Yang, Liu, and Lin 2012; Yao and Yang 2022), the slip generated at the hypocenter of the earthquake propagated and relieved where the strike of the NHF is dramatically bent (the red rectangle in Figure 8a). As a result, another slip might be caused at the shallow depths of the northern NHF.

Some coseismic surface displacements were observed at the east of Lake Hovsgol. In previous studies, Liu et al. (2022) presented that the Hovsgol Fault extended the Hovsgol Basin, and He, Wang, and Zhao (2022) reported that the 2-D (E–W and U–D) displacements from multiple DInSAR observations showed little subsidence and some eastern displacements in the east of the lake. Similarly, the slip distribution estimated in this study had a significant dip slip at depths of 8–20 km near the hypocenter, which suggests the possibility of eastern displacement with slight subsidence in the east of the lake.

4.4 Static Coulomb stress changes and the seismic hazard potential

Figure 9 shows the spatial distribution of ∆CFS at different depths calculated from the source parameters. The aftershock epicenters reported by the USGS (white dots in Figure 9) at all depth ranges were located in the region showing positive static ∆CFS (stress-increased region), except four aftershocks at a depth of 0–10 km and an aftershock at a depth of 15–20 km. We also observed that in the region close to the mainshock depth (10–15 km), the NHF released up to 3 MPa stress.

Figure 9. ∆CFS values at depths of (a) 0–10 km, (b) 10–15 km, and (c) 15–20 km caused by the January 2021 earthquake calculated using the optimal source parameters of Model 3 given in Table 4. The white dots indicate the aftershock epicenter reported by USGS, and the black lines represent the active fault traces.

Figure 10 shows the ∆CFS distribution of the faults around the mainshock epicenter. The ∆CFS distribution was calculated every 5 km for depths of 0–20 km and was segmented where the fault strikes changed. Investigations revealed that the Southern Hovsgol Fault (SHF) was stressed (positive ∆CFS) at all depths by >0.01 MPa. Stress increased at shallow depths (0–10 km) of the northernmost segments of the NHF. In the Yammatinskiy Fault (YF) and MF, stress increased at all depths only in segments close to the NHF, and it was released in the other fault segments. While stress relief was observed at all depths in the Northern Darkhad Fault (NDF), except at depths of 15–20 km, stress was accumulated in the Southern Darkhad Fault (SDF) at all depths.

Figure 10. ∆CFS values of the faults around the epicenter of the January 2021 earthquake.

We inferred that aftershocks could occur at depths of 0–10 km of the northernmost segment of the NHF while evaluating the seismic risk associated with a critical ∆CFS value of 0.01 MPa caused by the mainshock (Figure 10). Aftershocks of > Mw 4 occurred 46 times within 10 km of depth, and most of the accumulated stresses were released in the region. The northernmost region of the SHF, the easternmost region of the YF, and the westernmost region of the MF showed ∆CFS values of ~0.01 MPa. As these faults are assumed to be active faults (Arzhannikova et al. 2011) and high-magnitude earthquakes have not occurred for decades, the seismic risk potential of the fault zones can be interpreted as extremely high. Besides, at the SHF, while earthquakes of < Mw 4 have been occurring continuously since the 1900s, aftershocks of the January 2021 earthquake have rarely occurred. In the northern SHF segment, the ∆CFS values of ~0.05 MPa were calculated at deep depths. These findings show that the SHF has a high seismic risk potential.

5. Conclusions

The coseismic surface displacements caused by the January 2021 Hovsgol Mw 6.7 earthquake were measured using the Sentinel-1 and ALOS-2 DInSAR pairs at three different radar look directions. The multi-track DInSAR observations provided 3-D coseismic surface displacement fields, which showed that while the east of the NHF moved up to 18 cm eastward, 5 cm southward, and 33 cm subsidence, respectively, the west of the NHF shifted to the west by ~14 cm and was uplifted by ~3 cm. We estimated the slip mechanism by using inversions of the DInSAR-measured coseismic displacement fields based on the uniform slip and the distributed slip of a fault plane. The results indicated that the surface displacement fields simulated by the distributed slip model were consistent with the DInSAR-measured displacement fields and improved from the uniform slip model, especially in the main coseismic deformation region. In addition, the optimal source parameters and slip distribution derived from the best-fit model showed that the ruptured fault experienced a predominant right-lateral strike slip (~1.3 m) with a significant dip slip (~1.0 m). Based on the slip distribution, we infer that the slip was propagated from the hypocenter, and it was relieved where the NHF strike is dramatically bent, thereby triggering slips at the shallow depths of the northern NHF. Finally, the Coulomb stress distribution changes suggested that active faults in the Hovsgol Basin could have high seismic risk potential.

Overall, our results showed that the multi-track InSAR observations incorporating Sentinel-1 and ALOS-2 data could provide rich and precise information on coseismic deformations to ensure the reliability of the inversion solutions, which enables the identification of complex rupture mechanisms. In the future, we will study the seismological characteristics of the active faults in this region comprehensively by identifying and analyzing variations in micro-displacements through time-series InSAR techniques.

Acknowledgments

This study was supported by the National Research Foundation of Korea (NRF-2021R1C1C1009621 and No. 2019R1A6A1A03033167) and Ministry of the Interior and Safety as Earthquake Disaster Prevention Human resource development Project. This study was also supported by the Japanese Aerospace Agency (JAXA), under the 2^nd Earth Observation-Research Announcement (EO-RA2) collaboration with PI number of PER2A2N144 to obtain PALSAR-2 data. The authors would like to thank the European Space Agency (ESA) for providing Sentinel-1 SAR data and GLO-30 DEM and Centre for the Observation and Modelling of Earthquakes, Volcanoes and Tectonics (COMET) for sharing Geodetic Bayesian Inversion Software (GBIS).

Disclosure statement

No potential conflict of interest was reported by the authors.

Data availability statement

The data that support the findings of this study are available upon request from the corresponding author (H.H.), upon reasonable request.

Additional information

Funding

This study was funded by the National Research Foundation of Korea (NRF-2021R1C1C1009621 and No. 2019R1A6A1A03033167) and Ministry of the Interior and Safety as Earthquake Disaster Prevention Human resource development Project.