# Internal deformation of the North Andean Sliver in Ecuador-southern Colombia observed by InSAR

L. Marconato, L. Audin, M. -P. Doin, Jean-Mathieu Nocquet, P. Jarrin, F.

Rolandone, N. Harrichhausen, P. Mothes, H. Mora-Páez, D. Cisneros

► To cite this version:

L. Marconato, L. Audin, M. -P. Doin, Jean-Mathieu Nocquet, P. Jarrin, et al.. Internal deformation of the North Andean Sliver in Ecuador-southern Colombia observed by InSAR. Geophysical Journal International, 2024, 239 (3), pp.1557-1575. 10.1093/gji/ggae338. insu-04726372

HAL Id: insu-04726372

https://insu.hal.science/insu-04726372v1

Submitted on 20 Oct 2024

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

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

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

Distributed under a Creative Commons Attribution 4.0 International License
# Internal deformation of the North Andean Sliver in Ecuador and southern Colombia observed by InSAR

L. Marconato<sup>1</sup>, L. Audin,<sup>1</sup> M.-P. Doin,<sup>1</sup> J.-M. Nocquet<sup>2,3</sup> P. Jarrin<sup>2</sup>, F. Rolandone,<sup>4</sup> N. Harrichhausen<sup>1,5</sup> P. Mothes<sup>6</sup>, H. Mora-Páez<sup>7</sup> and D. Cisneros<sup>8</sup>

Accepted 2024 September 12. Received 2024 September 11; in original form 2024 July 4

## SUMMARY

In the northern Andes, partitioning of oblique subduction of the Nazca plate beneath the South American continent induces a northeastward motion of the North Andean Sliver. The strain resulting from this motion is absorbed by crustal faults, which have produced magnitude 7+ earthquakes historically in the Andean Cordillera of Ecuador and southern Colombia. In order to quantify the strain in that area, we derive a high-resolution surface velocity map using InSAR time-series processing. We analyzed 6–8 yr of Sentinel-1 data and combined different satellite line-of-sight directions to produce a reliable velocity map in the east direction. We use interpolated GNSS data to express the velocity map with respect to Stable South America and remove the long-wavelength pattern due to the postseismic deformation following the 2016 M<sub>w</sub> 7.8 Pedernales earthquake. The InSAR velocity map finds high east-west shortening strain rates along north-south trending structures within the Western Cordillera and the Interandean valley, with little deformation taking place east of them. This result strengthens the previous proposition of a ~350 km long Quito-Latacunga tectonic block, forming a restraining bend in the overall right-lateral strike-slip fault system accommodating the northeastward escape motion of the North Andean Sliver. However, the high spatial resolution provided by InSAR indicates that previously proposed boundaries for this block need to be revised. In particular, InSAR results highlight high strain rates (>300 nstrain yr<sup>-1</sup>) along undescribed active structures, south and west of the proposed limits for the Quito-Latacunga block, respectively, in Peltetec and Ibarra regions. Interestingly, the two areas with the largest strain rates spatially correlate with the proposed areas of large historical earthquakes. Modeling of the InSAR and GNSS velocities in these areas suggests shallow coupling and high slip rates on structures which, previously, were not identified as active. We also demonstrate a slow-down of the shallow aseismic slip on the Quito fault after the Pedernales earthquake, suggesting that stress changes following large megathrust events might trigger transient slip behaviors on crustal faults. The high-resolution strain map provided by this work provides a new basis for future tectonic models in the Ecuadorian and southern Colombian Andes, and will contribute to the seismic hazard assessment in this highly populated area of the Andes.

**Key words:** Radar interferometry; Satellite Geodesy; South America; Time-series analysis; Earthquake hazards; Continental neotectonics.

<footnote n="1"><sup>1</sup>Univ. Grenoble Alpes, Univ. Savoie Mont Blanc, CNRS, IRD, Univ. Gustave Eiffel, ISTerre, 38000 Grenoble, France. E-mail: leo.marconato@univ-grenoble-alpes.fr</footnote>

<footnote n="2"><sup>2</sup>Université Côte d'Azur, IRD, CNRS, Observatoire de la Côte d'Azur, Géoazur, 06560 Valbonne, France</footnote>

<footnote n="3"><sup>3</sup>Université Paris Cité, Institut de Physique du Globe de Paris, 75005 Paris, France</footnote>

<footnote n="4"><sup>4</sup>Sorbonne Université, Institut des Sciences de la Terre Paris, ISTE, UMR 7193, 75005 Paris, France</footnote>

<footnote n="5"><sup>5</sup>Department of Geological Sciences, University of Alaska Anchorage, Anchorage, 99508 AK, USA</footnote>

<footnote n="6"><sup>6</sup>Escuela Politécnica Nacional, Instituto Geofisico, Ladrón de Guevara E11-253, Apartado 2759 Quito, Ecuador</footnote>

<footnote n="7"><sup>7</sup>Colombian Geological Survey, Space Geodesy Research Group, Diagonal 5334-53 Bogotá D.C., Colombia</footnote>

<footnote n="8"><sup>8</sup>Instituto Geográfico Militar de Ecuador IGM-EC, Av. Seniergues E4-676 y Gral. Telmo Paz y Minõ, Sector El Dorado, Quito, Ecuador</footnote>
# 1. INTRODUCTION

Continental deformation in the northern Andes arises from the subduction of the Nazca plate along the western margin of Ecuador and Colombia, the subduction of the Caribbean plate beneath the northern margin of Colombia, and the collision of the Panama-Choco block with northwestern Colombia (Pennington 1981; Kellogg & Bonini 1982). Both widespread shallow seismicity spanning from the trench to the east of the Andean Cordillera (e.g. Guillier et al. 2001) and geomorphic markers (e.g. Winter et al. 1993; Lavenu et al. 1995; Audemard et al., 2003) indicate that active crustal deformation affects a broad area, including the Andean Cordillera, the sub-Andean thrust-and-fold belts, and parts of the coastal region.

Early work by Pennington (1981), later confirmed by direct GPS measurements (Trenkamp et al. 2002; Nocquet et al. 2014; Mora-Paéz et al. 2016), shows that most patterns of the present-day deformation can be explained by the northeastward escape motion of a continental domain known as the North Andean Sliver (NAS), squeezed between the Nazca and South America plates (Fig. 1). While the western boundary of the NAS frequently experiences large megathrust earthquakes (e.g. Kanamori & McNally 1982; Nocquet et al. 2017), crustal faults along its eastern boundary have also produced M7+ earthquakes over the past centuries, including some of the largest known in the Andes (e.g. Beauval et al. 2010; Audemard et al. 2021). Regional-scale GPS measurements have quantified the kinematics of the NAS, showing a northeastward movement at a rate of 8–10 mm yr⁻¹ relative to South America, consistent with a partial (~25 percent) partitioning of the oblique convergence vector of the Nazca plate towards South America (Nocquet et al. 2014; Alvarado et al. 2016; Mora-Páez et al. 2019). Besides its translational motion, the NAS also experiences considerable internal deformation, especially in northern Colombia, where the collision with the Panama Arc adds an eastward component to its overall motion, resulting in rapid contractional strain (Mora-Paéz et al. 2019; Jarrin et al. 2023). Another area of internal deformation is the Andean Cordillera in Ecuador and southern Colombia, which is the focus of this study.

In Ecuador, previous active tectonic studies and GPS data have proposed a boundary separating the NAS from the Inca or Peru sliver, which moves southeastward relative to South America at approximately 4–5 mm yr⁻¹ (Nocquet et al. 2014; Villegas-Lanza et al. 2016; Fig. 1). This boundary extends from the Gulf of Guayaquil, runs along the Puna and Pallatanga faults striking northeast, crosses the Andes, and then joins the Cosanga and Chingual faults along the eastern flank of the Eastern Cordillera in a north-south (N–S) direction (Alvarado et al. 2016; Fig. 1). This boundary is supported by evidence of rapid slip rates from morphological feature offsets along the Puna (Dumont et al. 2005) and Pallatanga faults (Winter et al. 1993; Baize et al. 2020), paleoseismological studies (Baize et al. 2015), sharp gradients across the Pallatanga fault observed in GPS data (Nocquet et al. 2014; Alvarado et al. 2016), the occurrence of large historical earthquakes such as the 1797 Riobamba earthquake (Giesecke et al. 2009; Beauval et al. 2010), and current seismicity (Beauval et al. 2013; Champenois et al. 2017; Vaca et al. 2019).

Concurrently, other observations indicate significant internal deformation within the Andean Cordillera, west of the proposed NAS eastern boundary. Lavenu et al. (1995) showed evidence of east-west (E–W) shortening since the early Pleistocene across the Interandean depression separating the Western and Eastern Cordilleras north of the Pallatanga fault. Farther north, geomorphological studies (Alvarado et al. 2014), GPS, and InSAR data

reveal active E–W shortening across the Quito fault at a rate of 4 mm yr⁻¹ (Marinière et al. 2020), about half the rate of the NAS translational motion. Historical records further show that significant earthquakes occurred within the NAS, like the Ibarra 1868 M ~7 earthquake (Beauval et al. 2010). This internal deformation of the NAS has been interpreted within the framework of a large restraining bend (Ego et al. 1996; Alvarado et al. 2016). Indeed, the northeastward motion of the NAS results in almost pure right-lateral strike-slip motion along the Puna-Pallatanga faults and the Afiladores-Sibundoy fault in southern Colombia (Fig. 1). Along the N–S trending Cosanga and Chingual faults, the change of strike results in a transpressive regime along a ~350 km long restraining bend. Alvarado et al. (2016) proposed that the restraining bend narrowed throughout time since 15 Myr leading to the current deformation, idealized by an individual Quito-Latacunga block (QLB; Fig. 1). The QLB is delimited to the west by a series of active faults along the western slope of the Interandean Valley. Using denser GPS data, Jarrin et al. (2023) further developed this proposition in an elastic block model at the scale of the northern Andes showing fragmentation of the NAS in several regions. For Ecuador and southernmost Colombia, their model suggests that 3–3.5 mm yr⁻¹ of E–W shortening is required west of the Interandean depression and 3 mm yr⁻¹ of right-lateral motion north of Quito. The strength of the elastic block modeling approach is that it ensures consistency of slip rates along major faults with the larger-scale kinematics. However, GPS block modelling has some limitations: (1) restraining bends are caused by inelastic deformation within the crust, potentially causing diffuse and spatially variable strain within areas; spatially variable strain is difficult to account for in the elastic block modelling approach and may result in biased fault slip rate estimates; (2) block boundaries are subjective, especially when GPS data are sparse and no active tectonics data are available; (3) GPS data are too sparse to estimate the locking depth along major faults, which was chosen fixed in their model, except for the Quito fault. In addition, many data used in the Cordillera were from survey mode observations, making it impossible to evaluate any time-dependent change of fault behavior through time.

InSAR from Sentinel-1 satellites, with a high spatial resolution (~100 m) and acquisition intervals of every 6–12 d, offers the ability to provide a continuous image of deformation, overcoming some limitations of GNSS data for crustal deformation studies. Sentinel-1 satellites launched in 2014 and 2016 now provide a wealth of data, allowing the computation of time-series in areas previously considered unsuitable for InSAR due to low coherence, steep slopes, or atmospheric noise (e.g. Mathey et al. 2022; Daout et al. 2023; Maubant et al. 2023). So far, an InSAR velocity map has not been produced in the equatorial Andes, possibly because of multiple challenges: vegetation growing up to 4000 m elevation results in poor phase coherence, strong topographic gradients complicate unwrapping, and a very turbulent troposphere increases noise in both space and time. Special care is therefore required during InSAR processing in this area to produce robust velocity maps, particularly regarding the optimization of coherence and mitigation of fading signals.

In this study, we produce the first high-resolution velocity maps of the Ecuador and southern Colombia Cordilleras, based on time-series of ~8 yr of Sentinel-1 data. The analysis period encompasses the Pedernales megathrust earthquake (2016, M<sub>w</sub> 7.8), which caused postseismic deformation that is superimposed on the secular interseismic strain of crustal faults in the Ecuadorian Cordillera. We describe our strategy for correcting the postseismic transient deformation to extract signals associated with crustal faults. We compare
<figure><img src="image_2.png" /><figcaption>**Figure 1.** Simplified tectonic map of the study area and data coverage. The block boundaries from the model of Jarrin *et al.* (2023) are shown in dashed black lines, with block names in capital italic: NAS: North Andean Sliver, QLB: Quito-Latacunga micro-block (named UIOL in Jarrin *et al.* 2023); SUB: Subandean, PS: Peruvian (or Inca) Sliver; SSA: Stable South America. Arrows show the average direction and velocity of motion of each block with respect to SSA from Jarrin *et al.* (2023). Slip model of the Pedernales earthquake from Nocquet *et al.* (2017). The coverage of Sentinel-1 tracks used in this study is shown by white polygons, and the area covered by both ascending and descending geometries is highlighted by a thick white line. Topographic data are from Copernicus GLO-30 DEM, and bathymetric data are from SRTM15+.</figcaption></figure>

our results with the elastic block model from Jarrin *et al.* (2023) and identify several areas of previously unknown high-strain rate. In these high-strain areas, we compute simple fault models to evaluate slip rate and locking depth. Finally, we discuss the localization of deformation around the Quito-Latacunga micro-block and the potential temporal variations in the interseismic behavior of crustal faults.

## 2. PROCESSING OF SENTINEL-1 INSAR TIME-SERIES

The region of interest is covered by Sentinel-1 Interferometric Wide Swath ascending tracks A018 and A120 and descending track D142. The processing of the D142 track is divided into two azimuth subsets

(D142 and D142_North), overlapping by one burst (Fig. 1). All available acquisitions since 2014 (2017 for track A120) to early-2022 or early-2023, depending on tracks are included (see Table 1 for details). We focus on the Andean Cordilleras and discard data from lowlands, which are covered by dense tropical vegetation where the interferometric coherence is extremely low. For the same reason, we also discard one sub-swath from track D142, located in the Amazonian Forest. The resulting area covered by both ascending and descending data is about 100 km wide and 400 km long (Fig. 1). It encompasses the Andean Cordillera from south of the Pallatanga fault, the whole area proposed for the Quito-Latacunga block and extends north of it in southern Colombia.

INSAR time-series are processed using the NSBAS chain (Doin *et al.* 2011), and following a SBAS approach (Ho Tong Minh *et al.* 2022). For each track, interferogram networks are first built by
<table><caption>**Table 1.** Sentinel-1 data used in this study.</caption><thead><tr><th>Track name</th><th>Geometry</th><th>Sub-swaths processed</th><th>First acquisition</th><th>Last acquisition</th><th>Number of acquisitions</th><th>Final number of interferograms</th></tr></thead><tbody><tr><td>A018</td><td>Ascending</td><td>iw1+iw2+iw3</td><td>October 6, 2014</td><td>April 15, 2022</td><td>252</td><td>1420</td></tr><tr><td>A120</td><td>Ascending</td><td>iw1+iw2+iw3</td><td>April 18, 2017</td><td>February 4, 2023</td><td>136</td><td>678</td></tr><tr><td>D142</td><td>Descending</td><td>iw2+iw3</td><td>October 15, 2014</td><td>November 14, 2022</td><td>267</td><td>1506</td></tr><tr><td>D142_North</td><td>Descending</td><td>iw2+iw3</td><td>October 27, 2014</td><td>March 26, 2023</td><td>287</td><td>1456</td></tr></tbody></table>

creating pairs with the next three acquisitions, and then by adding 3 months, 6 months, and 1-year pairs with minimal perpendicular baseline (Figs 2 and S1). This strategy aims at maximizing the global coherence of interferograms while limiting the influence of phase signals associated with soil moisture or vegetation changes and affecting coherently only the short-term interferograms (referred to in the following as fading signals; e.g. Ansari et al. 2021; Espín Bedón et al. 2022). Interferograms are computed from co-registered Single Look Complex images at a resolution of two looks in azimuth and eight looks in range. Geometric and topographic phases are removed using ESA precise orbits and the Copernicus 30 m DEM (European Spatial Agency 2022). Phase jumps between acquisition bursts are corrected using the spectral diversity approach (Grandin 2015) and finally inverted into time-series (Lemrabet et al. 2023). A tropospheric correction from the ECMWF ERA5 global reanalysis is applied (Doin et al. 2009; Jolivet et al. 2011).

The area of interest is primarily located at high elevation. Yet, because of its equatorial location, vegetation is present almost everywhere in the Cordillera. Thus, interferometric coherence is generally medium, but remains good in the urbanized areas. Coherence is also rapidly decreasing along the Cordillera slopes, where the vegetation becomes denser and thicker.

Spatial filtering operations such as multilooking or filtering, required to reduce the phase noise due to decorrelation, respectively, average or smooth the phase information of neighboring pixels, whose quality can vary over short distances. In order to optimize the coherence of interferograms, we implement a strategy to give higher weight to the most robust pixels within a multilooking cell, similar to the one described by Cheaib & Doin (2023). For this approach to be efficient, we must quantify adequately the reliability of each pixel and associate it with a weight.

To do so, prior to the multilooking step allowing to reach 8 looks in azimuth and 32 looks in range, we combine two measurements of the quality of each pixel. First, the colinearity $C^i$, computed for each interferogram $i$ at a 2-looks resolution, reflects the local self-consistency of a phase value with respect to its neighbors (Pinel-Puyssegur et al. 2012). Like the classical interferometric coherence, this indicator is lying between 0 (worst quality) and 1 (best quality), but it does not take the amplitude into account and represents how a single wrapped phase value is close the ones of neighboring pixels. Second, in order to account for the temporal stability of the phase in each pixel independently of its neighbors, we include a proxy for the temporal coherence $\omega$, computed from all subsequent triplets of 2-looks interferograms in the network (Thollard et al. 2021; Doin et al. 2023) and also lying between 0 (worst quality) and 1 (best quality). From these two metrics, we compute a weight used during complex multilooking:

$$W_{\text{2look}}^{i} = \begin{cases} C^{i} & \text{if } \omega < 0.6 \\ [80. (\omega - 0.6) + 1] \cdot C^{i} & \text{if } \omega \ge 0.6 \end{cases} \quad (1)$$

The coefficients of eq. (1) have been chosen empirically such that ~10 per cent of most stable pixels through time ($\omega > 0.6$) are weighted more heavily, by a factor proportional to $\omega$ (Fig. S2),

and such that a pixel with $\omega$ close to 1 has a weight sufficient to dominate the average within the multilooking cell containing 4 × 4 pixels. We thus set the weight of very stable pixels to be about 30 (~2 × 16) times larger than the one of unstable pixels ($\omega < 0.6$). This procedure allows to preserve the reliable phase information of isolated stable pixels.

The phase is then filtered by applying a moving average using a 6 by 6 pixels window weighted by $W_{\text{8look}}^i$, itself filtered to produce $W_{\text{8look-filt}}^i$, a measure of the quality of the filtered phase. Interferograms are then unwrapped, using a region-growing algorithm that expands the unwrapping from high to low values of $W_{\text{8look-filt}}^i$ (Doin et al. 2023). Espín Bedón et al. (2022) showed that using a large proportion of long-temporal baseline (>3 months) interferograms removes a large part of the fading signals in this area. Thus, we ensure that the long-temporal baseline interferograms, which represent about the half of all the interferograms included in the analysis (Fig. 2b), are unwrapped correctly in the Cordillera area. We define an unwrapping threshold as the optimal balance that enables the unwrapping of all interferograms in the Cordillera while avoiding phase unwrapping errors elsewhere. As a result, only high-elevation areas (>2500–3000 m) are unwrapped on all interferograms, although very short baseline interferograms (6 or 12 d) are unwrapped totally in most cases (Fig. S3).

Then, we fit linear phase ramps in the unwrapped interferograms and invert the parameters of these ramps according to the SBAS network. We correct interferograms for the phase ramps computed from inverted parameters, following Biggs et al. (2007). The flattened interferograms are then inverted into time-series (López-Quiroz et al. 2009; Doin et al. 2011; Ho Tong Minh et al. 2022), after having referenced the phase to an area common to all interferograms. We perform a decomposition of the time-series in order to separate a linear trend from seasonal signals, using a double-iteration approach to reduce the weight of outliers and mitigate the influence of acquisitions with large atmospheric phase screen, as described in Lemrabet et al. (2023). The ground displacement in line of sight (LOS) $d_{\text{LOS}}$ at each pixel $(x, y)$ and each time-step $t_k$ is modeled by

$$d_{\text{LOS}}(x, y, t_k) = A_0(x, y) + v_{\text{LOS}}(x, y) \cdot t_k + A_c(x, y) \cdot \cos(2\pi t_k) \\ + A_s(x, y) \cdot \sin(2\pi t_k) \quad (2)$$

where $A_0$ is a constant, $v_{\text{LOS}}$ the linear velocity, $A_c$ the cosine amplitude, and $A_s$ the sine amplitude. For tracks A120 and D142_North, a coseismic step term is added in order to model the offset due to the July 25, 2022 M$_w$ 5.7 earthquake near El Angel. The fit is restricted to a specific time-span that will be discussed in the next section. The resulting velocity map $v_{\text{LOS}}(x, y)$ is masked for pixels with RMS misclosure >0.6 rad, a criterion measuring the pixel-wise misclosure within the interferogram network after time-series inversion (López-Quiroz et al. 2009).
<figure><img src="image_3.png" /><figcaption>Figure 2. Example of interferogram network (a) and proportions of temporal baselines (b) for the track A018.</figcaption></figure>

## 3. POSTSEISMIC DEFORMATION IN GNSS DATA

The period of Sentinel-1 data analyzed in this study encompasses the April 16, 2016, $M_w$ 7.8 Pedernales megathrust earthquake, resulting from a ~100 km rupture of the subduction between 0.5° S and 0.5° N (Nocquet et al. 2017). Though the rupture area is located more than 150 km west of our study area, the Pedernales earthquake produced significative surface deformation within the Andean Cordillera, with a westward coseismic offset of ~6 cm in Quito for instance (e.g. QUEM station at 180 km from the coast; Fig. 3b). GNSS data also show that significant deformation occurred within the Cordillera during the weeks, months, and years following the earthquake, due to both afterslip and visco-elastic relaxation (Rolandone et al. 2018; Tian et al. 2023). More precisely, the first order pattern of postseismic deformation consists of westward displacements, with a magnitude increasing towards the rupture area. In our study area, this results in an average 5 mm yr⁻¹ mostly E–W extensional velocity gradient (Fig. 3a). This extensional deformation tends to hide the long-term interseismic deformation, which is expected to have a significant E–W shortening component because of the transpressive regime within the NAS and along its boundaries.

We thus develop a strategy aiming at reducing the impact of the post-Pedernales earthquake deformation on the InSAR velocity field (Fig. S4). For that purpose, we rely on the continuous GNSS networks deployed in Ecuador by the Instituto Geofísico (IG-EPN), in collaboration the French National Research Institute for Sustainable Development (IRD; Mothes et al. 2013; Alvarado et al. 2018), and the Military Geographical Institute (IGM), and in Colombia by the Colombian Geological Survey (SGC; Mora-Páez et al. 2018). We used displacement time-series derived from a processing using the GAMIT-GLOBK software, according to the strategy described in Jarrin et al. (2023), updated to early 2023. GNSS time-series are corrected for coseismic offsets associated with two deep intraslab earthquakes: $M_w$ 7.3 on 2019/02/22 in Ecuador and $M_w$ 8 in May 2019/05/26 in Peru (Vallée et al. 2023).

GNSS time-series are divided into three periods (Fig. 3b): (1) the pre-Pedernales secular or interseismic velocity, (2) the early postseismic phase where rapid displacement decelerates through time, and (3) the current long-term postseismic phase that can be approximated by a steady velocity, dominated by visco-elastic relaxation (Tian et al. 2023). We empirically determine the limit between the phases 2 and 3 as the earliest date from which the time-series can be well approximated with a constant velocity for the East component. To do so, we iteratively fit a linear trend in the time-series, varying the date at which we start the linear regression, and we compute the

average misfit provided by the weighted root mean square (wrms) for a selection of stations. 57 stations are selected in a box (from −79° to −76° in longitude and from −2° to 1° in latitude) surrounding the study area, but excluding the coast. We find that after July 2017 and for a period of 6 yr, a single velocity fit provides a misfit of wrms = 3.3 mm, similar to the one found for a linear trend fitted before the earthquake (wrms = 3.5 mm, Fig. S5). Our empirical approach is consistent with the results from Tian et al. (2023), who show that after 1 year following the Pedernales earthquake, viscoelastic relaxation becomes the main contribution to postseismic deformation. Using the viscosity inferred from their best-fitting model, i.e. $3 \times 10^{19}$ Pa s Maxwell viscosity, and assuming a 50 GPa shear modulus, the characteristic time for the viscoelastic effect of this earthquake equal to about 20 yr, using equations for a classical Maxwell rheology (Bürgmann & Dresen 2008). Such characteristic time predicts a negligible velocity change during our period 3 (2017.6–2023).

We then aim to compute a differential velocity field of period 3 with respect to period 1. For that purpose, we select stations within a box spanning −81° to −76° in longitude and −4° to 3° in latitude. As a robustness criterion, for each station, a velocity is extracted on the condition that the time-series contains at least 1000 daily measurements (equivalent to a ~3 yr duration) during the period of interest. For period 1, we use the horizontal velocity field from Jarrin et al. (2023), and complement it with the vertical velocities estimated during the 2007–2016.3 period (82 stations with enough data). For period 3, we extract both horizontal and vertical velocities between 2017.6 and 2023 (116 stations with enough data). All velocity fields are then rotated from ITRF2020 (Altamimi et al. 2023) to a Stable South America reference frame using the Euler pole from Jarrin et al. (2023). We eventually compute the difference between the velocities extracted from period 3 and the ones extracted from period 1 at stations where they are both available, i.e. 78 stations for the horizontal components and 72 stations for the vertical component. This differential velocity field (Fig. 3a), representing the ‘anomaly’ due to the postseismic transient during the period 3 with respect to period 1, will be used in the next section to correct InSAR velocity maps.

## 4. POSTPROCESSING OF InSAR VELOCITY MAPS

We estimate velocities from InSAR time-series for period 3 only (from July 2017 to December 2022), because Sentinel-1 acquisitions are very sparse during period 1, while period 2 is highly impacted by nonlinear postseismic deformation. We then tie the flattened LOS
<figure><img src="image_4.png" /><figcaption>**Figure 3.** (a) InSAR velocity map for the east component, computed on period 3. The LOS velocities have been tied to GNSS data and decomposed into east and up directions. This map is not corrected for the postseismic deformation due to the Pedernales earthquake. Arrows show the velocity field corresponding to the postseismic deformation (i.e. the difference between periods 1 and 3, see the text for details). Slip model of the Pedernales earthquake from Nocquet *et al.* (2017). (b) GNSS time-series of eastward displacement at the QUEM station, located in Quito. This allows to distinguish three time periods with distinct deformation behavior (see the text for details).</figcaption></figure>

velocity maps to the GNSS velocity field also computed over period 3 and expressed in the Stable South America reference frame using the pole from Jarrin *et al.* (2023). The density of GNSS stations is globally good but spatially heterogeneous, and some areas in southern Ecuador and southern Colombia have notably fewer stations than in northern Ecuador (Fig. 3a). To derive a homogeneous

velocity field for referencing InSAR velocity maps, we thus interpolate the GNSS velocity field (using parameters in Table S1), once projected to LOS using the 3-D components (Fig. S6). We then fit a surface to the difference between the InSAR velocity maps and the interpolated GNSS LOS map. The surface $d_{ramp}$ is parameterized as a 2-D plane linear in longitude *x* and latitude *y*, plus, for the longer
tracks A018 and D142, a cross term leading to a twisted plane:

$$d_{\text{ramp}} = \begin{cases} A + Bx + Cy + Dxy & \text{for A018 and D142} \\ A + Bx + Cy & \text{for A120 and D142..North} \end{cases} \qquad (3)$$

After referencing the LOS velocity maps by removing this ramp, we compute the RMS of the difference between InSAR and GNSS velocities at the stations (Fig. S7). The agreement is of the order of 2 mm yr⁻¹, consistent with the RMS of residuals found in overlapping regions between adjacent InSAR tracks.

The LOS velocities, now referenced with respect to Stable South America (Fig. S8), are then decomposed into East and Up components, taking advantage of the complementary ascending and descending LOS geometries (Fuhrmann & Garthwaite 2019). As the North component cannot be precisely constrained by the InSAR measurements, we first remove its effect from the velocity LOS maps. To do so, we interpolate the North component of GNSS velocity computed over the period 3, using parameters in Table S1, project it into LOS, and then subtract it from the InSAR velocity maps. The obtained 2017–2023 velocity maps for the East component with respect to Stable South America and Up components with respect to the ITRF2020, with a pixel size of 0.001° (~110 m), are displayed in Figs 3a and S9.

As a final step, we correct the obtained InSAR velocity map for postseismic deformation. The difference of velocity between period 3 (2017–2023) and period 1 (before Pedernales) represents the postseismic perturbation that is superimposed on the ‘secular’ interseismic velocity field. Removing this contribution from our velocity maps should produce a map of the interseismic deformation, as it would have presumably been observed before the earthquake, at least for long wavelength. Given that mantle flow is the main process driving the postseismic deformation during the 2017–2023 period considered here (Tian et al. 2023), the associated pattern must be spatially smooth throughout the Cordillera. We first determine 2017–2023 minus pre-Pedernales GNSS velocities as seen in previous section and perform a smooth interpolation for the East and Up components (Fig. S10), using parameters in Table S1. We finally subtract the interpolated differential velocity from the East and Up InSAR velocity maps to obtain a velocity field for the Ecuador and southern Colombia Cordillera, that is representative of the interseismic deformation.

# 5. RESULTS

The InSAR vertical velocity map (Fig. S11) is noisier than the horizontal one. The noise sources impacting the vertical component have the same effect in both the ascending and descending LOS maps (Figs 4a,d). Consequently, during the East and Vertical decomposition, these common signals preferably go in the vertical velocity map (Fig. 4b), while having a negligible impact on the East velocity map (Fig. 4e). For instance, areas showing negative sign in both ascending and descending LOS velocity maps (i.e. apparent motion away from the satellite) appear to be correlated with the land cover type (Fig. 4c). They especially have strong magnitudes in ‘áparo’, a tundra-like soil and ecosystem typical of the tropical Andes above 3000 m elevation (Fig. 4f). We thus think that Sentinel-1 velocities in these areas are actually dominated by residual fading signals (Ansari et al. 2021; Maghsoudi et al. 2022), and should not be interpreted in terms of tectonic deformation. By contrast, Fig. 4e shows that the East velocity map displays smooth variations, independent of the type of land cover.

For such reasons, in the following sections, we focus on the East velocity map, which allows quantification of the interseismic deformation. In Fig. 5(a), we show the East component of velocity predicted by the elastic block model from Jarrin et al. (2023), together with the GNSS velocity used as input of their model. Our East velocity map from InSAR is shown in Fig. 5(b). Fig. 6 shows the comparison of the East velocities for InSAR results, GNSS data and predictions from Jarrin et al. (2023) model along (40 km wide x 200 km long) profiles crossing the Andean Cordillera. The RMS of the difference between InSAR and block model prediction is 1.77 mm yr⁻¹. The main pattern highlighted by InSAR is a decrease of eastward velocity towards the east-southeast, consistent with transpressive deformation in this part of NAS. InSAR results further show significant variations of EW shortening rates from 2 mm yr⁻¹ in the southern part of the study area (profile 12) and up to 8 mm yr⁻¹ in southern Colombia (profile 1).

The substantial increase in spatial resolution brought by InSAR compared to GNSS allows us to identify local velocity gradients which were not previously resolved using GNSS data alone. In Fig. 6, some profiles show areas of enhanced velocity gradient, witnessing localized deformation on active faults. An eastward drop of East velocity can reflect either shortening across a N-S oriented fault, right-lateral strike-slip along a NE-SW oriented fault, or in between these two end-members, transpressional deformation. Several areas of localized velocity gradients are consistent with the prediction of the block model from Jarrin et al. (2023). However, InSAR also identify some areas that differ from the elastic block model prediction, either regarding the location of strain or the sharpness of the gradient. Hereafter, we review from north to south the new insights brought by our results on the main faults of the area, especially the strands defining the western boundary of the Quito-Latacunga micro-block (QLB) where InSAR provides a continuous image of deformation.

## 5.1 Evidence for activity along the Romeral fault in southern Colombia

In southernmost Colombia, the motion of the NAS is thought to be accommodated mostly by the Afiladores-Sibundoy fault, the main structure delimiting the eastern border of the NAS (Velandia et al. 2005; Tibaldi et al. 2007; Audemard et al. 2021). Our Sentinel-1 InSAR map shows little reliable data east of the Afiladores fault, but still suggests a ~4 mm yr⁻¹ velocity gradient across it (profile 1). In that area, the most significant pattern is a shortening strain west of the Afiladores fault, around 7–8 mm yr⁻¹ accommodated over a distance of 30 km. That area might correspond to the strain induced by the Romeral fault in southern Colombia (e.g. Paris et al. 2000; Audemard et al. 2003) or more diffuse strain in that area. Interestingly, Jarrin et al. (2023) tested two possible geometries for that area: one considering that the Romeral fault merges with the QLB (model C in Jarrin et al. 2023) and one discarding the Romeral fault as a significant active fault in southern Colombia (model D). Model C was discarded because it would predict 6–7 mm yr⁻¹ shortening rates to be accommodated by the Romeral fault, a value much higher than the ~1 mm yr⁻¹ previously proposed rates based on geology (Paris et al. 2000). The InSAR results here strongly support an active deformation along the southern segment of the Romeral fault, possibly merging with the western boundary of the QLB in Ecuador. The plausibility of Jarrin’s model C should thus be reassessed given the new information from InSAR.
<figure><img src="image_5.png" /><figcaption>**Figure 4.** Comparison of LOS velocity maps (descending and ascending in a and d, respectively) with decomposed velocity maps (vertical and East component in b and e, respectively), an optical satellite image (c) and ecosystem classification (f). In both LOS maps, residual fading signals are negative and located on 'páramo' ecosystems, they are thus mostly visible in the Vertical velocity map, leading to a less noisy East velocity map.</figcaption></figure>

## 5.2 Northern boundary of the Quito-Latacunga block: the El Angel fault

Directly to the south, the WSW-ENE segment labeled as *El Angel* in Fig. 5 defines the northern boundary of the QLB. No major active structure has been identified in this area, and a recent study proposes instead that it could be a distributed boundary (Harrichhausen *et al.* 2023b). That view is consistent with our results, which show no localized velocity gradient, but a rather constant velocity slope (Fig. 5b, *profile 2*). However, this area is close to the Chiles volcano, which has been undergoing inflation since 2014 (Fig. 3a) in correlation with seismo-tectonic activity (Ebmeier *et al.* 2016). Our velocity map reveals a maximum uplift rate of 2.9 cm yr⁻¹ during period 3, S-E of Chiles volcano (Fig. S11). Volcanic deformation most probably overprints part of the tectonic signal in this area.

minor evidence of active faulting, delimits the Western Cordillera from the Interandean valley. According to our InSAR results, this area shows a strong and localized velocity gradient. The offset in East velocity is about 4–5 mm yr⁻¹, which is larger than the prediction of the block model. This step is accommodated over a short distance of less than 20 km, suggesting a shallow locking depth. By contrast, the elastic block model predicts a much smoother velocity gradient, that does not fit well the GNSS East velocities either, especially on *profile 3*. As such, the observed pattern is similar to the observations of shallow locking reported by Marinière *et al.* (2020) for the Quito fault. The InSAR map also allows here a more precise location of the western boundary of the QLB, which is for example here offset by about 15 km to the west with respect to Jarrin’s block model on *profile 3*.

## 5.3 Localized deformation in the Ibarra area

The western boundary of the QLB around latitudes 0° to 0.5° N (*profiles 3* and *4*), along the segment labelled as *Ibarra* in Fig. 5, is not very well defined from a tectonic point of view. Very sharp scarps have been observed within the Billecocha plateau, with evidence of both normal, reverse and strike-slip faulting, potentially related to both gravitational and tectonics effects (Jomard *et al.* 2021). East of the Billecocha fault system, a large topographic step, associated with

## 5.4 The Quito fault system

This reverse fault system lies under the city of Quito (Alvarado *et al.* 2014). Marinière *et al.* (2020) found that it was mostly creeping aseismically using GNSS data acquired between 1994 and 2016, and PS-InSAR observations from ERS scenes acquired between 1993 and 2000. They find a horizontal shortening rate of 3–5 mm yr⁻¹ accommodated over a width of about ~10 km. This feature was taken into account in the block model of Jarrin *et al.* (2023), whose prediction shows a sharp gradient in the East velocity (*profile 5*).
<figure><img src="image_6.png" /><figcaption>**Figure 5.** Velocity comparison between InSAR, GNSS data, and block model prediction. (a) Prediction of the block model by Jarrin *et al.* (2023) for the East component. The GNSS velocity field (computed before 2016) used to constrain the model is shown by the arrows. The color of the dots represents the East component. Boundaries and names of the blocks are shown by dashed black lines and black labels (see the caption of Fig. 1 for abbreviations). Fault names are written in white. (b) Interseismic InSAR velocity map from our study, for the East component, with the same colorscale as in (a). The color dots are the same as in (a). We mask the area of the Chiles volcano because the east component there is not well constrained. The footprints of Figs 4, 7, and 8 are shown by white rectangles.</figcaption></figure>

However, the Sentinel-1 InSAR velocity map only displays a very smooth eastward decrease of the East velocity, with no clear gradient along the Quito fault. This result appears thus at odds with the studies using pre-2016 GNSS velocities. We comment on the possible origin of this discrepancy in the discussion section.

## 5.5 Latacunga area

South of the Quito fault system, the main tectonic structure proposed to act as the western boundary of the QLB is the Latacunga fold-and-thrust belt (Lavenu *et al.* 1995; Baize *et al.* 2020; Fig. 5). There,

our results show a very good agreement with the GNSS data and the elastic block model predictions. One can observe along *profiles* 6 and 7 that the deformation is distributed, with about 4–5 mm yr⁻¹ of E–W shortening accommodated by the fold-and-thrust belt across a ~80 km width within the QLB. This is consistent with morphological analysis, showing a series of parallel folds oriented ~N–S (Lavenu *et al.* 1995; Fiorini & Tibaldi 2012). A slight increase in the velocity gradient can be noted at the western boundary of the QLB block (at ~−60 km on *profiles* 6 and 7). This suggests strain localization on a major structure at the block boundary, but it cannot be thoroughly resolved due to the lack of InSAR data west of this boundary.
<figure><img src="image_7.png" /><figcaption>**Figure 6.** Swath-profiles numbered from 1 to 12, whose footprints are shown by black rectangles in Fig. 5. The point clouds and the large dots represent the East velocities from InSAR (whose median is shown by black dashed curve) and GNSS data respectively, both coded with same colorscale as in Fig. 5. The white envelopes show the East velocities predicted by the block model. Vertical dashed black lines show the boundaries of blocks, with the corresponding fault names in white.</figcaption></figure>

## 5.6 Pallatanga and Peltetec faults segments

Profiles 8 and 9 cross the proposed boundaries for the QLB to the north, for the NAS to the west, and for the Inca/Peru sliver (Nocquet *et al.* 2014; Villegas-Lanza *et al.* 2016) to the southeast. The main active structure in that area is the Pallatanga fault, which hosted the 1797 M<sub>w</sub> 7.6 Riobamba earthquake (Giesecke *et al.* 2009; Beauval *et al.* 2010; Baize *et al.* 2015; Alvarado *et al.* 2016; Baize *et al.* 2020). This fault has a clear geomorphological expression (e.g. Winter *et al.* 1993; Baize *et al.* 2015, 2020), running from the Guayaquil basin to the Eastern Cordillera with a SW-NE orientation and was chosen as a block boundary in the model of Jarrin *et al.* (2023). In agreement with this previous finding, profile 8 shows a clear velocity gradient following closely the GNSS data, centered on the surface trace of the Pallatanga fault. However, on profile 9, the gradient appears to be stronger than the one observed in GNSS data and the prediction from the elastic block model. Moreover, the maximum gradient is geographically offset, as it is located 20 km east of the Pallatanga fault and follows instead another active fault colocated with the Peltetec melange and suture zone (Alvarado *et al.* 2016; Fig. 6).

South of this area, on Fig. 6's profile 10, the velocity map does not show any significant localized gradient neither on the Pallatanga fault, nor on the Peltetec suture zone, but shows instead a rather smooth velocity variation. On profiles 11 and 12, the imaged area does not cross any boundary of the block model, and only

encompasses the Inca/Peru sliver. There also, the InSAR velocity shows no significant gradients. The velocity trend decreases towards the South, with no evidence of internal deformation in the Peruvian sliver away from its limit with the NAW block, in agreement with the block model prediction.

## 6. MODELING

Our analysis of InSAR results highlights two areas where previously unknown high velocity gradients take place, namely the Ibarra and Pallatanga-Peltetec areas, north and south of the proposed Quito-Latacunga block, respectively. While a full model is beyond the scope of this study, we aim to constrain the main fault characteristics (location, geometry, slip rate, and locking depth) using 2-D models using InSAR and GNSS velocities projected along profiles perpendicular to the main structures. Since the shallow locking depth expected from the observed localized strain might be transient, we perform the modeling using the InSAR and GNSS data without the correction of the long wavelength postseismic deformation. Consequently, the results presented here represent the strain around crustal faults during the 2017–2023 period, potentially impacted by stress changes due to faraway megathrust processes. We model local and relative motions within the studied profiles, and thus we believe that this smooth postseismic deformation pattern included in the data does not impact the modeling results.
## 6.1 Modeling approach

We model the surface deformation observed along a two-dimensional profile using a Bayesian inversion approach fully described in Daout *et al.* (2016a, b). This approach explores a fault geometry in a convergent setting, composed of a flat decollement (FD) and a dipping ramp (e.g. Marinière *et al.* 2020). The explored parameters are the depth of the decollement (*H*), the horizontal position of the decollement's junction with the ramp (*D*), the shortening rate (*short*) and the strike-slip rate (*FDss*) on the FD, the strike-slip rate on the ramp (*Rss*), the depth of the ramp tip (*h*) and the length of the ramp (*L*) (e.g. Fig. 7d). To be sure of consistent kinematics, the slip rate along the ramp is forced to short/cos(θ) where θ = sin⁻¹[(H − h)/L] is the dip of the ramp. The fault strike is set perpendicular to the profile. We constrain the model using both the InSAR East velocity map, without correction of the postseismic deformation associated with the Pedernales earthquake, and the horizontal GNSS velocities extracted from the post-2017.6 period. We choose uniform prior distributions for all parameters, and sample the parameter space with 10⁵ samples, rejecting the first half of them. Although the Bayesian approach results in numerous models sampling the posterior covariance, we choose the means of marginal posterior probability density functions as a set of preferred parameters (Figs S12 and S13) to compute synthetic displacement rates along the profile for plots (Figs 7 and 8). Uncertainties on parameters are given for a 95 per cent confidence level (two times the standard deviation).

## 6.2 Pallatanga-Peltetec profile

The profile is perpendicular to the main orientation of the velocity gradient, around N15° (Fig. 7). GNSS stations located close to Tungurahua volcano were avoided, since they show velocities inconsistent with the regional pattern, probably because of deformation induced by volcanic processes. The best fit fault geometry has a decollement depth around 7.3 ± 2.9 km depth, a shallow dipping ramp (13°) slipping up to 4.8 ± 3.7 km below the surface. The decollement accommodates 7.0 ± 1.2 mm yr⁻¹ of shortening, and 6.1 ± 0.9 mm yr⁻¹ of right-lateral slip, while the right-lateral rate on the ramp is 2.0 ± 2.1 mm yr⁻¹.

In our preferred model, the decollement-ramp junction is located roughly on the trace of the Pallatanga fault, and the upper tip of the ramp is very close to the Peltetec suture (Fig. 7a). While the Pallatanga fault is considered as the main and most studied active structure in this area (Baize *et al.* 2015, 2020; Harrichhausen *et al.* 2023a), the Peltetec segment corresponds to a major suture zone and ophiolite belt delimiting the Interandean Valley from the Eastern Cordillera from a geological and tectonic point of view (Litherland & Aspden 1992). Although the reactivation of sutures in Ecuador was suggested based on the seismicity (Guillier *et al.* 2001), there exists no study showing the seismotectonic activity of the Peltetec segment. Our InSAR data suggests that the current strain accumulation south of the Quito-Latancunga block occurs not only along the Pallatanga fault, but also at shallow depth close to the Peltetec suture. The amount of right-lateral slip predicted by our model on the ramp, close to the Peltetec suture, is not precisely constrained but it appears to be significantly lower than the slip rate accommodated on the decollement. This difference suggests that right-lateral slip occurs on the Pallatanga fault, and thus a possible partitioning of deformation between the two structures. The Pallatanga fault accommodating most of strike-slip motion would be locked at a depth of about 7 km, while the Peltetec tectonic belt, would accommodate

shortening, with a shallower locking depth (5 km). However, the ratio between shortening and strike-slip rates is hard to constrain, because it strongly depends on the fault strike. Since the Pallatanga and Peltetec faults are not parallel, it would require more complex 3-D multifaults models to constrain the exact rates accounting for slip-partitioning.

In that area, the very different geometries used for analyzing this area in this work and in the block model by Jarrin *et al.* (2023) do not allow a direct comparison of the proposed slip rates. Nevertheless, the net amount of slip predicted by our best model (~9 mm yr⁻¹) is higher than the relative motion of the NAS with respect to the PFW block estimated in their in their study (7–8 mm yr⁻¹). The pre-Pedernales GNSS-derived slip rates were already higher than the geological slip rates on the northern segment of the Pallatanga fault compiled by Harrichhausen *et al.* (2023a), accounting for strike-slip motion only (2.4–6.6 mm yr⁻¹ according to geomorphology and 0.3–4 mm yr⁻¹ according to paleoseismology). Our joint InSAR and GNSS modeling confirms that geodetic slip rates are at odds with the geological ones, and attests to even higher slip rates during the observation period of our study (2017–2023) than the ones inferred from pre-2016 GNSS data.

The discrepancy between our results and the block model may be associated with the difference of scale between both models, since the one from Jarrin *et al.* (2023) is a large-scale model, fitting the GNSS velocities of entire blocks, whereas we are here focusing on a small portion of fault. In this case, our dataset better constrains the slip rates locally than the block model. Indeed, the modeled profile includes InSAR data at least 30 km away from the maximum gradient on each side, and quite evenly distributed GNSS values (Fig. 7). An alternative explanation for such discrepancy could be a transient deformation occurring in the observation period of our study, and more specifically an increase of slip rate in the aftermath of the Pedernales earthquake. One way to discriminate between these two options, would be to build a model using high-resolution InSAR data from before the Pedernales earthquake (2016), and compare the results with the block model, which uses geodetic data from the same time span. Such datasets exist, derived from Envisat data (Champenois *et al.* 2014) or ALOS-1 data (Marconato *et al.* 2024), but the velocities in this area are dominated by a large-scale radial inflation induced by the Tungurahua volcano activity, preventing a reliable interpretation of the observed gradients.

In conclusion, InSAR data emphasizes the major role played by the Peltetec suture zone in the accommodation of the shortening south of the QLB. This suggests the possibility of a QLB boundary actually located south of the Pallatanga fault, or distributed between two sub-parallel structures. Also, since InSAR data covers only a short time span, possibly impacted by transient stress changes due to the Pedernales earthquake, it is possible that our results witness a transient regime of strain, whose localization may be varying through time within this complex system.

## 6.3 Ibarra profile

For this area, the lack of data in the far field west of the location of the maximum velocity gradient prevents us from constraining the depth and slip rate of the flat decollement. Therefore, we fixed the H parameter to 5 km, and explored the other parameters. The best geometry is a shallow dipping ramp (17°), slipping up to 2.8 ± 1.1 km below the surface (Fig. 8c). The junction between the flat decollement and the ramp is located below the Billecocha plateau (Fig. 8a),
<figure><img src="image_8.png" /><figcaption>**Figure 7.** Fault modeling results for a profile crossing the Pallatanga-Peltetec high-strain area. (a) InSAR East velocity map and GNSS velocity field used to constrain the model, both computed on the post-2017.6 period. The location of the decollement-ramp intersection and up-dip limit of the ramp are shown by the thick and thin red lines, respectively. Black dashed lines show the block limits from Jarrin *et al.* (2023). White dashed lines show the sutures from Jaillard *et al.* (2022). Thin black lineaments show the active structures (Alvarado *et al.* 2014, 2016; Baize *et al.* 2020). (b) Swath-profile showing the East velocity for InSAR (dots), GNSS (crosses) and the best-fitting model (dashed red line). (c) Modeled displacement rates for fault parallel (blue) and fault perpendicular (green) components, with horizontal GNSS velocities projected in the corresponding components (circles). (d) Topographic profile showing the model parameterization. The best fit dip-slip (ds) and strike-slip (ss) slip rates are given for each fault segment.</figcaption></figure>

in which strong evidence of active faulting was observed, but the mechanisms explaining the mixed kinematics in this area are not clear yet (Jomard *et al.* 2021). The fault tip is located below the topographic front limiting the Interandean Valley at west, where again, Jaillard (2022) draw a main suture zone. Yet, no major active structure has been mapped in that area, although small reverse

faults were observed affecting the Quaternary volcano sedimentary landscape (Alvarado *et al.* 2014, 2016).

The best fit model predicts 8.0 ± 1.6 mm yr⁻¹ of shortening, and 6.6 ± 1.9 mm yr⁻¹ of right-lateral slip, but they are poorly constrained due to the lack of data NW of the decollement tip. In contrast, for comparison, the equivalent fault segment (El Angel
<figure><img src="image_9.png" /><figcaption>**Figure 8.** Same as Fig. 7 for a profile crossing the Ibarra high-strain area.</figcaption></figure>

fault) included in the block model of Jarrin *et al.* (2023) shows 2.4 mm yr⁻¹ of shortening and 2.9 mm yr⁻¹ of right-lateral slip, with a locking depth of 12 km and 50 per cent of coupling. Our model thus predicts slip rates more than three time larger than the one inferred from the elastic block model based on pre-2016 GNSS. First, the absence of data further west than −25 km probably yields overestimated slip rates by about 30 per cent. However, the InSAR and GNSS data both points towards slip rates larger than ~5 mm yr⁻¹ (Fig. 8c), which is still higher than the block model prediction. The remaining discrepancy may be associated with the

difference of scale between models, as mentioned for the Pallatanga-Peltetec profile. A last explanation could be the occurrence of a transient increase of slip rates of this section of the fault system over the observation period (2017–2023), potentially linked with the Pedernales earthquake.

Overall, the Ibarra area exhibits a strong deformation gradient along poorly mapped structures striking ~N30° and following the main topographic front between the Western Cordillera and the Interandean Valley. The accommodation of ~6 mm yr⁻¹ of East displacement over a 15 km distance requires a very shallow locking
depth (2.5–3 km). Due to lack of data in the far-field west of the main gradient, it is not possible to constrain accurately the slip rates for this specific area. The ranges of fault slip rates predicted by our fault models and block modeling (Jarrin et al. 2023) are 2.4–8.0 mm yr⁻¹ for shortening and 2.9–6.6 mm yr⁻¹ for right-lateral strike-slip.

# 7. DISCUSSION

## 7.1 Deformation mode within a restraining bend

Deformation within the Ecuadorian Cordillera has been interpreted as a large restraining bend, marking the transition between strike-slip segments: the Pallatanga fault in the south, which delimits the boundary with the Inca/Peru sliver, and the Afiladores-Sibundoy fault in the north, which accommodates the sliver motion with respect to South America. In a restraining bend, part of the strike slip component converts into a shortening component accommodated at reverse faults oblique to the sense of the strike-slip motion (e.g. Ego et al. 1996). In the case of Ecuador, previous models have proposed that the restraining bend is actually made of a rigid block localizing deformation along its boundaries (Ego et al. 1996; Alvarado et al. 2016). Jarrin et al. (2023) followed this proposition in the elastic block modeling. Their model predicts that 2–3 mm yr⁻¹ of shortening should be accommodated across the western boundary of the QLB along the Latacunga-Quito-Ibarra faults and 2 mm yr⁻¹ along the eastern boundary along the Cosanga fault, the latter also accommodating in addition 6 mm yr⁻¹ of right-lateral strike-slip motion.

InSAR, by providing a continuous image of the velocity field allows a test of whether deformation is localized on major structures or is diffuse within the restraining bend. In order to specifically answer this question, Figs 9(b)–(d) show the norm of the gradient vector for the InSAR East component velocity map. For that purpose, the velocity map is previously smoothed using a Gaussian filter with a critical distance of 5.5 km to smooth out short-scale correlated noise. For comparison, the same quantity is calculated from the prediction of Jarrin et al. (2023) elastic block model (Fig. 9a). The elastic block model predicts strain rate of the order of 150–300 nstrain yr⁻¹ at the block boundaries and <100 nstrain yr⁻¹ inside the QLB. Fig. 9(b) shows that the InSAR strain rate map reveals several areas of high strain rate, with N–S elongated shape, almost continuous over a length of 450 km from latitude 2.5° S to 1.3° N. Aside the two areas of Peltetec and Ibarra studied in the previous section where strain rates reach 300 nstrain yr⁻¹, elongated areas are found south of the Peltetec fault, west of the Interandean depression near Latacunga, Quito, and then parallel to the Chingual and Afiladores-Sibundoy faults in southern Colombia. Such elongated shapes of the strain rate field suggest en-echelon fault segments. East of those linear segments, strain rate inside the QLB only occurs at specific areas. High strain along the Pisayambo segment (Fig. 9b) delineates the extent of the Pallatanga fault in the Eastern Cordillera, where a shallow strike-slip earthquake was studied using InSAR (Champenois et al. 2017). Significant strain rate is found near the Cotopaxi and Cayambe volcanoes, both known to have seismic activity during the period analyzed here (https://www.igepn.edu.ec/red-de-observatorios-vulcanologicos-rovig). Otherwise, low strain rates (<100 nstrain yr⁻¹) are observed within the QLB block, similarly to the prediction from the elastic block model. Hence, linear high strain rate areas surrounding lower strained areas favor the view of localized deformation along the boundaries of relatively rigid QLB.

## 7.2 New constraints on the Quito-Latacunga block boundaries

Jarrin at. (2023) defined block boundaries based on a review of available knowledge of geodesy, active tectonics, geomorphology, paleoseismology and seismicity. Our results indicate that some of the proposed boundary locations need to be revised in the light of the InSAR strain map. For instance, in the Ibarra region, the boundary proposed in Jarrin et al. (2023) is located slightly east of the maximum gradient that we demonstrate. This higher strain rate is located in an area with few recognized active faults and showing a low level of instrumental seismicity (Figs 9b,d). Nevertheless, this area is located along a major tectonic suture zone (Figs 8 and 9; Jaillard 2022) and coincides with the epicentral area of the 1868 M~7 earthquake (Fig. 9c; Beauval et al. 2013). This example illustrates the benefit of incorporating data from multiple sources, including InSAR, to define accurately limits of tectonic blocks.

Similarly, the second high-strain area, close to the southern boundary of the QLB, aligns more with the Peltetec suture zone than with the Pallatanga fault system (Figs 7 and 9). This strain pattern suggests a more complex tectonic setting than previously thought at this triple junction area, implying a potential partitioning of the deformation between the Pallatanga and the Peltetec faults. The Pallatanga fault was considered until now as the main structure of the region (Baize et al. 2015, 2020), notably in the seismic hazard models (Beauval et al. 2018). Our results suggests that the Peltetec fault also accommodates a significant part of the strain (at least during the 2017–2023 period), possibly joining another strike-slip fault near Guamote, south of Pallatanga (Fig. 9b). Together with shortening across the newly identified thrust, partitioning of the strike-strike slip motion between the Guamote-Peltetec and Pallatanga fault would reconcile the lower slip rates found from geomorphic markers compared to geodesy estimates, as suggested by Audemard et al. (2021). Such result has important implications on the slip budget in this area, and thus on seismic hazard assessment.

Although scarce focal mechanisms from global (Ekström et al. 2012) and local (Vaca et al. 2019) catalogs do not allow to draw the continuous limit of the QLB, they still provide some useful indications on the orientation of strain. In particular, several almost purely reverse mechanisms can be observed close to the Peltetec suture area (Fig. 9d), supporting a view of partitioning between E-W shortening and NE-SW right-lateral motion in this area.

## 7.3 Transient behaviors on crustal faults

The Quito fault system benefits from a dense geodetic dataset, both before and after the 2016 Pedernales earthquake. Using 1994–2015 GPS data and ERS 1993–2000 InSAR data, Marinière et al. (2020) documented shallow aseismic slip along the Quito fault. A surprising result of our analysis is that only a moderate velocity gradient is observed in our 2017–2023 Sentinel-1 results (Fig. 6, profile 5). We further checked that our results are consistent with GNSS time-series for our analyzed time period. GNSS time series indeed show a similar change in behavior: relative velocity in the East component between stations QUEM and SALF, located 38 km apart on each side of the Quito Fault System, dropped by about 3.2 mm yr⁻¹ after the Pedernales earthquake (Fig. S14).

In Fig. 10, we compiled the pre-2016 GNSS velocities from Marinière et al. (2020) and Jarrin et al. (2023), and compare them with our GNSS and InSAR data for the 2017–2023 period. The GNSS velocities around the Quito fault tend to diverge towards east, showing a clear radial pattern (Fig. 10a). While it cannot be
<figure><img src="image_10.png" /><figcaption>**Figure 9.** Norm of the East velocity gradient of (a) the block model prediction and (b-d) the InSAR results of this study. Strain maps are compared with (a) boundaries of the block model from Jarrin *et al.* (2023), (b) active faults (Alvarado *et al.* 2014, 2016; Baize *et al.* 2020), (c) sutures mapped by Jaillard (2022) and 95 per cent confidence contour (in blue) for the intensity centers of the two largest historical earthquakes from Beauval *et al.* (2010) and (d) historical seismicity from Beauval *et al.* (2013) (in blue) and shallow (<35 km) focal mechanisms from gCMT catalog (Ekström *et al.* 2012) (in grey) and from Vaca *et al.* (2019) (in black).</figcaption></figure>

<figure><img src="image_11.png" /><figcaption>**Figure 10.** Velocities across the Quito Fault System before and after the Pedernales earthquake. (a) Velocity map from our study, for the East component, without correction of the postseismic effect. Arrows show GNSS velocities from Marinière *et al.* (2020) (in blue), from Jarrin *et al.* (2023) (in red), both computed before 2016, and from this study (in black), computed on the period 3 (2017.6–2023). (b) Swath-profile across the QFS. InSAR data are shown by the colored point cloud, and the dashed black line is a sliding median. East components of GNSS velocity are shown with the same colors as in (a). GNSS velocities are expressed in a local reference-frame defined by stations ANTN and SALF (in dashed circles in panel a). The InSAR profile was re-referenced to the post-Pedernales GNSS East velocities by adjusting a linear trend along the profile.</figcaption></figure>

ruled out that this deformation is due to a far-away and deep vol-
canic source, we favor, just like Marinière *et al.* (2020), that these
radial velocities have a tectonic origin, possibly linked with the
curved shape of the Quito fault. We extract a profile perpendicular
to the Quito fault in order to reduce the potential influence of this
diverging pattern (Fig. 10b). This profile shows that the gradient

observed before 2016 is not seen by the GNSS velocity field com-
puted after 2017.6 (Fig. 10). This change cannot be explained by
the postseismic deformation of the Pedernales earthquake itself, as
it is spatially smooth and could not account for a velocity change of
several mm yr⁻¹ over a distance as short as 40 km. More precisely,
using our interpolated differential East velocity field (Fig. S10a),
we estimate as 1.7 mm yr⁻¹ the contribution of the large-scale post-seismic effect to the decrease of relative velocity between QUEM and SALF, which stands for only half of the observed variation.

A tentative explanation for this change through time could be the stress perturbation induced by the Pedernales earthquake and its postseismic deformation: both induced extensional strain, resulting in a decrease of the shear stress to normal stress ratio along the Quito thrust. That change could explain a decrease of the creep rate or relocking of the fault plan. This hypothesis would be consistent with slow-slip modulations following an earthquake documented in different tectonic contexts (e.g. Lienkaemper et al. 1997; Schmidt & Bürgmann 2008; Copley & Jolivet 2016).

The change of behavior of the Quito fault after the Pedernales earthquake also questions the surprisingly shallow locking observed from InSAR for the other fault segments in the area. We showed in Section 6 that our GNSS and InSAR dataset implies that slip rates in Ibarra and Peltetec areas are higher on the 2017–2023 period than the ones inferred from block modeling based on GNSS observations on the pre-Pedernales earthquake period (Jarrin et al. 2023). This suggests that part of the current (2017–2023) strain might actually be transient. Monitoring those fault system over long periods might provide answer to this question in future.

## 8. CONCLUSIONS

This study demonstrates that an optimized processing strategy of Sentinel-1 SAR data can produce a robust regional scale surface velocity map in a slowly deforming equatorial region, despite turbulent troposphere, large topographic gradients and abundant vegetation. Comparison with available GNSS data indicate a 2 mm yr⁻¹ consistency in LOS, similar to the consistency among InSAR overlapping tracks. Taking advantage of the GNSS data available in the Ecuadorian and southern Colombian Andes, we could derive a velocity field for the East component, expressed in a Stable South-America reference-frame suitable for tectonic analysis. The obtained continuous velocity field confirms the known overall northeastward translational motion of the NAS. In Ecuador, north of 2° S, active N–S trending faults documented along the eastern front of the Andean Cordillera and along the Western Cordillera, have been interpreted as deformation of a restraining bend. Although the loss of coherence along the eastern flank of the Andean Cordillera prevents constraining the deformation pattern there, the high spatial resolution provided by InSAR shed light on the active internal deformation pattern within the sliver and especially along the western Cordillera. Higher strain rates appear to be primarily localized along active faults in the Western Cordillera and within the Interandean Valley. The velocity gradient varies from north to south, reflecting along strike variation of the locking depth. Two areas show large strain rates reaching 300 nstrain yr⁻¹ equivalent to 6 mm yr⁻¹ of shortening accommodated over a distance of 20 km, consistent with shallow (<5 km) locking depth. Interestingly, the two areas of high strain rate correlate with the location proposed for two large historical earthquakes. In contrast, the Quito fault, which showed evidence of shallow creep during the 1993–2015 period, now exhibits a lower velocity gradient during the 2017–2023 period analyzed here. This result suggests that stress changes induced by the 2016 Mw 7.8 megathrust earthquake and its postseismic deformation have changed the behavior of crustal faults. InSAR data thus reveals both spatial and temporal complexity of the tectonics of this transpressive continental boundary, providing new information for seismic hazard assessment in this highly populated area of the Andes.

## DATA AVAILABILITY

Sentinel-1 data have been downloaded from the PEPS website: peps.cnes.fr. Copernicus DEMs have been downloaded from space-data.copernicus.eu. ERA5 climate reanalysis has been downloaded from the ECMWF website: ecmwf.int. Ecosystems classification map for Colombia has been accessed from https://www.icde.gov.co/geoservicio/iavh-proyecto-pacbao-ecopetrol. Seismic catalogs from the Instituto Geofísico de la Escuela Politécnica Nacional have been downloaded from https://www.igepn.edu.ec/catalogos-sismicos/formulario-catalogos-sismicos (last accessed in September 2023). InSAR and GNSS data provided in this study are available at https://zenodo.org/records/10262497.

## ACKNOWLEDGMENTS

This research was supported by the Centre National d'Études Spatiales (CNES) through the MONA project, the IRD International LMI SVAN project, and the ENS de Lyon (through the L. Marconato PhD funding). All the computations presented in this article were performed using the GRICAD infrastructure (https://gricad.univ-grenoble-alpes.fr), which is supported by Grenoble research communities. The acquisition of GPS data and their analysis in this study were funded by the ANR (Agence Nationale de la Recherche) under the S5 project, contract number ANR-19-CE31-0003-01. We thank the Editor and the two reviewers Franck Audeارد and Jeffrey Freymueller for their constructive and insightful comments.

## SUPPORTING INFORMATION

Supplementary data are available at *GJIRAS* online.

**Figure S1.** Interferogram networks used to compute the time-series. (a) Track A018. (b) Track A120. (c) Track D142 (South segment). (d) Track D142 (North segment). The dashed red line shows the date of the 2016 Mw 7.8 Pedernales earthquake.

**Figure S2.** Histogram of the 2-looks temporal coherence for track A018 (red curve) and weights attributed to each pixel as a function of ω. C is the colinearity.

**Figure S3.** Number of interferogram used in the time-series inversion for each pixel for the track A018. In the cordilleras, all interferograms are unwrapped, while in lowlands only short temporal baseline interferograms are unwrapped (< 1 month).

**Figure S4.** Diagram showing the workflow used to postprocess InSAR time-series, using two GNSS velocity fields, one extracted from the period 1 (pre-2016) and one extracted from the period 3 (post-2017.6).

**Figure S5.** RMS residual of a linear fit of the East component of GNSS time-series, during the postseismic period of the Pedernales earthquake (grey lines). We vary the start date *tmin* of the time-span on which we compute the fit, for 57 stations in the area of our InSAR dataset. The average and median for all stations are shown in red and blue, respectively. The dashed line shows the average RMS for a linear fit on the period 2012-2016, i.e. before the earthquake, and represents the amount of noise in the time-series. The arrow shows the date (July 2017) after which the postseismic deformation is found linear, i.e. the beginning of period 3 (see the main text for details).

**Figure S6.** GNSS velocities projected into line of sight using the ENU local LOS vectors. (a) Track A018. (b) Track A120. (c) Track
D142 (South segment). (c) Track D142 (North segment). Filled circles: GNSS data points. Maps: smoothed and interpolated GNSS velocity maps. Velocities are computed over period 3 (starting in July 2017).

**Figure S7.** Comparison between GNSS and InSAR velocity after referencing. (a) Track A018. (b) Track A120. (c) Track D142 (South segment). (c) Track D142 (North segment). Both velocities are computed on period 3. GNSS velocities are projected into the line of sight of each track using the 3 ENU components. InSAR velocity is taken by averaging velocity values in a radius of 12 pixels (~1200 m) around each GNSS station. RMS of the difference between GNSS and InSAR velocities is shown for each track in the lower-right corner.

**Figure S8.** LOS velocity maps tied to GNSS data for (a) track A018, (b) track A120, (c) track D142 (south segment), (c) track D142 (north segment). The dots show the GNSS velocity, projected into the line of sight, with the same colorscale as the InSAR velocity field, shown in background. Both InSAR and GNSS velocities are computed over period 3.

**Figure S9.** East (a) and vertical (b) InSAR velocity maps before correction of the interpolated large-scale GNSS velocities reflecting the postseismic effect. The dots show the GNSS velocity fields extracted on the same period as the InSAR velocity (period 3), with the same colorscale. Dashed lines are the block boundaries from Jarrin *et al.* (2023).

**Figure S10.** Smoothed and interpolated postseismic GNSS velocity field for the east (a) and up (b) components. The arrows show the horizontal velocity field associated with the postseismic effect of the Pedernales earthquake (period 3–period 1, see the main text for details). The circles show the interpolated component: east (a) or up (b) velocity components, with the same colorscale as the smoothed and interpolated velocity, in background. Dashed lines are the block boundaries from Jarrin *et al.* (2023). The GNSS velocity field is highly smoothed to retain only the large-scale postseismic effect originating from distant sources (postseismic relaxation at depth and/or after slip on the subduction interface).

**Figure S11.** Vertical velocity map and swath profiles. (a) Vertical InSAR velocity map corrected from the interpolated large-scale GNSS velocities reflecting the postseismic effect. The circles show the Up component of the GNSS velocity field computed over period 1 (pre-Pedernales). Boundaries and names of the blocks are shown by dashed black lines and black labels (see Fig. 1 for abbreviations). Faults names are written in white. (b) Swath-profiles numbered from 1 to 12, whose footprint is shown by black rectangles in panel a (same profiles as in Figs 5 and 6). The point clouds and the large dots represent InSAR and GNSS data respectively, both coded with same colorscale as in panel (a).

**Figure S12.** Posterior probability density functions (PDFs) for the Peltetec-Billecocha high-strain area (see Fig. 7). PDFs using GPS data only or InSAR data only are shown by blue and orange unfilled histograms. PDF using GPS + InSAR data are shown by grey filling. A vertical dashed line indicates the retained value for each parameters. See the main text and Fig. 7 for the definition of parameters.

**Figure S13.** Same as Fig. S12 for a profile crossing the Ibarra high-strain area (see Fig. 8).

**Figure S14.** East component of relative time-series between stations SALF and QUEM. Relative velocities are computed on period 1 (in blue) and period 3 (in yellow). The dashed red line shows the date of the 2016 M<sub>w</sub> 7.8 Pedernales earthquake.

**Table S1.** Parameters used for the different GNSS data interpolation steps. See Fig. S4 for the context of each interpolation.

Please note: Oxford University Press is not responsible for the content or functionality of any supporting materials supplied by the authors. Any queries (other than missing material) should be directed to the corresponding author for the paper.

## REFERENCES

Altamimi, Z., Rebischung, P., Collilieux, X., Métivier, L. & Chanard, K., 2023. ITRF2020: an augmented reference frame refining the modeling of nonlinear station motions, *J. Geod.*, **97**, 47. https://doi.org/10.1007/s00190-023-01738-w

Alvarado, A. et al. 2014. Active tectonics in Quito, Ecuador, assessed by geomorphological studies, GPS data, and crustal seismicity, *Tectonics*, **33**, 67–83. https://doi.org/10.1002/2012TC003224

Alvarado, A. et al. 2016. Partitioning of oblique convergence in the northern Andes subduction zone: migration history and the present-day boundary of the North Andean Sliver in Ecuador, *Tectonics*, **35**, 1048–1065. https://doi.org/10.1002/2016TC004117

Alvarado, A. et al. 2018. Seismic, volcanic, and geodetic networks in Ecuador: building capacity for monitoring and research, *Seismol. Res. Lett.*, **89**, 432–439. https://doi.org/10.1785/0220170229

Ansari, H., De Zan, F. & Parizzi, A., 2021. Study of systematic bias in measuring surface deformation with SAR interferometry, *IEEE Trans. Geosci. Remote Sensing*, **59**, 1285–1301. https://doi.org/10.1109/TGRS.2020.3003421

Audemard, F.A., 2003. Geomorphic and geologic evidence of ongoing uplift and deformation in the Mérida Andes, Venezuela, *Quat. Int.*, **101–102**, 43–65. https://doi.org/10.1016/S1040-6182(02)00128-3

Audemard, F.A., Mora-Páez, H. & Fonseca P, H.A., 2021. Net right-lateral slip of the Eastern Frontal Fault System, North Andes Sliver, northwestern South America, *J. South Amer. Earth Sci.*, **109**, 103286. https://doi.org/10.1016/j.jsames.2021.103286

Baize, S. et al. 2015. Paleoseismology and tectonic geomorphology of the Pallatanga fault (Central Ecuador), a major structure of the South-American crust, *Geomorphology*, **237**, 14–28. https://doi.org/10.1016/j.geomorph.2014.02.030

Baize, S. et al. 2020. Active tectonics and earthquake geology along the Pallatanga Fault, *Central Andes of Ecuador. Front. Earth Sci.*, **8**, 193. https://doi.org/10.3389/feart.2020.00193

Beauval, C. et al. 2013. An earthquake catalog for seismic hazard assessment in Ecuador, *Bull. Seismol. Soc. Am.*, **103**, 773–786. https://doi.org/10.1785/0120120270

Beauval, C. et al. 2018. A new seismic hazard model for Ecuador, *Bull. Seismol. Soc. Am.*, **108**, 1443–1464. https://doi.org/10.1785/0120170259

Beauval, C., Yepes, H., Bakun, W.H., Egred, J., Alvarado, A. & Singaucho, J.-C., 2010. Locations and magnitudes of historical earthquakes in the Sierra of Ecuador (1587-1996), *Geophys. J. Int.*, **181**, 1613–1633. https://doi.org/10.1111/j.1365-246X.2010.04569.x

Biggs, J., Wright, T., Lu, Z. & Parsons, B., 2007. Multi-interferogram method for measuring interseismic deformation: denali Fault, *Alaska. Geophys. J. Int.*, **170**, 1165–1179. https://doi.org/10.1111/j.1365-246X.2007.03415.x

Bürgmann, R. & Dresen, G., 2008. Rheology of the lower crust and upper mantle: evidence from rock mechanics, geodesy, and field observations, *Annu. Rev. Earth Planet. Sci.*, **36**, 531–567. https://doi.org/10.1146/annu rev.earth.36.031207.124326

Champenois, J., Baize, S., Vallee, M., Jomard, H., Alvarado, A., Espin, P., Ekström, G. & Audin, L., 2017. Evidences of surface rupture associated with a low-magnitude (Mw5.0) shallow earthquake in the Ecuadorian Andes, *J. Geophys. Res. Solid Earth*, **122**, 8446–8458. https://doi.org/10.1002/2017JB013928

Champenois, J., Pinel, V., Baize, S., Audin, L., Jomard, H., Hooper, A., Alvarado, A. & Yepes, H., 2014. Large-scale inflation of Tungurahua volcano (Ecuador) revealed by Persistent Scatterers SAR interferometry, *Geophys. Res. Lett.*, **41**, 5821–5828. https://doi.org/10.1002/2014GL060956
Cheaib, A. & Doin, M.-P., 2023. Spatial unmixing of pixels for more accurate displacement time series obtained with a small baseline strategy: application on France, in: *IGARSS 2023–2023 IEEE International Geoscience and Remote Sensing Symposium. Presented at the IGARSS 2023–2023 IEEE International Geoscience and Remote Sensing Symposium*, IEEE, Pasadena, CA, USA, pp. 7875–7877. https://doi.org/10.1109/IGARSS52 108.2023.10282268

Copley, A. & Jolivet, R., 2016. Fault rheology in an aseismic fold-thrust belt (Shahdad, eastern Iran), *JGR Solid Earth*, **121**, 412–431. https://doi.org/ 10.1002/2015JB012431

Daout, S. et al. 2016b. Along-strike variations of the partitioning of conver-
gence across the Haiyuan fault system detected by InSAR, *Geophys. J.
Int.*, **205**, 536–547. https://doi.org/10.1093/gji/ggw028

Daout, S., Barbot, S., Peltzer, G., Doin, M.-P., Liu, Z. & Jolivet, R., 2016a. Constraining the kinematics of metropolitan Los Angeles faults with a slip-partitioning model, *Geophys. Res. Lett.*, **43**, 11,192–11,201. https://doi.org/10.1002/2016GL071061

Daout, S., D'Agostino, N., Pathier, E., Socquet, A., Lavé, J., Doin, M.-P., Riesner, M. & Benedetti, L., 2023. Along-strike variations of strain partitioning within the Apennines determined from large-scale multi-temporal InSAR analysis, *Tectonophysics*, **867**, 230076. https://doi.org/10.1016/j.tecto.2023.230076

Doin, M.-P., Cheiab, A. & Tollard, F., 2023. Strategy used for phase unwrapping in the NSBAS MT-InSAR chain, in: *IGARSS 2023–2023 IEEE International Geoscience and Remote Sensing Symposium. Presented at the IGARSS 2023–2023 IEEE International Geoscience and Remote Sensing Symposium*, IEEE, Pasadena, CA, USA, pp. 8210–8213.
https://doi.org/10.1109/IGARSS52108.2023.10283158

Doin, M.-P., Guillaaso, S., Jolivet, R., Lasserre, C., Lodge, F., Ducret, G. & Grandin, R., 2011. Presentation of the small baseline NSBAS processing chain on a case example: the Etna deformation monitoring from 2003 to 2010 using Envisat data, *Presented at the Proceedings of the Fringe Symposium, ESA SP-697*, Frascati, Italy, pp. 3434–3437.

Doin, M.-P., Lasserre, C., Peltzer, G., Cavalié, O. & Doubre, C., 2009. Corrections of stratified tropospheric delays in SAR interferometry: validation with global atmospheric models, *J. Appl. Geophys.*, **69**, 35–50.
https://doi.org/10.1016/j.jappgeo.2009.03.010

Dumont, J.F., Santana, E., Vilema, W., Pedoja, K., Ordóñez, M., Cruz, M., Jiménez, N. & Zambrano, I., 2005. Morphological and microtectonic analysis of quaternary deformation from Puná and Santa Clara Islands, Gulf of Guayaquil, Ecuador (South America), *Tectonophysics*, **399**, 331–350.
https://doi.org/10.1016/j.tecto.2004.12.029

Ebmeier, S.K. et al. 2016. Shallow earthquake inhibits unrest near Chiles-Cerro Negro volcanoes, Ecuador-Colombian border, *Earth Planet. Sci. Lett.*, **450**, 283–291. https://doi.org/10.1016/j.epsl.2016.06.046

Ego, F., Sóbrier, M., Lavenu, A., Yepes, H. & Egues, A., 1996. Quaternary state of stress in the Northern Andes and the restraining bend model for the Ecuadorian Andes, *Tectonophysics*, **259**, 101–116. https://doi.org/10.1016/j.tecto.1995.09.00075-5

Ekström, G., Nettles, M. & Dziewoński, A.M., 2012. The global CMT project 2004–2010: centroid-moment tensors for 13,017 earthquakes, *Phys. Earth Planet. Inter.*, **200–201**, 1–9. https://doi.org/10.1016/j.pepi.2012.04.002

Espin Bedón, P.A. et al. 2022. Unrest at Cayambe Volcano revealed by SAR imagery and seismic activity after the Pedernales subduction earthquake, Ecuador (2016), *J. Volcanol. Geotherm. Res.*, **428**, 107577. https://doi.org/10.1016/j.jvolgeores.2022.107577

European Space Agency, Airbus, 2022. Copernicus DEM. https://doi.org/10.5270/ESA-c5d3d65

Fiorini, E. & Tibaldi, A., 2012. Quaternary tectonics in the central In-
terandean Valley, Ecuador: fault-propagation folds, transfer faults and
the Cotopaxi Volcano, *Global Planet. Change*, **90–91**, 87–103. https:
//doi.org/10.1016/j.gloplacha.2011.06.002

Fuhrmann, T. & Garthwaite, M.C., 2019. Resolving three-dimensional sur-
face motion with InSAR: constraints from multi-geometry data fusion,
*Remote. Sens.*, **11**, 241. https://doi.org/10.3390/rs11030241

Giesecke, A., Capera, A.A.G., Leschiutta, I., Migliorini, E. & Valverde, L.R., 2009. The CERESIS earthquake catalogue and database of the Andean Region: background, characteristics and examples of use, *Ann. Geophys.*, **47**, 421–435. https://doi.org/10.4401/ag-3310

Grandin, R., 2015. Interferometric processing of SLC Sentinel-1 TOPS data, in: *Proceedings of Fringe 2015: Advances in the Science and Applications of SAR Interferometry and Sentinel-1 InSAR Workshop. Presented at the Fringe2015: Advances in the Science and Applications of SAR Interferometry and Sentinel-1 InSAR Workshop*, European Space Agency. https://doi.org/10.5270/Fringe2015.pp116

Guillier, B., Chatelain, J.-L., Jaillard, É., Yepes, H., Poupinet, G. & Fels, J.-F., 2001. Seismological evidence on the geometry of the Orogenic System in central-northern Ecuador (South America), *Geophys. Res. Lett.*, **28**, 3749–3752. https://doi.org/10.1029/2001GL013257

Harrichhausen, N. et al. 2023a. Fault source models show slip rates measured across the width of the entire Fault zone best represent the observed seismicity of the Pallatanga-Puna Fault, Ecuador, *Seismol. Res. Lett.*, **95**, 95–112. https://doi.org/10.1785/0220230217

Harrichhausen, N. et al. 2023b. Distributed right-lateral faults accommodating strain at the northern boundary of the Quito-Latacunga microblock of the Northern Andean Sliver, https://doi.org/10.5194/egusphere-egu2 3-7834

Ho Tong Minh, D., Hanssen, R., Doin, M. & Pathier, E., 2022. Advanced methods for time-series InSAR, in: *Surface Displacement Measurement from Remote Sensing Images*. Wiley Online Library, pp. 125–153. https://doi.org/10.1002/9781119986843.ch5

Jaillard, E., 2022. Late cretaceous-paleogene orogenic build-up of the Ecuadorian Andes: review and discussion, *Earth Sci. Rev.*, **230**, 104033. https://doi.org/10.1016/j.earscirev.2022.104033

Jarrin, P. et al. 2023. Continental block motion in the Northern Andes from GPS measurements, *Geophys. J. Int.*, **235**, 1434–1464. https://doi.org/10.1093/gji/ggad294

Jolivet, R., Grandin, R., Lasserre, C., Doin, M.-P. & Peltzer, G., 2011. Systematic InSAR tropospheric phase delay corrections from global meteorological reanalysis data, *Geophys. Res. Lett.*, **38**, L17311. https://doi.org/10.1029/2011GL048757

Jomard, H. et al. 2021. Interactions between active tectonics and gravitational deformation along the Billecocha fault system (Northern Ecuador): insights from morphological and paleoseismological investigations, *J. South Amer. Earth Sci.*, **111**, 103406. https://doi.org/10.1016/j.jsames.2021.103406

Kanamori, H. & McNally, K.C., 1982. Variable rupture mode of the subduction zone along the Ecuador-Colombia coast, *Bull. Seismol. Soc. Am.*, **72**, 1241–1253.

Kellogg, J.N. & Bonini, W.E., 1982. Subduction of the Caribbean plate and basement uplifts in the overriding South American plate, *Tectonics*, **1**, 251–276. https://doi.org/10.1029/TC001i003p00251

Lavenu, A., Winter, T. & Dávila, F., 1995. A pliocene-quaternary compres-
sional basin in the Interandean Depression, *Central Ecuador. Geophys. J.
Int.*, **121**, 279–300. https://doi.org/10.1111/j.1365-246X.1995.tb03527.x

Lemrabet, L., Doin, M., Lasserre, C. & Durand, P., 2023. Referencing of continental-scale InSAR-derived velocity fields: case study of the eastern Tibetan plateau, *JGR Solid Earth*, **128**, e2022JB026251. https://doi.org/10.1029/2022JB026251

Lienkaemper, J.J., Galehouse, J.S. & Simpson, R.W., 1997. Creep response of the Hayward fault to stress changes caused by the Loma Prieta earthquake, *Science*, **276**, 2014–2016. https://doi.org/10.1126/science.276.53 21.2014

Litherland, M. & Aspden, J.A., 1992. Terrane-boundary reactivation: a con-
trol on the evolution of the Northern Andes, *J. South Amer. Earth Sci.*, **5**,
71–76. https://doi.org/10.1016/0895-9811(92)90060-C

López-Quiroz, P., Doin, M.-P., Tupin, F., Briole, P. & Nicolas, J.-M., 2009. Time series analysis of Mexico City subsidence constrained by radar interferometry, *J. Appl. Geophys.*, **69**, 1–15. https://doi.org/10.1016/j.ja ppgeo.2009.02.006

Maghsoudi, Y., Hooper, A.J., Wright, T.J., Lazecky, M. & Ansari, H., 2022. Characterizing and correcting phase biases in short-term, multilooked
interferograms, *Remote Sens. Environ.*, **275**, 113022. https://doi.org/10.1016/j.rse.2022.113022

Marconato, L., Doin, M.-P., Audin, L. & Pathier, E., 2024. Ionospheric compensation in L-band InSAR time-series: performance evaluation for slow deformation contexts in equatorial regions, *Sci. Remote Sensing*, **9**, 100113. https://doi.org/10.1016/j.srs.2023.100113

Mariniere, J., Nocquet, J.-M., Beauval, C., Champenois, J., Audin, L., Alvarado, A., Baize, S. & Socquet, A., 2020. Geodetic evidence for shallow creep along the Quito fault, Ecuador, *Geophys. J. Int.*, **220**, 2039–2055. https://doi.org/10.1093/gji/ggz564

Mathey, M., Doin, M.-P., André, P., Walpersdorf, A., Baize, S. & Sue, C., 2022. Spatial heterogeneity of uplift pattern in the western European alps revealed by InSAR time-series analysis, *Geophys. Res. Lett.*, **49**, e2021GL095744. https://doi.org/10.1029/2021GL095744

Maubant, L., Frank, W.B., Wallace, L.M., Williams, C.A. & Hamling, I., 2023. Imaging the spatiotemporal evolution of plate coupling with interferometric radar (InSAR) in the Hikurangi Subduction Zone, *Geophys. Res. Lett.*, **50**, e2023GL105388. https://doi.org/10.1029/2023GL105388

Mora-Páez, H. et al. 2019. Crustal deformation in the northern Andes—a new GPS velocity field, *J. South Amer. Earth Sci.*, **89**, 76–91. https://doi.org/10.1016/j.jsames.2018.11.002

Mora-Páez, H., Mencin, D.J., Molnar, P., Diederix, H., Cardona-Piedrahita, L., Peláez-Gaviria, J. & Corchuelo-Cuervo, Y., 2016. GPS velocities and the construction of the Eastern Cordillera of the Colombian Andes, *Geo-phys. Res. Lett.*, **43**, 8407–8416. https://doi.org/10.1002/2016GL069795

Mora-Páez, H., Peláez-Gaviria, J., Diederix, H., Bohórquez-Orozco, O., Cardona-Piedrahita, L., Corchuelo-Cuervo, Y., Ramírez-Cadena, J. & Díaz-Mila, F., 2018. Space geodesy infrastructure in Colombia for geodynamics research, *Seismol. Res. Lett.*, **89**, 446–451. https://doi.org/10.1785/0220170185

Mothes, P.A., Nocquet, J.-M. & Jarrín, P., 2013. Continuous GPS network operating throughout Ecuador, *Eos Trans. AGU*, **94**, 229–231. https://doi.org/10.1002/2013EO260002

Nocquet, J.-M. et al. 2014. Motion of continental slivers and creeping subduction in the northern Andes, *Nature Geosci*, **7**, 287–291. https://doi.org/10.1038/ngeo2099

Nocquet, J.-M. et al. 2017. Supercycle at the Ecuadorian subduction zone revealed after the 2016 Pedernales earthquake, *Nature Geosci*, **10**, 145–149. https://doi.org/10.1038/ngeo2864

Paris, G., Machette, M.N., Dart, R.L. & Haller, K.M., 2000. Map and database of Quaternary faults and folds in Colombia and its offshore regions (Report No. 2000–284), Open-File Report. https://doi.org/10.3133/ofr00284

Pennington, W.D., 1981. Subduction of the Eastern Panama Basin and seismotectonics of northwestern South America, *J. Geophys. Res.*, **86**, 10753–10770. https://doi.org/10.1029/JB086iB11p10753

Pinel-Puyssegur, B., Michel, R. & Avouac, J.-P., 2012. Multi-link InSAR Time Series: enhancement of a wrapped interferometric database, *IEEE*

J. Sel. Top. Appl. Earth Observations Remote Sensing, **5**, 784–794. https://doi.org/10.1109/JSTARS.2012.2196758

Rolandone, F. et al. 2018. Areas prone to slow slip events impede earthquake rupture propagation and promote afterslip, *Sci. Adv.*, **4**, eaao6596. https://doi.org/10.1126/sciadv.aao6596

Schmidt, D.A. & Bürgmann, R., 2008. Predicted reversal and recovery of surface creep on the Hayward fault following the 1906 San Francisco earthquake, *Geophys. Res. Lett.*, **35**, L19305. https://doi.org/10.1029/2008GL035270

Thollard, F. et al. 2021. FLATSIM: the ForM@Ter LArge-scale multi-temporal Sentinel-1 InterferoMetricy service, *Remote. Sens.*, **13**, 3734. https://doi.org/10.3390/rs13183734

Tian, Z., Freymueller, J.T., Yang, Z., Li, Z. & Sun, H., 2023. Frictional properties and rheological structure at the Ecuadorian subduction zone revealed by the postseismic deformation due to the 2016 Mw7.8 Pedernales (Ecuador) earthquake, *JGR Solid Earth*, **128**, e2022JB025043. https://doi.org/10.1029/2022JB025043

Tibaldi, A., Rovida, A. & Corazzato, C., 2007. Late quaternary kinematics, slip-rate and segmentation of a major Cordillera-parallel transcurrent fault: the Cayambe-Afiladores-Sibundoy system, NW South America, *J. Struct. Geol.*, **29**, 664–680. https://doi.org/10.1016/j.jsg.2006.11.008

Trenkamp, R., Kellogg, J.N., Freymueller, J.T. & Mora, H.P., 2002. Wide plate margin deformation, southern Central America and northwestern South America, CASA GPS observations, *J. South Amer. Earth Sci.*, **15**, 157–171. https://doi.org/10.1016/S0895-9811(02)00018-4

Vaca, S., Vallée, M., Nocquet, J.-M. & Alvarado, A., 2019. Active deformation in Ecuador enlightened by a new waveform-based catalog of earthquake focal mechanisms, *J. South Amer. Earth Sci.*, **93**, 449–461. https://doi.org/10.1016/j.jsames.2019.05.017

Vallée, M. et al. 2023. Self-reactivated rupture during the 2019 M = 8 northern Peru intraslab earthquake, *Earth Planet. Sci. Lett.*, **601**, 117886. https://doi.org/10.1016/j.epsl.2022.117886

Velandia, F., Acosta, J., Terraza, R. & Villegas, H., 2005. The current tectonic motion of the Northern Andes along the Algeciras Fault System in SW Colombia, *Tectonophysics*, **399**, 313–329. https://doi.org/10.1016/j.tecto.2004.12.028

Villegas-Lanza, J.C., Chlieh, M., Cavalié, O., Tavera, H., Baby, P., Chire-Chira, J. & Nocquet, J.-M., 2016. Active tectonics of Peru: heterogeneous interseismic coupling along the Nazca megathrust, rigid motion of the Peruvian sliver, and Subandean shortening accommodation: active tectonics of Peru, *J. Geophys. Res. Solid Earth*, **121**, 7371–7394. https://doi.org/10.1002/2016JB013080

Winter, T., Avouac, J.-P. & Lavenu, A., 1993. Late quaternary kinematics of the Pallatanga strike-slip fault (Central Ecuador) from topographic measurements of displaced morphological features, *Geophys. J. Int.*, **115**, 905–920. https://doi.org/10.1111/j.1365-246X.1993.tb01500.x