# Satellite imaging of the 2015 M7.2 earthquake in the Central Pamir, Tajikistan, elucidates a sequence of shallow strike-slip ruptures of the Sarez-Karakul fault

A Elliott, J Elliott, J Hollingsworth, G Kulikova, B Parsons, R Walker

► To cite this version:

A Elliott, J Elliott, J Hollingsworth, G Kulikova, B Parsons, et al.. Satellite imaging of the 2015 M7.2 earthquake in the Central Pamir, Tajikistan, elucidates a sequence of shallow strike-slip ruptures of the Sarez-Karakul fault. Geophysical Journal International, 2020, 221 (3), pp.1696-1718. 10.1093/gji/ggaa090. hal-03062254

HAL Id: hal-03062254

https://hal.science/hal-03062254

Submitted on 14 Oct 2021

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
<figure><img src="image_2.png" /></figure>

# Satellite imaging of the 2015 M7.2 earthquake in the Central Pamir, Tajikistan, elucidates a sequence of shallow strike-slip ruptures of the Sarez-Karakul fault

A. Elliott<sup>1,5,*</sup> J. Elliott<sup>2,5</sup> J. Hollingsworth<sup>3</sup> G. Kulikova<sup>4</sup> B. Parsons<sup>1,5</sup> and R. Walker<sup>1,5</sup>

Accepted 2020 February 14. Received 2020 February 7; in original form 2019 May 10

## SUMMARY

On 7 December 2015, a shallow $M_w$ 7.2 strike-slip earthquake struck the Murghab River Valley in the Central Pamirs of Tajikistan. Seismologically this event was similar to a large seismic event in 1911 whose causative fault has never been identified. We measure the displacement field of the 2015 event from satellite observations using Sentinel-1 radar interferometry, Landsat-8 optical pixel-tracking, and surface rupture mapping from high resolution SPOT-6/7 imagery to characterize the role this earthquake rupture plays in the accommodation of strain on its causative structure, the Sarez-Karakul fault. We present geomorphic mapping and interpretations of other Quaternary-active reaches of this fault system, which highlight variable rupture history of the different sections. These sections appear to be separated by inherited bedrock structural boundaries. Significantly, the reaches of the fault northeast and southwest of the 2015 rupture exhibit the freshest morphology prior to 2015, indicative of a more recent rupture than elsewhere. Using new high resolution imagery we map fresh scarps at the northern and southern ends of the Sarez-Karakul fault which may represent this 1911 rupture. To test which of these reaches could have been the source of the elusive 1911 event, we compare synthetic seismograms from three plausible fault sources determined from geomorphology, with observed seismic traces from 1911 at early recording stations throughout Europe. We find that the best fitting fault source is in fact southwest of the 2015 rupture, meaning that we have a record of three distinct recent events on the Sarez-Karakul fault system—two of them instrumentally recorded. Our mapping of these separate events reveals a correlation between their boundaries and the active and inherited thrust and suture systems that intersect the northeast striking left-lateral fault, suggesting structural control over the extents of individual ruptures on the active strike-slip fault.

**Key words:** earthquake source observations; seismic cycle; satellite geodesy; continental neotectonics; dynamics and mechanics of faulting; palaeoseismology.

## 1 INTRODUCTION

The $M_w$ 7.2 Murghob, Tajikistan earthquake of 7 December 2015, occurred along a particularly high-elevation mountainous reach of the left-lateral Sarez-Karakul fault (SKF), a distinctive northeast-striking structure that cuts across several of the large east-west

convergent fault zones and tectonic sutures of the North-Central Pamir mountains (Fig. 1) (Metzger et al. 2017; Sangha et al. 2017; Gan et al. 2018). The 2015 event is similar in location and kinematics to a seismic event in 1911 whose causative fault has been debated and has evaded recognition despite repeated investigations spanning over a century (Klotz & Galitzin 1915; Jeffreys 1923; Oldham 1923; Macelwane 1926; Ambraseys & Bilham 2012; Kulikova et al. 2016). Remote sensing of the 2015 ground displacements using optical and radar imaging allows us to measure the fault slip

<footnote n="1">¹Department of Earth Sciences, University of Oxford, UK. Email: ajelliott@usgs.gov</footnote>

<footnote n="2">²School of Earth and Environment, University of Leeds, UK</footnote>

<footnote n="3">³ISTerre, Université Grenoble-Alpes, France</footnote>

<footnote n="4">⁴Institute of Geosciences, University of Potsdam, Germany</footnote>

<footnote n="5">⁵Centre for the Observation and Modelling of Earthquakes, Volcanoes and Tectonics (COMET), UK</footnote>

<footnote>* Now at: the U.S. Geological Survey Earthquake Science Center, Moffett Field, CA, USA.</footnote>
<figure><img src="image_3.png" /><figcaption>**Figure 1.** Setting of the 2015 and 1911 earthquakes in the Central Pamir, showing topography, seismicity (circles), active faults (thin black lines) and coverage (white boxes) of satellite data sets acquired to image deformation due to the 7 December 2015 event. Inset shows location of figure centred on Tajikistan. Earthquake hypocentres from the Advanced National Seismic System (ANSS) catalogue 1933–2018. Focal mechanisms of the 1911 and 2015 events from Kulikova *et al.* (2016). Faults drawn from HimaTibet database of Taylor & Yin (2009) and Styron *et al.* (2010). MF, Muji Fault; PFT, Pamir Frontal Thrust; SKF, Sarez-Karakul Fault; AV, Alai Valley; AM, Alai Mountains (Range). GPS vectors and $1\sigma$ uncertainty ellipses relative to stable Eurasia from Ischuk *et al.* (2013).</figcaption></figure>

distribution during the 7 December 2015 event and map in detail the previously unreported surface ruptures of this event. Our observations reveal comprehensively the extent of discrete surface ruptures versus blind faulting, and further mapping of the system shows the relationship of coseismic slip with other ruptures along strike as well as with regional bedrock geology.

Documentation of sequential ruptures along a fault system can be used to address questions of earthquake nucleation, termination and triggering; structural and topographic development; and strain accumulation and moment release. Having a record of multiple ruptures along a single fault affords a rare opportunity to understand the processes by which a whole fault system accommodates crustal strain through several earthquake cycles, linking tectonic structure and landscape morphology to seismic hazard (e.g. Stein *et al.* 1997; Duross *et al.* 2016; Schwartz 2018). As geometric irregularities along faults are both theorized and observed to control rupture propagation and earthquake ground motions, understanding rupture behaviour requires investigating the structures and processes involved in the intersection zones between rupture tips, where one earthquake may end and another begin (Sibson 1985; Wesnousky 2006; Lozos *et al.* 2011, 2013; Ong *et al.* 2019). We require

observations of multiple, varied ruptures along a particular fault system to begin to formulate a picture of the representative suite of events that may occur in different systems (Sharp 1982; Sharp *et al.* 1989; Crone & Haller 1991; Depolo *et al.* 1991; Stein *et al.* 1997; Berberian *et al.* 1999; Sibson *et al.* 2011; Walker *et al.* 2011; Philibosian *et al.* 2014; Duross *et al.* 2016; Bilham 2019).

In this paper, we present the first detailed map of surface ruptures from the 2015 Murghob earthquake along the SKF, as well as a near-fault coseismic deformation field from Landsat-8 image pixel-tracking that we use in conjunction with interferometric synthetic aperture radar (InSAR) line-of-sight (LOS) displacements to model slip at depth along the ~80-km-long fault rupture. Using these near-field observations, we construct a coseismic slip model with significant detail that allows us to interrogate the role of complexities along the fault which appear to play a role in modulating its slip history. We then present geomorphic evidence for prior surface ruptures along the reaches of the SKF northeast and southwest of the 2015 event.

To test whether these northern or southern reaches of the SKF could have been the source of the elusive 1911 event (e.g. Fig. 2),
<figure><img src="image_4.png" /><figcaption>**Figure 2.** Geomorphic and structural context of the Sarez-Karakul fault ruptures showing Cenozoic faults and Quaternary-active features. The 1911 event focal mechanism and location, in yellow, are from Kulikova *et al.* (2016). Red lines represent newly mapped 2015 surface ruptures reported in this study. SKF trace mapped from Quaternary scarps in 1.5-m SPOT-6/7 imagery and occasionally higher resolutions available in Google Earth. Queried fault traces come from previously published figures in, for example Stübner *et al.* (2013). Other structures from Schwab *et al.* (2004), Stübner *et al.* (2013) and Rutte *et al.* (2017b) with vergence and orientation preference given to most recent and most comprehensive cross sectional interpretation, Rutte *et al.* (2017b). White boxes indicate footprints of subsequent figures in this paper.</figcaption></figure>
we expand the analysis of analog seismograms conducted by Kulikova et al. (2016). From our mapping of coseismic displacement and regional fault scarps from satellite imaging, we derive four candidate sources for the 1911 Sarez earthquake and test the misfit of synthetic seismograms from each of these candidate sources against the observed waveforms analysed by Kulikova et al. (2016).

Our mapping of rupture extents on the SKF, coupled with revised analogue seismological analysis, indicates a history of rupture of three distinct reaches of this fault system. Relating the mapped extents and bounds of these respective ruptures to surrounding tectonic structures, we suggest that the regional geometry of suture zones including contemporary thrust faults and inherited extensional faults, may control the extent of coseismic ruptures along the Sarez-Karakul fault.

## 2 BACKGROUND: SURFACE FAULTING IN THE PAMIRS

Active shallow faulting in the Pamirs is broadly distributed among a set of conjugate strike-slip faults which cross-cut and/or merge with arcuate shear zones that are inherited from an alternating history of convergence and extension (e.g. Searle & Hacker 2018). Naturally this network of inherited structures may modulate where and how active faulting occurs, making it important to understand the relationships between contemporary earthquakes and regional tectonic structures.

### 2.1 Tectonic setting

The Murghob earthquake occurred along the transtensional Sarez-Karakul fault zone (Figs 1 and 2), separating the northward-convergent eastern Pamir from the western Pamir. The south-southwest-striking SKF stands out among the largely east-west tectonic structures of the Central Pamir, which predominantly reflect a long history of north-south-directed crustal shortening and thickening (e.g. Burtman & Molnar 1993; Strecker et al. 1995; Schwab et al. 2004; Stübner et al. 2013; Rutte et al. 2017a). The characteristic shape of the Pamir is defined by a series of arcuate suture zones (Fig. 1) that correspond to the major terrane-accretion sutures of Tibet (Burtman & Molnar 1993; Schwab et al. 2004; Yin 2010; Searle & Hacker 2018).

This zone of extreme elevations (most of the Pamir exceeds 4 km above sea level) reflects crustal thickening in a rare zone of continental subduction (e.g. Billington et al. 1977; Burtman & Molnar 1993; Kufner et al. 2016, 2017). Beneath the Pamir and Hindu Kush, a subvertical curtain of earthquake hypocentres extending to depths of ~250 km illuminates the subducted lithospheric slabs of both Asian and Indian crust (Billington et al. 1977; Pegler & Das 1998; Sippl et al. 2013; Kufner et al. 2016, 2017). Tomographic imaging further outlines the structure of the oppositely vergent subducted slabs (Koulakov & Sobolev 2006; Negredo et al. 2007; Schneider et al. 2013). These data sets, combined with modelling of lithospheric dynamics, suggest that this continental subduction is attended by lithospheric delamination, as well as stagnation and break-off of the slabs (Negredo et al. 2007; Kufner et al. 2017). Collectively these processes have induced a varied history of convergence and extension in the overlying crust as the Pamir plateau was built.

Based on timing relationships between extensional and convergent faulting in the Central Pamir, Rutte et al. (2017b, a) conclude that the alternation between crustal thickening and exhumation of deep crustal metamorphic rocks reflects the break-off of the Indian

slab followed by the delamination of Asian mantle lithosphere upon arrival of the Indian craton at the convergent boundary. This alternation of extensional and compressional tectonics has left a series of arcuate sutures and shear zones that trend roughly east-west through the central Pamir (rotating around either side of the apparent indenter; Fig. 1). Several of these have a long tectonic history involving repeated reactivation of continental sutures (e.g. the extensive Kunlun and Jinsha—locally called Tanymas—sutures; Robinson et al. 2012; Robinson 2015; Jepson et al. 2018). Some remain active today such as the Sarez-Murghab thrust system along the Pshart suture, and the Aksu-Murghab shear zone (Fig. 2), which are expressed with gouge zones and scarps through Quaternary deposits, as well as contemporary shallow seismicity (Verma & Chandra Sekhar 1985; Strecker et al. 1995; Pegler & Das 1998; Zhang et al. 2011). Contemporary conjugate strike-slip fault systems, such as the Sarez-Karakul fault and the Aksu-Murghab and Aksu-Rangkul fault systems (Fig. 2), occupy and cut across these convergent shear zones, accommodating east-west deformation of the Pamir as convergence propagates westward into the Afghan-Tajik Basin (Pozzi & Feinberg 1991; Strecker et al. 1995; Bourgeois et al. 1997; Mohadjer et al. 2010).

### 2.2 Modern geodetic deformation

India–Eurasia convergence continues at an average rate of 30–34 mm yr⁻¹ across the Pamirs as measured by the Global Positioning System (GPS, Mohadjer et al. 2010; Zubovich et al. 2010; Ischuk et al. 2013; Schurr et al. 2014). The largest share of this deformation is presently taken up across a single structure at the northern topographic, structural and seismotectonic boundary of the Pamir, the Pamir Frontal Thrust (PFT) in the Alai Valley (Coutand et al. 2002), across which there is crustal convergence of 10–15 mm yr⁻¹ (Mohadjer et al. 2010; Zubovich et al. 2010; Ischuk et al. 2013; Zubovich et al. 2016). The remainder of the shortening is taken up within the rest of the high Northern and Central Pamir, although low gradients in northward GPS velocity demonstrate slow and distributed modern internal shortening of the range (Mohadjer et al. 2010; Ischuk et al. 2013).

To the west, GPS vectors rotate counterclockwise with increasing westward velocities attributed to east-west extension of 5–10 mm yr⁻¹ within the mountain range (Ischuk et al. 2013). This extension is balanced by 10 mm yr⁻¹ of east-west convergence farther west within the Tajik depression, and because it exceeds the internal northward shortening of the range, the resulting strain field is interpreted to arise from gravitational collapse of the high Pamir rather than pure shear (Ischuk et al. 2013; Schurr et al. 2014). Indeed, Jay et al. (2017) calculate that west-directed extension in the Pamir can be entirely explained by gravitational collapse of the high plateau. This east-west extension is in part accomplished via slip on the oblique left-normal Sarez-Karakul fault, which runs between the Alai mountains in the north and the Gunt shear zone in the south (Figs 1 and 2).

### 2.3 Role of the Sarez-Karakul Fault

Both diverging GPS vectors (fig. 1; Mohadjer et al. 2010; Zubovich et al. 2010, 2016; Ischuk et al. 2013) and the large-scale geomorphology (Fig. 2) attest to contemporary activity of the Sarez-Karakul fault. Near its northern intersection with the west-east trending, right-lateral Muji fault (Fig. 1) and other active strike-slip and thrust structures within the Alai mountains (Fig. 1), the SKF
system forms a conspicuous series of horsts and grabens, which con-
tain and shape Karakul ('black lake' in Kyrgyz; Fig. 2). South of the
lake, at the conjunction of the Muzkol and Muksu rivers, the several
active strands of the Karakul graben system come together to form
a localized left-lateral fault that continues south-southwestward
(Fig. 2). This structural change in the fault zone from an array of
oblique left-normal faults to a single, largely left-lateral fault cor-
responds spatially to its intersection with the active Aksu-Rangkul
dextral/thrust fault system (Strecker et al. 1995; Rutte et al. 2017b),
probably accounting for the southward change from transtension to
more purely strike-slip as east-west motion is compensated by slip
on the Aksu-Rangkul fault.

Before 2015 the apparently youngest scarps on the northern SKF
occupied this 12-km reach between the Muzkol River and the north
end of the Kokuibel Valley where the SKF crosses (and appears to
offset) the Tanymas suture system (Strecker et al. 1995). Intersec-
tion with the Tanymas suture system (Rutte et al. 2017b) coincides
with a left-step in the SKF and another change in kinematics back to
a steeper rake reflecting oblique transtension. To the south, Strecker
et al. (1995), Schurr et al. (2014) and Rutte et al. (2017b) reported
young (late Quaternary) normal faulting scarps along the west-
ern front of the Officer range, marking the boundary of the large
Kokuibel valley graben—the most conspicuous large-scale geomor-
phic feature of the northern Pamir after the Karakul grabens, with
2 km of relief between hangingwall and footwall (Fig. 2).

At the southern end of the Kokuibel valley, the SKF takes a right-
ward double-bend as it cuts across the North Sarez Shear Zone, a
substantial extensional shear zone that bounds the gneissic Sarez
Dome to the north (fig. 2; Stübner et al. 2013; Rutte et al. 2017b;
Searle & Hacker 2018). Within the crystalline bedrock of the Sarez
Dome the SKF occupies a rugged mountainous regime with few
Quaternary deposits and landforms to offset. Nonetheless it is ex-
pressed as a series of shutter ridges, perched valleys and left-laterally
offset ridgelines. Largely, however, the active trace is obscured by
steep scree slopes and hummocky rock glaciers or young moraines.
The southern boundary of the Sarez Dome is the Rushan-Pshart
suture, demarked locally by the active Pshart thrust (Schurr et al.
2014), the trace of which coincides with the valley-filling Lake
Sarez and the epicentre of the 2015 earthquake. South of the Pshart
thrust and Lake Sarez, the fault is roughly traceable for a further 90–
100 km as a substantial topographic lineament with discontinuous
exposures of young offsets in rare Quaternary landforms. It can be
traced to around its general intersection with the Gunt shear zone.
Beyond the Gunt shear zone and eponymous river, the SKF trace
becomes unclear and active left-lateral fault traces are absent from
regional tectonic maps (Stübner et al. 2013; Schurr et al. 2014).
Deformation may be absorbed along an array of inherited faults,
or may be transformed westward onto the lateral boundary of the
Pamir along bedrock fault boundaries within the Gunt shear zone.

## 2.4 The enigmatic 1911 Sarez earthquake

On 18 February 1911 a $M_S$ 7.3 seismic event centred in the Pamir was recorded at nascent seismological observatories around the world (Klotz & Galitzin 1915; Kulikova et al. 2016). The earthquake coincided with the largest single non-volcanic landslide in recorded history (600 m high × ~3 km diameter), the Usoi landslide—named for the village it buried—which dammed the Murghab river, eventually impounding 17 cubic km of water upstream in Sarez Lake—named for the village it submerged (Schuster 2005; Ambraseys & Bilham 2012). Early on, the coincidence of the massive

landslide and the observed earthquake drove significant develop-
ment in the physical interpretation of seismic records: a series of
papers debated whether a tectonic earthquake had occurred at all, or
whether the recorded teleseismic signal could be entirely explained
by the impact of the sliding/falling mass on the valley floor (Klotz
& Galitzin 1915; Jeffreys 1923; Oldham 1923; Macelwane 1926).
Although more recent studies have conclusively established that the
seismic event had a double-couple source at finite depth, and there-
fore must have been a strike-slip earthquake (Ambraseys & Bilham
2012; Kulikova et al. 2016), the fault that produced it has never
been determined and thus researchers have only speculated about
the source (Ambraseys & Bilham 2012; Bindi et al. 2014; Schurr
et al. 2014; Kulikova et al. 2016; Metzger et al. 2017). The more
recent 2015 earthquake, with a similar magnitude, epicentre and
focal mechanism, provides an opportunity to revisit the 1911 earth-
quake with a comparative study with the aim of finally identifying
the fault that slipped in the earlier event.

## 2.5 Prior studies of the 2015 Murghob ('Murghab') earthquake

The 7 December 2015 earthquake originated at 07:50:05 UTC (11:50 a.m. local time) beneath Lake Sarez, with an epicentre of 38.2107°N 72.7797°E (U.S. Geological Survey 2017). Sangha et al. (2017) and Metzger et al. (2017) have shown through remote sensing data sets that nearly the entire fault rupture lies north of this epicentre, and of the lake, on the Sarez-Karakul fault. Teleseismic backprojection of high frequency seismic energy conducted by Sangha et al. (2017) shows that the event indeed propagated unilaterally toward the northeast, in agreement with the remote sensing geodetic analyses. Comparing their time-history of seismic radiation with the fault geometry evident in ALOS and Sentinel-1 radar interferograms and Landsat-8 optical pixel offsets, Sangha et al. (2017) conclude that the rupture reached supershear velocities between 10–16 and 20–28 km along the fault, then paused before proceeding substantially more slowly through the restraining double-bend of 22°, then picking up speed again as it propagated further north.

InSAR and pixel-tracking measurements by both Sangha et al. (2017) and Metzger et al. (2017) show 35–40 km of shallow fault rupture, but require an ~80 km length of slip at depth, extending north beneath the Kokuibel valley as blind strike-slip rupture. All fault models for this rupture require at least three fault segments, representing a restraining double-bend with a ~15 km central reach oriented clockwise of the regional fault strike. Sangha et al. (2017) analyse a single pair of Landsat-8 images from October 2015 and October 2016, respectively, to measure a representative slip distribution, which reaches a maximum of 4.1 m just before the bend and drops off northward. Their line-of-sight displacement maps from ALOS and Sentinel-1 radar interferograms require significant slip at depth for a further 40 km north beyond the bend. Metzger et al. (2017) visited the northernmost surface ruptures from the earthquake at 38.472°N, where they measured 2 m lateral offsets as well as a wide array of ground failures associated with shallow ground failures and landsliding. The narrow band of deformation in their geodetic results suggests that shallow slip extends between there and the epicentre, but only the northern 4 km of surface rupture were accessible due to the extreme elevation and relief of the Muzkol range (Fig. 2).

No previous studies have reported the full extent of surface rupture from the 2015 earthquake along the Sarez-Karakul fault. Here
we present a comprehensive map of the discrete tectonics surface breaks visible in 1.5-m scale SPOT-6/7 imagery and in higher resolution (submetre) DigitalGlobe imagery available in Google Earth. We improve the spatial coverage of the optically derived horizontal displacement field by averaging the pixel-tracking offsets between year-spanning Landsat-8 scene pairs, presenting a robust near-field slip distribution derived from the resulting 30-m resolution map of surface displacements. We corroborate our near-field slip distribution by jointly inverting InSAR phase displacement and Landsat-8 pixel-tracking deformation fields to model fault slip at depth. We incorporate new mapping and seismological analysis to place this earthquake in the context of rupture history along the SKF.

## 3 METHODS: IMAGING THE 2015 RUPTURE

We measure the deformation field of the Murghob earthquake rupture using pre- and post-event imagery from a combination of radar and optical satellite sensors. ESA's Sentinel-1A&1B satellites captured radar images with a temporal baseline of 24 days spanning the earthquake, permitting us to construct the most coherent possible interferograms on both descending (southward) and ascending (northward) passes of the satellite. These provide deformation fields resolved in two different satellite line-of-sight directions. The NASA/USGS Landsat-8 optical imaging satellite has a maximum resolution of 15 m per pixel in the panchromatic band and makes a repeat pass every 16 d, permitting detailed tracking of surface changes. We conduct pixel-tracking on pairs of these optical images using the COSI-Corr program developed by Leprince et al. (2007) and Ayoub et al. (2009), mapping horizontal displacements in relatively high resolution without losing coherence in the extreme near-field of the fault rupture. We purchased higher resolution optical SPOT-6/7 images) from CNES (Centre National d'Etudes Spatiales) covering dates before (October 2013, September 2014, September 2015) and after (January 2016) the 7 December 2015 earthquake to identify, map, and measure discrete surface ruptures. Collectively these form a detailed picture of the ground deformation caused by the 2015 earthquake.

### 3.1 Radar interferometry

Sentinel-1 satellite C-band (~5.545 cm wavelength) synthetic aperture radar (SAR) images of the epicentral area were acquired at 24-d intervals spanning the 7 December earthquake, during descending (track 5) passes on 18 November and 12 December 2015, and on ascending (track 100) passes on 6 December and 30 December (Figs 1 and 3). We computed differential interferograms using the GAMMA software package (Wegmüller & Werner 1997) to derive line-of-sight (LOS) ground surface displacement from both descending-track and ascending-track satellite vantages. Detailed processing steps may be found in Fig. S4.

Due to the rugged relief resulting in layover and shadowing in the side-looking SAR images, (as well as snow cover changes and landsliding) coherence is low and unwrapping challenging in the mountainous areas through which the southern half of the fault rupture passes (Figs 2 and 3). Where the nature of unwrapping errors were readily identified, they are corrected with integer two-pi additions, and where this is not possible the ambiguous data are removed (e.g. Fig. 3).

Additionally, high phase gradients due to the large displacements in the near-field mean that interferometric SAR (InSAR) displacements cannot be derived close to the fault. With these effects combined, no InSAR displacement data is available within 10 km on either side of the fault along its southern reaches; rather, the available data provides constraints on faulting in the mid and far-field (10–100 km).

The look direction on the ascending track is more aligned along the strike-direction of the rupture and therefore more sensitive to fault parallel lateral displacements than on the descending track (Fig. 3). Over 50 cm of motion towards the satellite is observed on the northwestern side of the rupture, and over 50 cm of motion away on the southeastern side (Fig. 3), consistent with left-lateral strike-slip faulting.

### 3.2 Optical pixel-tracking with Landsat-8

In the near-field of the fault rupture, high displacements make successive radar images lose coherence, and high InSAR fringes (strong displacement gradients) cause phase unwrapping errors, masking the displacement field where it contains the greatest share of deformation. To image the horizontal displacement field in the near-field of the fault rupture, we perform pixel-tracking on successive 15-m optical Landsat-8 Level 1 scenes to measure subpixel offsets using the COSI-Corr 'Coregistration of Optically Sensed Images and Correlation' software of Leprince et al. (2007, 2008); Ayoub et al. (2009). COSI-Corr uses the shift in phase of the images' Fourier transforms to map out the spatial translation of coherent windows of correlative pixels between two raster images, outputting two grids of orthogonal (Eastward and Northward, respectively) displacement components. We use these high resolution displacement maps to calculate the along-fault distribution of shallow strike-parallel slip, and to constrain the near-field in our inversions for slip at depth using InSAR LOS displacements.

The last Landsat-8 image of the area before the earthquake was captured a week before, on 30 November 2015; the next image was captured on 16 December 2015, just over a week after the earthquake. Comparing these two images would capitalize on the narrowest possible time span including the earthquake, but changing environmental conditions reduce the inherent correlation of these images. In general, natural environmental changes throughout the year produce discrepancies among optical imaging data, and these are pronounced in regions of high relief.

Seasonal signals introduce unwanted artefacts into the displacement field measured between any two images taken at different times in the year. Shifting sun angles cause apparent and spatially non-uniform north-south shifts throughout the scenes because the edges of shadows represent strong correlative features that change position as sun elevation migrates. Changing reflectance values due to snow cover, soil moisture, vegetation, etc. lead to temporal decorrelation between images, resulting in gaps in the displacement field. The rugged, high-altitude landscape of the Pamir exacerbates these seasonal variations. To reduce these seasonal errors in ground displacement, we correlate pairs of images that span the earthquake but are separated by approximately one full year (368 d, with the 16-d revisit interval of Landsat-8), taken on as close as possible to the same date in successive years so that sun angle and other seasonal environmental variables (e.g. snowcover) are similar (e.g. 29 September 2016 referenced to 27 September 2015, etc.—note 2016 is a leap-year with 366 d; Table 1; Fig. S1). Each 'path' (or-
<figure><img src="image_5.png" /><figcaption>**Figure 3.** Phase-unwrapped interferograms showing net Line-of-Sight range change for the observations, model, and misfit. Contours (thin black lines) are at 10 cm intervals. Arrows indicate orbit (long) and look direction (short) of satellite for the descending (track 005) and ascending (track 100) geometries. Thick black line represents simple fault used for initial uniform-slip inversion; thick dark red line represents 9-segment fault used for final slip-varying inversion. Phase-wrapped interferogram shown in Fig. S4.</figcaption></figure>

bital track) of Landsat scenes is captured every 16 d, totaling 22 acquisitions per year. However, only a subset of these scenes are sufficiently cloud-free to image the displacement of the solid earth surface. Ultimately a subset of 10 pre- and post-event pairs were usable from among the 22 pairs captured in the 2-yr period that spans the earthquake (Table 1 and Fig. S1).

We employed COSI-Corr's frequential correlator on the 15 m resolution panchromatic band 8 of Landsat-8's OLI (Operational Land Imager) sensor (Level-1 geotiff, UTM zone 43N, WGS84 datum). We used a correlation window size of 32×32 pixels (an area of 480 m per side) and a step size of 2 pixels (30 m) to obtain a displacement field as close to the raw image pixel resolution as
<table><caption>**Table 1.** Imaging dates of correlative pairs of Landsat-8 scenes (Worldwide Reference System path=151 row=033). Dates given in conventional Julian day-of-the-year (DoY) format, and relative to the event date (rel.), 7 December, that is, 341, 2015.</caption><thead><tr><th colspan="3">Pre-event scene</th><th colspan="3">Post-event scene</th></tr><tr><th>Calendar date</th><th>DoY<sup>a</sup></th><th>rel.</th><th>Calendar date</th><th>DoY</th><th>rel.</th></tr></thead><tbody><tr><td>13 Dec 2014</td><td>347</td><td>-359</td><td>16 Dec 2015</td><td>350</td><td>+9</td></tr><tr><td>29 Dec 2014</td><td>363</td><td>-343</td><td>1 Jan 2016</td><td>001</td><td>+25</td></tr><tr><td>14 Jan 2015</td><td>014</td><td>-327</td><td>17 Jan 2016</td><td>017</td><td>+41</td></tr><tr><td>19 Mar 2015</td><td>078</td><td>-263</td><td>21 Mar 2016</td><td>081</td><td>+106</td></tr><tr><td>6 May 2015</td><td>126</td><td>-215</td><td>8 May 2016</td><td>129</td><td>+154</td></tr><tr><td>23 Jun 2015</td><td>174</td><td>-167</td><td>25 Jun 2016</td><td>177</td><td>+202</td></tr><tr><td>11 Sep 2015</td><td>254</td><td>-87</td><td>13 Sep 2016</td><td>257</td><td>+282</td></tr><tr><td>27 Sep 2015</td><td>270</td><td>-71</td><td>29 Sep 2016</td><td>273</td><td>+298</td></tr><tr><td>13 Oct 2015</td><td>286</td><td>-55</td><td>15 Oct 2016</td><td>289</td><td>+314</td></tr><tr><td>30 Nov 2015</td><td>334</td><td>-7</td><td>2 Dec 2016</td><td>337</td><td>+362</td></tr></tbody></table>

is recommended for meaningful correlations (Ayoub *et al.* 2009). From each raster we remove statistical outlier pixels using COSI-Corr's non-local means filter, then further mask pixels that exhibit a signal-to-noise ratio of <0.98 (Leprince *et al.* 2007).

Each displacement map derived in this manner shows a discontinuity of similar shape and magnitude across the preexisting trace of the Sarez-Karakul fault, attesting to a robust signal of slip from this earthquake (Figs 4 and S1). However, slight geographical misalignment of pixel strips resulting from errors in satellite line-of-sight corrections lead to long wavelength artifacts of varying orientation and order (Landsat 8 Data Users Handbook, Department of the Interior U.S. Geological Survey 2016). These errors may stem from minor miscalculations or rounding within the 10 coordinate transformations among the satellite's attitude and the internal orientations of its optical components and detectors, and so they vary randomly from pass to pass amounting to non-systematic noise in the apparent temporal evolution of the pixel-tracking displacement field. Furthermore, the spatial coverage within individual scene comparisons varies among the pairs where particular pre- and post-event images have pixels that were unable to be correlated, leading to patches of no data and effectively lower resolutions. To improve spatial coverage of the displacement field and to reduce noise, we take the median value of each pixel among the 10 displacement maps, resulting in the comprehensively sampled horizontal displacement field in Fig. 4(b).

There are further signals in the displacement field imparted by real geomorphic processes that are not strictly coseismic tectonic deformation. These include landsliding and mistaken correlations among the largely decorrelated braided river channels, but they predominantly represent the downhill flow of alpine glaciers throughout the high Pamir. These signals are strong and coherent, and are thus not filtered out based on their signal-to-noise ratio. We mask them separately using Landsat-8's spectral bands for the classification of snow and ice. Because ice is highly reflective in the red and near-infrared (NIR, 0.85-0.89 mm), but almost completely absorptive in the short-wave infrared (SWIR, 1.56-1.66 mm, Paul *et al.* 2015, 2016; Winsvold *et al.* 2016), we automate the mapping of glaciers by detecting pixels with Red:SWIR ratios above a threshold value of 1 using Landsat-8's band 4 (Red) and band 6 (SWIR) as suggested in Paul *et al.* (2016). This excludes most high-displacement regions on the surface of icy glaciers, but rock-covered glacier tongues remain evident as small (10⁰ km) lobes of extremely high displacement at the mouths of valleys (Fig. 4).

After masking noise by removing these outliers from the results for each pair, and taking the median value for each pixel, we obtain

comprehensive maps of eastward and northward displacement ($U_E$ and $U_N$) at 30 m resolution (Fig. 1). We calculate the fault-parallel and fault-perpendicular components of the deformation field ($U_{\parallel}$ and $U_{\perp}$, respectively, with $U_{\parallel}$ positive northeast and $U_{\perp}$ positive northwest, in the dip direction of the SKF). We take a representative fault strike ($\theta$) of 034°, defined by the regional trend of the Sarez-Karakul fault and by the left-lateral nodal plane of the coseismic moment tensor from the USGS. The fault-oriented components are calculated as:

$$U_{\parallel} = U_N \cdot \cos(\theta) + U_E \cdot \sin(\theta)$$

$$U_{\perp} = U_N \cdot \sin(\theta) - U_E \cdot \cos(\theta).$$

Maps of these values are shown in Fig. 4. The smoothed, 30-m resolution horizontal displacement field from Landsat is combined with descending and ascending radar line-of-sight interferograms to improve the resolution, constrain the near-field deformation and resolve deformation in 3 dimensions for the inversion of slip distribution along the Sarez-Karakul fault (Section 4.1).

We also extract 1-km-wide fault-perpendicular swath profiles through the displacement field map in order to calculate a near-field surface slip distribution for comparison with geologic measurements of this and prior surface rupturing earthquakes (Figs 6 and 7). Offsets are calculated as the difference between median displacements within the 250-1000 m range on either side of the fault with uncertainties calculated as standard error (root mean square of the standard deviations on either side of the fault) (Figs 6 and 7; Supplementary spreadsheet). Values within 240 m of the fault are contaminated by slip in the opposite direction on the other side of the fault because of the correlation window size of 32×32 pixels (i.e. 480×480 m at pixel size 15 m), so we evaluate displacements from only pixels whose correlation windows did not include the other side of the fault, that is, greater than half a correlation window away (Fig. 6). We thus may lack local values of higher slip nearer the fault, however most investigations of near-fault coseismic deformation fields have found that discrete slip measurements of smaller apertures systematically miss distributed deformation within 10⁻¹-10² metres of a fault zone (Zachariasen & Sieh 1995; Shelef & Oskin 2010; Oskin *et al.* 2012; Zinke *et al.* 2014; Milliner *et al.* 2015; Vallage *et al.* 2015; Milliner *et al.* 2016; Gombert *et al.* 2018). As such we consider the omission of displacements within this zone fair in the representation of shallow (i.e. in the uppermost kilometre) coseismic fault slip, although not strictly comparable with geological measurements on discrete fault traces. To compare with modelled slip and to further reduce noise that may not reflect

<footnote><sup>a</sup>Landsat-8 scene names are indexed 'path' 'row' 'year' 'day of year', and so can be found as 'LC8151033YYYYDDDD'.</footnote>
<figure><img src="image_6.png" /><figcaption>**Figure 4.** Horizontal displacement fields spanning the 2015 earthquake calculated by pixel-tracking of 15-m Landsat-8 panchromatic imagery, presented in fault-orthogonal coordinate system. (a) East and north components of displacement as measured on respective seasonally matching pairs of Landsat-8 scenes that span the 7 December 2015 earthquake, labelled by day of the year. (b) Median stack of all year-spanning displacement maps, resolved into the fault-parallel direction (azimuth 034°E of N) to show left-lateral strike slip motion. Fault lines (black) represent geomorphic scarps visible in pre-event 1.5-m SPOT imagery; nine-fault model for source slip inversion, based on near-field discontinuity, shown as thin red line; 4-fault model used for initial uniform-slip inversion is shown in grey. (c) As in (b) but for the fault-normal direction—note the absence of a strong step in the signal across the predominantly strike-slip fault.</figcaption></figure>

real gradients in seismogenic slip (e.g. Gold *et al.* 2013), we group
the data into 1-km swath profiles, also presented in Fig. 7.

3.3 Geodetic inversion for slip at depth

We use both the downsampled InSAR LOS displacements and
the Landsat horizontal displacements jointly to invert for slip on

the fault plane at depth. Fault geometry is initially constrained by
performing a uniform-slip inversion using the smallest number of
planar fault segments that fairly represent the major bends and
steps of the SKF evident in prior mapping and in the Landsat-8
displacement map (the 4-segment model shown in Fig. 4). This
enables us to determine the fault dips and rakes of the major fault
segments based upon the far-field InSAR data while also honouring
<figure><img src="image_7.png" /><figcaption>**Figure 5.** Joint InSAR–Landsat slip distribution of the 2015 M7.2 Murghob earthquake on the 9-segment representation of the Sarez-Karakul fault that honors the surface discontinuity in the Landsat-8 displacement field. Note the highest values of surface slip in the southern half of the rupture, and the absence of shallow/surface slip in the northern reaches.</figcaption></figure>

<table><caption>**Table 2.** Modelled segments of the SKF defined from north to south for joint inversion of InSAR & Landsat-8 displacement fields. These faults correspond to the red traces in Fig. 1 and the planes displayed in cutaway perspective view in Fig. 5.</caption><thead><tr><th rowspan="2">Segment #</th><th rowspan="2">Strike</th><th rowspan="2">Dip</th><th rowspan="2">Rake</th><th colspan="2">Fault centre<sup>a</sup></th><th rowspan="2">Length(km)</th></tr><tr><th>Easting</th><th>Northing</th></tr></thead><tbody><tr><td>1</td><td>213.2</td><td>87.8</td><td>1.0</td><td>347.189</td><td>4283.78</td><td>15.5</td></tr><tr><td>2</td><td>182.2</td><td>87.8</td><td>1.0</td><td>342.812</td><td>4273.87</td><td>6.86</td></tr><tr><td>3</td><td>219.1</td><td>78.1</td><td>18.0</td><td>339.300</td><td>4266.29</td><td>10.7</td></tr><tr><td>4</td><td>237.1</td><td>67.7</td><td>-3.0</td><td>329.036</td><td>4257.68</td><td>16.4</td></tr><tr><td>5</td><td>219.4</td><td>79.7</td><td>-4.0</td><td>320.995</td><td>4251.82</td><td>3.64</td></tr><tr><td>6</td><td>204.1</td><td>79.7</td><td>-4.0</td><td>318.713</td><td>4247.89</td><td>5.52</td></tr><tr><td>7</td><td>226.2</td><td>79.7</td><td>-4.0</td><td>314.158</td><td>4242.09</td><td>9.50</td></tr><tr><td>8</td><td>199.7</td><td>79.7</td><td>-4.0</td><td>309.575</td><td>4235.57</td><td>6.86</td></tr><tr><td>9</td><td>208.8</td><td>79.7</td><td>-4.0</td><td>304.590</td><td>4225.37</td><td>15.9</td></tr></tbody></table>

the surface location of the fault as imaged in the Landsat-8 offsets. These four planes comprise the southern and northern reaches that have roughly the regional strike; a reach between the Muzkol range and the Kokuibel valley that has a more easterly strike; and the Officer range mountain front fault, which is set back eastward from the others by the bend and a left-step at the northern end of the Kokuibel valley (Fig. 2).

To solve for the fault geometry, we downsampled the InSAR deformation fields to 3156 and 2898 points for the ascending and descending tracks respectively using a nested uniform sampling approach with 2.5 and 5 km spacing in the medium and far-field respectively. We equally weight the two InSAR data sets as we only have a single interferogram in each of the ascending and descending look directions which are each differently sensitive to surface motion. Due to the level of the far-field noise in the Landsat-8 offsets that may affect the inversion for dip geometry in particular, we do not use this data set in the inversion for the uniform fault geometry (although we do use it to help define the surface trace location of the faults). Instead we formally use the Landsat-8-derived offsets for the subsequent variable slip model where the near-field, larger magnitude of the Landsat-8 offsets are above the noise and are more useful in constraining the along-strike distribution of shallow slip.

To model the InSAR displacements, we use the formulations of Okada (1985) for the surface displacements due to rectangular fault

planes dislocations at depth, assuming a uniform elastic half space and a rigidity of $3.2 \times 10^{10}$ Pa. To invert the downsampled data, we use a non-linear Powell's algorithm (Press et al. 1992) with Monte Carlo restarts to avoid local minima (Wright et al. 1999) to solve for the dip, rake and bottom depth of the four major fault segments. We fix the fault location, strike and lengths to the geomorphic fault trace and discontinuity in the Landsat-8 displacement field (Fig. 4). We fix the fault to the surface for the southern three segments and solve for the top depth of the fault for the most northern segment, where no sharp surface discontinuity is observed in either data set. We also solve for nuisance parameters in the two interferograms that take into account the arbitrary unwrapping point in the InSAR technique and also any long wavelength signals. The ascending data fit very well in the far-field, with some greater misfits (~25 cm) on the northern side of the rupture in the centre, and also to the north in the Kokuibel valley. The descending data are noisier in the far-field with 10 cm misfits in the Kokuibel valley and around the northermost segment. All data, models, and misfits are shown in Fig. S5.

Understanding the mechanics and dynamics of the fault rupture requires identifying variation in slip along strike, so we perform a variable-slip inversion using a further refined fault geometry with 9 segments (Figs 4 and 5; Table 1). Our near-field displacement map from Landsat pixel-tracking enables us to define a more detailed fault geometry that accurately captures undulations in the fault

<footnote><sup>a</sup> in UTM kilometres, Zone 43N.</footnote>
<figure><img src="image_8.png" /><figcaption>**Figure 6.** (a) Swath profiles through the Landsat-8 fault-parallel displacement field taken perpendicular to the simplified 4-segment fault trace at 1 km intervals, illustrating net motion of either side, which declines northward. (b) Representative swath profile through the Landsat-8 displacement field, illustrating the spread in the raw pixel data and various possible methods of deriving an offset measurement (explained further in Figs S2 and S3). Dots represent individual pixel values within the swath, coloured by distance across the swath. Thin black line represents swath median value at pixel resolution. Additional colored lines represent offsets derived from a variety of conventional approaches, for comparison. We prefer differencing the median values within 1 km of the fault rather than fitting lines to either side.</figcaption></figure>

plane that may influence the distribution of slip at a length scale of a few kilometres, that is, a scale over which surface slip distributions reflect mechanical and dynamical phenomena that control the propagation of rupture and the release of most seismic moment, which can modulate ground motion (Zachariasen & Sieh 1995; Oglesby 2008; Elliott et al. 2009; Gold et al. 2013; Lozos et al. 2013, 2015). We define a 9-segment model that is the simplest fault trace that fits within the discontinuity observed in Landsat displacements (Figs 4 and 5). Dips and rakes are taken as those calculated from the uniform-slip inversion for the corresponding fault sections—slightly west-dipping but nearly subvertical (78–88°) in all but the ‘bent’ section which requires a dip of 68° NW (Table 2). We extend the northern and southernmost fault segments beyond the extent of the uniform slip model and mapped surface ruptures to allow for subsurface slip beneath Lake Sarez and along the continuation of the mapped Sarez-Karakul Fault towards Lake Karakul. We downsample the Landsat-8 data as for the InSAR using a nested uniform approach (Figs S7 and S8). To determine the distribution of slip, we subdivide each of the nine fault segments into 1 km × 1 km patches and solve for the magnitude of slip on each subfault, as well as nuisance parameters (interferogram offset and planar long-wavelength tilt), in a linear inversion on each subfault following the method of Funning et al. (2005). A Laplacian smoothing constraint is used to regularize the inverse problem, as well as a positivity constraint. The Laplacian smoothing is applied across the fault segments and has boundary condition tapers to zero slip at the northern, southern and deepest segments, but is unconstrained for the surface slip patches. We select a smoothing value that balances the roughness of the slip solution with a meaningful reduction in the root mean square misfit to the data compared to the noise (Fig. S9).

The respective contributions to the solution from InSAR and from Landsat were weighted based on a minimization of residuals in the far-field. If the Landsat-8 offsets are equally weighted with

the InSAR data, large misfits are observed in the interferograms (Fig. S7). Reducing the weight of the Landsat-8 measurements to one fifth of those of the InSAR does not increase the misfit to the Landsat-8 markedly (10 per cent), but does greatly improve the InSAR fits (Fig. S8), halving the weighted misfit (Fig. S9). The resulting slip distribution from the joint inversion is shown in Fig. 5, and the model and misfits to the data in Fig. 3. The fit to the downsampled data is shown in Fig. S8. The variable-slip model and residuals for the InSAR-only solution (i.e. Landsat-8 with zero weight) are shown in Fig. S6 for comparison. The fit to the InSAR data improves over the uniform slip model (Fig. S5) with the greater number of free parameters from the increased fault slip patches and longer faults.

We calculate the uncertainty on the slip distribution (Fig. S10b) with a Monte Carlo approach in which we generate 100 perturbed data sets and invert for slip on each of these to make an estimate of the impact of the spatial noise in the data on the retrieved slip inversion (where the measure of uncertainty is the standard deviation of slip on each fault patch across the 100 inversion runs). We perturb the input data set using spatially correlated noise (Funning et al. 2005) estimated from the far-field region away from the deformation in each of the input data sets. We radially average the 2-D autocorrelation function of this selected region and fit an exponential function to estimate the maximum variance (0.06, 0.43, 1770 cm²) and e-folding wavelength (4.8, 3.9, 0.36 km) for the descending InSAR, ascending InSAR and Landsat-8 offsets, respectively.

The median stack of Landsat-8 displacements contains on average 6 months of post-seismic motion given the one year stepping intervals used, and will therefore potentially contain a greater magnitude of slip from post-seismic deformation (we assume there is no pre-seismic slip). There is a greater degree of systematic misfit for the Landsat-8 offsets in the north on the northeastern side of the
<figure><img src="image_9.png" /><figcaption>**Figure 7.** Slip distribution along the 2015 Sarez-Karakul fault rupture plotted to correspond spatially with surface rupture & horizontal displacement map below (coordinates in UTM Zone 43N—grid spacing 10 km). Fault offset is measured as the difference between median displacements within the 250–1000 m range away from either side of the fault along profiles through the fault-parallel Landsat-8 pixel-tracking map. Dots and shaded region represent displacements and standard deviations, respectively, from individual profiles. Continuous black line represents differential displacement across fault in 1-km-binned swath profiles (white lines on map). Red circles represent upper-kilometre offsets in joint InSAR-Landsat geodetic slip inversion. Pre-2015 fault scarps in black; 2015 surface ruptures in red; thick white line represents simplified fault model used in this study for slip inversion at depth. Lettered inset boxes refer to panels in Fig. 8 showing 1.5-m resolution satellite images of fault rupture.</figcaption></figure>

rupture in the Kokuibel valley (Fig. S8) that may be indicative of post-sesismic motion, consistent with the coseismic rupture being blind.

## 3.4 High resolution fault and rupture mapping

In order to map surface ruptures in detail and identify any topographic changes, we acquired 1.5-m resolution tri-stereo SPOT-6/7 optical imagery covering 700 sq. km around the Sarez-Karakul fault from before and after the earthquake. Pre-event tri-stereo imagery was available from acquisitions in October 2013, September 2014 and September 2015. Once we had defined the rupture extent using the first Sentinel-1 interferogram within a week of the earthquake,

we tasked Airbus to collect post-event tri-stereo imagery in a 5-km-wide swath along the coseismically active reach of the SKF (Fig. 1).

Digital elevation models constructed from the tri-stereo images were of insufficient resolution and precision to identify and trace the coseismic surface faulting in detail; thus, we relied on the orthorectified images themselves to locate and map fresh fault ruptures. We identified coseismic surface breaks in two ways: where displacements are large enough, surface ruptures appear as ‘moletacks’ or fracture arrays (e.g. Fig. 8a), or as sharp, dark scarps (e.g. Figs 8b and c) that are not present on pre-event imagery; in other places surface breaks are identified by subtle shadows and a narrow discontinuity across which features differentially change position when toggling between pre- and post-event imagery.
<figure><img src="image_10.png" /><figcaption>**Figure 8.** Pre- and post-event SPOT-6/7 imagery showing the formation of surface scarps along the mountainous portion of the SKF. Locations of these sites shown in Fig. 7. Coorindate grid of UTM Zone 43N shown for reference. 'Before' imagery dates from October 2013 or September 2014, before the snowy season; 'after' images are from 23 to 24 January, 2016, hence snowcover substantially alters surface reflectance. Panels (a) and (b) show fresh surface breaks in the steep slopes through the Muzkol Range. Panel (c) shows a newly formed scarp in the bed of a river between offset moraines on the SKF. Panel (d) shows the oblique-normal portion of the fault adjacent to the Officer range, where no new (i.e. 2015) surface breaks are observed in undeformed alluvial fans along the projection of the fault trace.</figcaption></figure>

# 4 RESULTS: RUPTURES ALONG THE SAREZ-KARAKUL FAULT

We first document our new analysis of the rupture parameters (length, slip, geometry, location) of the 2015 Murghob earthquake, then we present additional mapping and observations of prior ruptures along other reaches of the Sarez-Karakul fault, before seismologically evaluating these as possible sources of the 1911 Sarez earthquake.

## 4.1 Slip in the 2015 earthquake

We present new, comprehensive mapping of fresh surface ruptures from the 2015 earthquake visible in 1.5-m SPOT-6/7 images along the Sarez-Karakul fault and in higher resolution DigitalGlobe imagery available in Google Earth (Fig. 7). In conjunction with the map of discrete surface breaks, we present a 30-m resolution map of the near-field horizontal ground deformation that delineates shallow fault dislocation down to a zone of ~250 m width. We use this map to derive a slip distribution that may be inferred to represent shallow coseismic slip, that is, fault offset in the upper half-kilometre of the crust, which we may compare with our slip distribution at depth jointly inverted from Landsat and InSAR.

### 4.1.1 Slip distribution at depth

The slip distribution derived from inversion of InSAR and Landsat displacements (Fig. 5 and Supplementary file 'MurghabQuakeSlip-Inversion.csv') shows five notable features: (1) fault slip is almost entirely north of the hypocentre, with only a small deep patch of relatively minor slip south of Sarez Lake required by the far-field displacements from the interferograms, which is consistent with

the deep hypocentre location from the U.S. Geological Survey; (2) the greatest share of the slip is concentrated on the Sarez Dome section of the fault, immediately northeast of the epicentre, where it traverses crystalline bedrock and its highest elevations, achieving magnitudes of up to 5 m; (3) slip declines steadily northward along the 'bent' portion of the fault; (4) a further patch of relatively high slip occupies the deep portion of the fault beneath Kokuibel Valley, with little to no slip in the upper kilometre; (5) finite slip of <1 m continues on the northern section of the fault, although it is both blind and only a fraction (<25 per cent) of the slip on the Sarez Dome section. The deep slip beneath Kokuibel valley as modelled amounts to 2–3 m with an uncertainty of about 20 per cent (Supplementary Fig. S10). However, the relatively deep nature means that with surface geodetic data constraints alone, the resolution at 15-20 km depth is likely to be around 6 km [based on other checkerboard test of a similarly deep InSAR-constrained slip model, for example, Elliott & Copley (2013)]. Therefore, the slip imaged here is smoothed over a length-scale of this order of magnitude and cannot be resolved more finely. The distribution of slip we have resolved compares favorably with slip distributions inverted on simpler 3-segment fault models by Sangha et al. (2017), Metzger et al. (2017) and Gan et al. (2018), corroborating the findings that slip occurred in two main (supershear) pulses, interrupted by the major 22° rightward bend in the fault (Sangha et al. 2017).

### 4.1.2 Rupture length

The surface ruptures run discontinuously for 37 km through the Sarez Dome from the scree slope bounding Sarez Lake to the southern end of the Kokuibel valley (Fig. 7, Supplementary file '2015ruptures.kml'). The discontinuous nature of the fault traces
can be attributed to active ice and rock glaciers that overlie the fault, steep talus slopes that lack the competence to sustain discrete rupture and/or coseismic avalanches that have buried the surface rupture (Fig. 8a). No surface rupture nor substantial differential offset is evident on the southern shore of Lake Sarez (Figs 3 and 4), meaning the southern extent of surface rupture coincides with the quake's epicentre in the middle of the narrow lake. At the other end, remotely visible surface rupture terminates northward at 4260000N (38.472°N 72.082°E). This matches the northernmost ruptures observed in the field by Metzger et al. (2017), implying that satellite observations comprehensively capture the fresh scarps from this earthquake.

Importantly, conspicuous normal-faulting scarps apparent in the pre-event imagery north of 34.472° exhibit no new surface breaks in the 2016 images, with intervening alluvial fans appearing undeformed by faulting from the 7 December earthquake (Fig. 8d). These remote observations were corroborated in the field by Metzger et al. (2017) who identified the northernmost continuous ground ruptures at 38.472°N and observed that the normal-fault scarps along the Officer range had not been reactivated. (They identified a short zone of fractures with several cm left-lateral offset 10 km north of this site, but reported no fresh surface ruptures on the intervening fault scarps.) The absence of discrete coseismic surface faulting along this reach of the fault bounding the Officer range is consistent with greater width of the fault zone discontinuity in the Landsat-derived displacement field (Figs 4 and 7), and with the dearth of shallow slip in the joint Landsat-InsAR slip inversion (Fig. 5). Thus together detailed rupture mapping and the geodetic deformation field confirm that the northern third of the rupture remained blind. As such, the northward extent of faulting at depth is more difficult to define than the length of surface rupture, but our inversion requires slip along approximately 74 km of the fault.

### 4.1.3 Shallow slip distribution

The slip distribution presented in Fig. 7 represents fault offsets measured from the Landsat displacement map by differencing median values of fault-parallel displacement within a window extending 250–1000 m on either side of the fault (Fig. 6, further explanation in Figs S2 and S3). Results are shown from both individual pixel-width profiles, and from 1-km swath profiles to reduce short wavelength noise in the slip distribution which may not reflect the seismogenic fault slip process (Gold et al. 2013). The southern ~30 km of the rupture exhibit a relatively uniform 3–4 m of slip (possibly locally reaching up to 5 m), increasing gently northward until abruptly beginning to decline where the fault takes a sharp bend 22 degrees to the right (clockwise) at 28 km. This sudden decline in slip coincides spatially with the pause in rupture propagation inferred by Sangha et al. (2017). The slip decreases rather constantly along this 15-km misoriented reach of the fault, and drops more steeply from 2 to <1 m at exactly the terminus of surface rupture as identified in imagery. North of this, along the normal-faulted margin of the Kokuibel Valley and the Officer Range, a differential displacement of <1 m remains evident across the fault, but is no longer a localized discontinuity; rather, the deformed zone here has a width exceeding 1 km (Fig. 7). This is consistent with a deeply buried patch of slip as imaged in our joint inversion. In the far north, beyond the leftward stepover, a mild differential displacement persists for ~14 km, attesting to a small amount of slip at depth.

## 4.2 Pre-2015 rupture to north

To the north of the deformation field of the 2015 earthquake, the 1.5-m resolution pre-event imagery captures the active trace of the Sarez-Karakul fault as it continues to the southern shore of Karakul. Normal faulting scarps offset abandoned alluvial surfaces along the Kokuibel Valley margin (Fig. 8d), beneath which our slip inversion images ~2.5 m of slip in 2015 buried >15 km depth but no significant surface slip (Figs 4 and 5). These scarps have been noted by most previous studies as the conspicuous expression of the active, oblique SKF (Strecker et al. 1995; Rutte et al. 2017b). However, the surface faulting disappears through the active alluvial fans of the major rivers. A more northwesterly striking scarp emerges within the stepover zone at the north end of the Kokuibel Valley, striking discontinuously towards bedrock scarps of the fault, which outcrop from 38.632°N 73.196°E northeastward, as detailed in Strecker et al. (1995).

We map this reach of the fault, termed the Muzkol section, in new detail revealed by its clear manifestation in the SPOT imagery (Fig. 9, Supplementary file 'SKF_traces.kml'). The small bedrock scarp—suggestive of a single event, prior to 2015—continues northward up to the low drainage divide between the Kokuibel and Karakul basins, where it enters the sediments of the Muksu River (figs 2 and 9; Strecker et al. 1995). Here, the fault exhibits relatively fresh scarps and moletracks that can be traced across all but the youngest (presently active) alluvial surfaces (Fig. 9). While the absolute age of these scarps cannot be determined remotely, the fact that they cross the braided gravels of active alluvial channels here (e.g. Fig. 9c) indicates that they are younger than scarps in the Kokuibel Valley to the south and have not yet been eroded nor buried by surface processes, that is, they represent a more recent event than the penultimate rupture of the fault section that ruptured in 2015. In total, the fresh scarps on the northern SKF extend ~35 km from the north end of the Kokuibel Valley to the shore of Karakul. The fault lies beneath the lake for 28 more km before emerging on the north shore with clear surface expression offsetting palaeoshorelines. The 35 km rupture length implied by the continuous fresh scarps we can trace support an earthquake in excess of M=6.5 (e.g. Wells & Coppersmith 1994; Stirling et al. 2013).

The emergence of these scarps (i.e. the absence of similar surface fault scarps farther south) at around the left step in the SKF suggests that this structure may have served as the terminus of a prior rupture on the northern section of the fault. Furthermore, the difference in apparent scarp morphology before 2015 between this reach of the fault and the reach that exhibits surface rupture in 2015 indicates that the 2015 earthquake occurred along the less recently ruptured portion of the SKF, compatible with our expectations for slip recurrence. The fresh fault morphology in the north naturally makes that reach of the fault a strong candidate for the 1911 earthquake, which had a similar focal mechanism to the 2015 event (Kulikova et al. 2016).

## 4.3 SKF scarps south of Lake Sarez

The active trace of the SKF continues linearly southward from Lake Sarez and its intersection with the Sarez-Murghab (Pshart) Thrust system on the lake's southern bank. Previously published maps differ in the precise tracing of the SKF south of Lake Sarez (e.g. Strecker et al. 1995; Schurr et al. 2014), and no specific active trace is present for this reach in the most recent compilation of active faults by Mohadjer et al. (2016). In this high terrain (elevations uniformly >4 km), the fault is commonly obscured by active ice
<figure><img src="image_11.png" /><figcaption>**Figure 9.** SPOT-6/7 Imagery and interpretation diagrams of the Muksu-Muzkol section of the Sarez-Karakul fault, north of the Kokuibel Valley and the 2015 rupture. Area is up and right of centre in Fig. 2. Panel (a) serves as an index map for the remaining panels. Note the abundance of strike-slip fault morphology in young sedimentary deposits (panels b and c), including ponding of fine sediments and preserved fissures in the youngest abandoned surface (panel d). Small, likely single-event scarp is evident in bedrock hillslopes in (panel d). Coordinate grid in panels (b)–(d) is UTM Zone 43N.</figcaption></figure>

and rock glaciers, but occasionally crops out as scarps on bedrock hillsides and in Quaternary deposits (Fig. 10). We present mapping of these scarps as well as left-laterally offset ridges and valleys continuing for at least 70 km (and possibly 90 km) further south along a linear trend of 222° until they roughly merge with the Gunt shear zone (Figs 2 and 11, Supplementary file 'SKF_traces.kml'). Newly documented evidence of the active fault trace here includes a continuous linear scarp running along a fault line valley, which cuts the entire cross-section of a potentially left-laterally offset glaciated valley, evident in Google Earth imagery as shown in Fig. 11. Through much of this reach fault scarps are similarly present but discontinuous, largely due to poor preservation in the glaciated terrain. However, from 37.825°N 72.389°E for 30–50 km southwestward, we have mapped a more continuous fresh trace identified by Strom (2014) and Strom & Abdrakhmatov (2018) (Fig. 11a, Supplementary file '1911Source.kml') formed of intricate en echelon fractures (Fig. 11b), moletracks (Fig. 11c) and small scarps that pond groundwater and sediments (Fig. 11d). Based on the small dimensions of these individual features and the absence of offset abandoned Quaternary landforms, this scarp appears to show offset of only a few metres. Along this reach, we interpret the fine detail that is preserved and the relatively small apparent offsets as evidence of a single surface rupture more recent than on any other reach of the SKF prior to 2015.

## 4.4 Evidence of faulting in the Bartang/Murghab River Valley

Previous estimates of the 1911 Sarez earthquake have all placed the epicentre west of the Sarez-Karakul fault (albeit with uncertainties

of 10¹⁻¹⁰² km; Klotz & Galitzin 1915; Oldham 1923; Ambraseys & Bilham 2012; Schurr *et al.* 2014; Kulikova *et al.* 2016), which would appear to imply the existence of a yet-undocumented active left-lateral fault somewhere in the Bartang river valley, downstream of the Usoi landslide dam and the junction of the Kudara and Murghab tributaries, where they become the Bartang (Fig. 2). The bedrock geology surrounding this valley contains no major strike-slip faults (Stübner *et al.* 2013) so the fault would have to be a youthful structure expressed in the morphology and Quaternary sediments of the valley floor. Like previous researchers (e.g. Schurr *et al.* 2014; Metzger *et al.* 2017; Strom & Abdrakhmatov 2018), we find no strong evidence of young (historical) strike-slip faulting in the ranges surrounding Kulikova *et al.* (2016)'s epicentre. Within the Bartang river valley, however, we identify an assortment of ambiguous geomorphic features suggestive of active faulting, including linear scarps that pond sediments or align natural springs (Fig. 12). The anomalous width of the valley floor (ranging from 2 to 6 km, in contrast to other valleys surrounding it whose widths are limited to a few hundred meters, e.g. Figs 2 and 12) and the narrow low-relief bedrock ridge that runs up the middle do suggest that the valley shape may be modified by crustal faulting, perhaps in an oblique-normal sense as in the Kokuibel valley along the SKF, or possibly as active thrusts parallel to the Pshart thrust. Note that these ridges have specifically not been interpreted as landslide deposits in the thorough inventory of (Strom & Abdrakhmatov 2018, section 9.2, p. 285). Despite these potential indicators of faulting, highlighted in Fig. 12, there exist no scarps with the fresh morphology as seen along every reach of the SKF itself. We nonetheless cannot rule out the existence of a strike-slip fault parallel to the SKF running along the Bartang valley, where it would conceivably
<figure><img src="image_12.png" /><figcaption>**Figure 10.** Quaternary expression of the SKF south of Lake Sarez. Upper panel highlights large-scale geomorphic expression of the left-lateral fault, including a linear fault-line valley (indicated by small black arrows) and possible left offset of at least one glacially carved valley. Elevation is SRTM 30-m topography displayed in red-relief to emphasize elevation, slope, and topographic openness (Chiba *et al.* 2008; Kaneda & Chiba 2019). White lines show vantage/footprint of oblique view used in the lower panel. Lower panel shows perspective view from Google Earth indicating large scarps of the SKF (black arrows) including a small uphill-facing fault scarp damming Quaternary sediments (red arrows).</figcaption></figure>

have much lower preservation potential in the steep scree slopes or narrow river gorge further downstream.

## 5 RELATION TO THE 1911 EARTHQUAKE

Detailed mapping of the 2015 deformation field and pre-2015 scarps in the landscape using satellite imaging allows us to interrogate the succession of earthquakes on this fault system. In particular, the relative positions of the 1911 and 2015 earthquakes may indicate recurrence characteristics or a history of stress changes along this fault system. Previous studies have uniformly located the epicentre of the 1911 earthquake to the west or southwest of the 2015 epicentre (Klotz & Galitzin 1915; Oldham 1923; Ambraseys & Bilham 2012; Schurr *et al.* 2014; Kulikova *et al.* 2016), making them incompatible with rupture on the northern section of the fault that otherwise appears to be recently ruptured. However, these previously catalogued epicentres for the 1911 earthquake are based on either collocation

with the Usoi landslide (originally thought to be the source of the teleseismic signal, Klotz & Galitzin 1915; Jeffreys 1923; Oldham 1923), on sparse macroseismic data (Ambraseys & Bilham 2012) or on relatively crude early seismograms with significant timing errors (Kulikova *et al.* 2016). Thus the reported epicentres either do not refer to the tectonic source, or have location uncertainties exceeding 50 km. Furthermore, all of them predate detailed neotectonic mapping of the SKF system. With our new detailed mapping of the SKF system based on contemporary remote sensing, we can revisit the historical seismic recordings to test various realistically plausible sources of the 1911 earthquake constrained by observations of faulting in addition to seismological data.

### 5.1 Possible sources of the 1911 Sarez earthquake

Based on our mapping of the SKF on imagery from both before and after the 2015 earthquake, we suggest three plausible strike-slip sources for the 1911 seismic event (Fig. 13, Table 3). We list them here in reverse-order of their hypothesized plausibility based on available evidence for a recent earthquake.

(i) **Bartang/Murghab River Valley** – Section 4.4 – Because all previous assessments of the 1911 Sarez earthquake have located it west of the 2015 epicentre, our first candidate fault runs along the axis of the Bartang River Valley, downstream of the Usoi landslide dam and the junction between the Kudara and Murghab tributaries.

(ii) **Northern SKF (Muksu-Muzkol Section)** – Section 4.2 – This section presents strong continuous evidence of Quaternary faulting along its entirety, including small scarps and scarp arrays in the youngest abandoned alluvial deposits. However, it is the farthest of the candidates from any prior assessment of the earlier earthquake's epicentre.

(iii) **Southern SKF (Gunt-Sarez Section)** – Section 4.3 – While evidence of Quaternary faulting exists sporadically along the 80 km of the SKF identifiable south of Sarez Lake, the southern 30 km of this reach, previously termed the 'Pathur-Nemos Fault' (Strom 2014; Strom & Abdrakhmatov 2018) exhibit the strongest evidence in the region of a single recent surface rupture before 2015.

### 5.2 Seismological reevaluation of 1911 Earthquake source

In order to test which among these fault sections most likely produced the 1911 earthquake, we compare the observed seismograms presented in Kulikova *et al.* (2016) with synthetic seismograms for each respective fault source. The original digital reevaluation of them located their origin somewhere west of the 2015 earthquake (Figs 1, 2 and 13), but large timing errors led to >60 km uncertainties in the epicentre (Kulikova *et al.* 2016). A supplementary approach to constrain the epicentre is to test which of the plausible sources we have identified produces seismograms that best fit the digitized traces from 1911.

The original 1911 seismic records were digitized and processed from European stations as described in Kulikova *et al.* (2016). Synthetic seismograms were calculated for each station and component based on known instrument parameters, a layered global velocity model, and the respective source locations listed in Table 3 and described in Section 5.1. Based on the similarity in magnitude and mechanism between 1911 and 2015, the initial source parameters (source depth and rupture width and length) were taken to be similar, with the 1911 Sarez event slightly larger than 2015. Further, we let the source parameters vary in a grid search looking for the best waveform fit between observed and synthetic seismograms. Our
<figure><img src="image_13.png" /><figcaption>**Figure 11.** Detailed mapping of fresh surface rupture scarps 50 km south of the 2015 rupture, identified as the Pathur-Nemos Fault in Strom (2014) and Strom & Abdrakhmatov (2018). Area occupies lower left of Fig. 2. Hillshade elevation in panel (a) is 30-m SRTM (Shuttle Radar Topography Mission). Panel (b) shows a perspective view in Google Earth from the vantage point indicated on the index map in panel (a). Panels (c) and (d) show Digital Globe orthoimagery from Google Earth, bare and annotated, respectively, to highlight evidence of faulting. Note the exceptionally fresh scarp morphology with small offsets and fine detail preserved, and still damming stream channels to form ponds.</figcaption></figure>

detailed fault mapping permits us to constrain the strike of each potential source, so their azimuths align along the respective mapped traces as shown in Fig. 13. The best waveform fit was found for the uniform 80 km long × 20 km wide fault with rupture propagating almost unilaterally northwards with rupture velocity of 2 km per second (for all three test sources). Fault centroids as well as epicentres are listed in Table 3, with faults and epicentres shown in Fig. 13.

The comparison among the sources was accomplished by cross correlation between observed and synthetic waveforms, and a misfit is calculated as the inverse of the correlation coefficient between the relevant waveforms (Table 3). Because the digitized analog records have substantial uncorrected instrumental response, timing irregularity, and digitization artifacts, correlation of the full traces cannot provide a suitable metric for comparison. Instead we evaluate only the portions of the traces with clear body-wave phase arrivals. We evaluated misfits among synthetic and observed seismograms for *P*, *PP*, *S* and *SS* phase arrivals, as was done in Kulikova *et al.* (2016). The best-fitting source is the *Bartang valley* fault plane (Source 1) and the poorest fit results from the *Northern SKF* source (Source 2), however there is little differentiation between the misfits of the *Bartang valley* (Source 1) and the *Southern SKF* source (Source 3, Table 3; Fig. 14). The ambiguity between the latter two sources may arise from compounded timing errors when using multiple phase arrivals per record. Thus, we re-evaluate the misfits using only the clearest and strongest arrivals—the first-arriving *P*

and *S* waves. Using only these two phases, we observe stronger differentiation in the compatibility of various sources (Table 3).

For both tests, the northern fault source remains the least compatible with observed waveforms (Fig. 14). In the test using all four body wave phase arrivals, there is insufficient difference in correlation (misfit) to distinguish whether a *Bartang river valley* source or a *Southern SKF* source fits the observations better. However when only the *P* and *S* arrivals are used there is better correlation between the observed and synthetic seismic waveforms of the *Southern SKF* source than the *Bartang river valley* source. Visually inspecting the waveforms (e.g. Fig. 14) we observe better waveform fit for Source 3 in comparison to Source 1, particularly evident for the *E* and *Z* components of the *S* wave.

Thus, we propose that the southernmost mapped section of the SKF—that identified as the Pathur-Nemos Fault by Strom & Abdrakhmatov (2018)—from near at least the Gunt river in the south (37.738°N 72.202°E) to potentially as far north as Lake Sarez itself, a length of up to 75 km, was the source of the 1911 earthquake that damaged the region and caused the catastrophic Usoi landslide at the northern end of its extent.

## 6 ACTIVE FAULTING AND BEDROCK STRUCTURE

The overall structure of the Sarez-Karakul fault can be defined by sections whose boundaries coincide with the intersections of major
<figure><img src="image_14.png" /></figure>

**Figure 12.** Perspective view southwestward down the Bartang river valley from approximately 38.33°N 72.42°E in Google Earth; area located at centre left in Fig. 2. Upper panel highlights geomorphic features of interest including the Usoi landslide dam, the anomalous width of the Bartang valley, and the medial ridge that cannot be explained by landsliding (Strom & Abdrakhmatov 2018). Elevation is SRTM 30-m topography presented in red-relief to highlight elevation, slope, and openness (Chiba et al. 2008; Kaneda & Chiba 2019). Annotations in lower panel indicate potential evidence of tectonic faulting, including linear scarps, alignment of natural springs, and the anomalous topography of the valley floor. Note that although this is the site of the relocated 1911 epicentre from Kulikova et al. (2016), geomorphic scarps like those seen all along the SKF are absent here.

bedrock shear zones and other, roughly perpendicular active faults throughout the Central and Northern Pamir (Section 2.3, Fig. 2). Indeed, the coseismic rupture of the fault in 2015 was bounded by two of those major structures: the Pshart (Sarez-Murghab) thrust in the south, and the stepover coincident with the Tanymas thrust in the north. Furthermore coseismic slip appears to have been modulated by a third intersecting fault zone, insofar as the bend along which surface slip declines from its maximum to zero coincides with the reach of the SKF that runs along or intersects obliquely the North Sarez shear zone, an early Miocene detachment fault bounding the Sarez Dome (Rutte et al. 2017a; Sangha et al. 2017).

Similarly, fresh ruptures on the southern SKF can be traced northward from roughly the fault's junction with the Gunt shear zone. While the freshest scarps appear at the southern end of this reach,

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

**Figure 13.** Plausible fault sources for the 1911 earthquake determined based on observed 2015 deformation field and detailed fault mapping on satellite imagery. Seismological origins of the two events and $2\sigma$ error are plotted as calculated by Kulikova et al. (2016), along with our corresponding maps of surface rupture. Geometric details of each test source (centre coordinates and azimuth/strike) are listed in Table 3, as well as the misfit between synthetic seismograms for each and those observed in 1911. Only the southern SKF and the Bartang Valley are reasonable sources of the 1911 earthquake, and only the southern SKF shows clear evidence of a recent rupture.

evidence of Quaternary offsets persists northward to the Pshart Thrust, but are largely absent beyond, in the Sarez Dome section where the majority of slip in the 2015 event took place. Thus the distribution of timing of slip along the SKF appears to be divided by the structural boundaries of bedrock suture zones. Far from unexpected, this behaviour may reflect the influence on dynamic rupture propagation exerted by geometric changes (steps and bends) and rheological changes (effective friction and cohesion). This suggests observable parameters that may help interpret the lengths of past ruptures and forecast the lengths of future events. These qualitative observations of correlation between rupture boundaries and bedrock tectonic structure will be important to interrogate rigorously through dynamic, physics-based modelling of fault ruptures, as well as more detailed field investigations informed by these observations of coseismic slip.

Along-strike changes in geometry and rheology inherently change the parameters that permit rupture propagation (Harris & Day 1997; Duan & Oglesby 2005, 2006; Lozos et al. 2011), and the development of a 'boundary' by the termination of one rupture may induce an apparent segmentation by leading to differing rupture histories along different reaches of the fault (Duan & Oglesby 2005; Schwartz et al. 2012; Duross et al. 2016; Duan et al. 2019). The degree to which such boundaries become entrenched as persistent barriers to rupture is of paramount importance in forecasting the suite of ruptures that may occur along a fault system.
<table><caption>**Table 3.** Results of synthetic seismogram test for plausible fault sources of 1911 seismic records.</caption><thead><tr><th rowspan="2">Source #</th><th rowspan="2">Fault</th><th colspan="3">Centroid</th><th colspan="2">Misfits<sup>a</sup></th></tr><tr><th>Lat.</th><th>Lon.</th><th>Az.</th><th>P, PP, S, SS</th><th>P&S only</th></tr></thead><tbody><tr><td>1</td><td>Bartang R.</td><td>38.2123</td><td>72.2612</td><td>050</td><td>1.51</td><td>1.46</td></tr><tr><td>2</td><td>SKF north</td><td>38.8066</td><td>73.3041</td><td>020</td><td>1.73</td><td>1.67</td></tr><tr><td>3</td><td>SKF south</td><td>37.9798</td><td>72.5363</td><td>042</td><td>1.53</td><td>1.20</td></tr><tr><td>-</td><td>2015</td><td>38.4007</td><td>72.9910</td><td>036</td><td>-</td><td>-</td></tr></tbody></table>

<figure><img src="image_16.png" /><figcaption>**Figure 14.** Synthetic and observed waveforms for the 1911 Sarez-Pamir earthquake. The plot shows the synthetic (red lines) and observed (black lines) data for vertical Z and horizontal N and E components of GTT (Göttingen, Germany) station for three different test sources (Fig. 13).</figcaption></figure>

Now that its rupture behaviour has been illuminated by two historic earthquakes, the specific geometric relationships between the oblique left-normal SKF and these bedrock shear zones warrant further investigation through detailed field mapping, as the former may be shaped by exploitation of these preexisting structures. Recognition of such systematic behavior may improve estimates of the suite of earthquakes possible on a given fault system, and enhance our ability to recognize features that may control that distribution.

## 7 CONCLUSIONS

The 7 December 2015 $M_w$ 7.2 Murghob earthquake occurred along the transtensional Sarez-Karakul fault, separating the northward-convergent eastern Pamir from the western Pamir, which is deforming under both northward convergence due to the collision of

India with Eurasia and westward gravitational collapse into the relatively low-lying Tajik-Afghan Basin. We use a novel approach of stacking year-spanning optical-image pixel displacements to derive a comprehensive near-field horizontal deformation map from this earthquake, which we jointly invert with Sentinel-1 interferograms to compute the deep and shallow distributions of slip along the ~75–80 km rupture. Along with our slip distribution we present the first comprehensive map of surface ruptures from this earthquake. Prior to this study the only surface ruptures of the 2015 earthquake that had been confirmed were the northernmost 4–5 km. We present evidence of clear but discontinuous scarps formed coseismically reaching from the shore of Lake Sarez 37 km northward to the terminus of the surface rupture identified by Metzger *et al.* (2017) in the Kokuibel valley.

The boundaries of the rupture and the locations of changes in the magnitude of coseismic slip coincide with the intersections of major Cenozoic suture zones and/or crustal detachment faults with the SKF, suggesting that these inherited tectonic structures

<footnote><sup>a</sup>Misfit ~ 1/correlation (scaled to number of stations).</footnote>
not only compartmentalize the geological domains of the Pamir but also play a role in controlling earthquake rupture extents along this cross-cutting structure in the shallow crust.

Our mapping of the neotectonic expression of the SKF beyond the 2015 rupture zone reveals coherent reaches of prior surface rupture that may be interpreted as previous seismic events. These observations provide plausible sources of the 1911 earthquake, for each of which we compute synthetic seismograms for comparison with analog records. Poor misfits of the northernmost potential source indicate that fresh ruptures observed there represent a third earthquake prior to the known historic events on this fault system. Favourably low misfit values for a more southern source for 1911 suggest that newly mapped fresh ruptures south of Sarez Lake may finally reveal the culprit of the enigmatic earthquake. The variable apparent age of surface rupture along the trace of the SKF represents record of multiple distinct earthquakes, extending our recognition of slip events on this fault system and revealing an out-of-sequence series of earthquake ruptures along the strike of this fault.

## ACKNOWLEDGEMENTS

This work was supported by the UK National Environmental Research Council (NERC)'s Looking Inside the Continents from Space (LiCS) grant (NE/K011006/1), the NERC-ERSC Earthquakes without Frontiers (EwF) consortium (EwF/NE/J02001X/1-1) and the Centre for the Observation and Modelling of Earthquakes, Volcanoes and Tectonics (COMET, GA/13/M/031, http://comet.nerc.ac.uk). SPOT-6/7 imagery was purchased from Airbus through Apollo Mapping; Sentinel-1 radar images were obtained from ESA and this work contains modified Copernicus Sentinel data 2015-2016; Landsat-8 scenes are freely available through the USGS Earth Explorer portal. JRE is supported by a Royal Society University Research fellowship (UF150282). Two anonymous reviews guided improvements to the presentation of content and figures in this manuscript.

## REFERENCES

* Ambraseys, N.N. & Bilham, R.G., 2012. The Sarez-Pamir earthquake and landslide of 18 February 1911, *Seismol. Res. Lett.*, **83**(2), 294–314.

* Ayoub, F., Leprince, S. & Avouac, J.P., 2009. Co-registration and correlation of aerial photographs for ground deformation measurements, *ISPRS J. Photogram. Remote Sens.*, **64**(6), 551–560.

* Berberian, M., Jackson, J., Qorashi, M., Khatib, M.M., Priestley, K., Talebian, T. & Ghafuri-Ashtiani, M., 1999. The 1997 May 10 Zirkuh (Qa'enat) earthquake (M-w 7.2): faulting along the Sistan suture zone of eastern Iran, *Geophys. J. Int.*, **136**(3), 671–694.

* Bilham, R., 2019. Himalayan earthquakes: a review of historical seismicity and early 21st century slip potential, *Geol. Soc., Lond., Spec. Publ.*, **483**, 423–482.

* Billington, S., Isacks, B.L. & Barazangi, M., 1977. Spatial distribution and focal mechanisms of mantle earthquakes in the Hindu Kush-Pamir region: a contorted Benioff zone, *Geology*, **5**, 699–704.

* Bindi, D., Parolai, S., Gómez-Capera, A., Locati, M., Kalmetyeva, Z. & Mikhailova, N., 2014. Locations and magnitudes of earthquakes in Central Asia from seismic intensity data, *J. Seismol.*, **18**, 1–21.

* Bourgeois, O., Cobbold, P.R., Rouby, D., Thomas, J.C. & Shein, V., 1997. Least squares restoration of Tertiary thrust sheets in map view, Tajik depression, central Asia, *J. geophys. Res.*, **102**(B12), 27 553–27 573.

* Burtman, V.S. & Molnar, P., 1993. Geological and geophysical evidence for deep subduction of continental crust beneath the Pamir, *Geol. Soc. Am. Spec. Paper*, **281**, 76.

* Chiba, T., Kaneta, S.I. & Suzuki, Y., 2008. Red relief image map: new visualization method for three dimensional data, in *The International Archives of the Photogrammetry, Remote Sensing and Spatial Information Sciences*, Vol. XXXVII, Part B2, pp. 1071–1076.

* Coutand, I., Strecker, M.R., Arrowsmith, J.R., Hilley, G., Thiede, R.C., Korjenkov, A. & Omuraliev, M., 2002. Late Cenozoic tectonic development of the intramontane Alai Valley, (Pamir-Tien Shan region, central Asia): an example of intracontinental deformation due to the Indo-Eurasia collision, *Tectonics*, **21**(6), 3–1–3-19.

* Crone, A.J. & Haller, K.M., 1991. Segmentation and the coseismic behavior of Basin and Range normal faults: examples from east-central Idaho and southwestern Montana, U.S.A., *J. Struct. Geol.*, **13**(2), 151–164.

* Depolo, C.M., Clark, D.G., Slemmons, D.B. & Lli, A.R.R.A.M., 1991. Historical surface faulting in the Basin and Range province, western North America: implications for fault segmentation, *J. Struct. Geol.*, **13**(2), 123–136.

* Duan, B., Liu, Z. & Elliott, A.J., 2019. Multicycle dynamics of the Ackay bend along the Altyn Tagh Fault in Northwest China, Part 2: the realistically complex fault geometry, *Tectonics*, **38**(3), 1120–1137.

* Duan, B. & Oglesby, D.D., 2005. Multicycle dynamics of nonplanar strike-slip faults, *J. geophys. Res.*, **110**(B3), doi:10.1029/2004JB003298.

* Duan, B. & Oglesby, D.D., 2006. Heterogeneous fault stresses from previous earthquakes and the effect on dynamics of parallel strike-slip faults, *J. geophys. Res.*, **111**(B5), doi:10.1029/2005JB004138.

* Duross, C.B., Personius, S.F., Crone, A.J., Olig, S.S., Hylland, M.D., Lund, W.R. & Schwartz, D.P., 2016. Fault segmentation: new concepts from the Wasatch Fault Zone, Utah, USA, *J. geophys. Res.*, **121**(2), 1131–1157.

* Elliott, A.J., Dolan, J.F. & Oglesby, D.D., 2009. Evidence from coseismic slip gradients for dynamic control on rupture propagation and arrest through stepovers, *J. geophys. Res.*, **114**(B2), B02313.

* Elliott, J.R. & Copley, A., 2013. The 2011 Mw 7.1 Van (Eastern Turkey) earthquake, *J. geophys. Res.*, **118**, 1–19.

* Funning, G.J., Barke, R.M., Lamb, S.H., Minaya, E., Parsons, B. & Wright, T.J., 2005. The 1998 Aiquile, Bolivia earthquake: a seismically active fault revealed with InSAR, *Earth planet. Sci. Lett.*, **232**(1–2), 39–49.

* Gan, J., Hu, J., Li, Z., Yang, C., Liu, J., Sun, Q. & Zheng, W., 2018. Mapping three-dimensional co-seismic surface deformations associated with the 2015 MW7.2 Murghab earthquake based on InSAR and characteristics of crustal strain, *Sci. China Earth Sci.*, **61**(10), 1451–1466.

* Gold, P.O., Oskin, M.E., Elliott, A.J., Hinojosa-Corona, A., Taylor, M.H., Kreylos, O. & Cowgill, E., 2013. Coseismic slip variation assessed from terrestrial lidar scans of the El Mayor-Cucapah surface rupture, *Earth planet. Sci. Lett.*, **366**, 151–162.

* Goldstein, Richard M. & Werner, Charles L., 1998. Radar interferogram filtering for geophysical applications, *Geophys. Res. Lett.*, **25**(21), 4035–4038, doi.org/10.1029/1998GL00033.

* Gombert, B., Duputel, Z., Jolivet, R., Doubre, C., Rivera, L. & Simons, M., 2018. Revisiting the 1992 Landers earthquake: a Bayesian exploration of co-seismic slip and off-fault damage, *Geophys. J. Int.*, **212**, 839–852.

* Harris, R.A. & Day, S.M., 1997. Effects of a low-velocity zone on a dynamic rupture, *Bull. seism. Soc. Am.*, **87**(5), 1267–1280.

* Ischuk, A. et al., 2013. Kinematics of the Pamir and Hindu Kush regions from GPS geodesy, *J. geophys. Res.*, **118**(5), 2408–2416.

* Jay, C.N., Flesch, L.M. & Bendick, R.O., 2017. Kinematics and dynamics of the Pamir, Central Asia: Quantifying surface deformation and force balance in an intracontinental subduction zone, *J. geophys. Res.*, **122**(6), 4741–4762.

* Jeffreys, H., 1923. The Pamir Earthquake of 1911 February 18, in relation to the depths of earthquake foci, *Mon. Not. R. Astron. Soc. Geophys. Suppl.*, **1**(s2), 22–31.

* Jepson, G., Glorie, S., Konopelko, D., Gillespie, J., Danišik, M., Evans, N.J., Mamadjanov, Y. & Collins, A.S., 2018. Thermochronological insights into the structural contact between the Tian Shan and Pamirs, Tajikistan, *Terra Nova*, **30**(2), 95–104.

* Kaneda, H. & Chiba, T., 2019. Stereopaired morphometric protection index red relief image maps (Stereo MPI-RRIMs): effective visualization of high-resolution digital elevation models for interpreting and mapping small tectonic geomorphic features, *Stereo MPI-RRIMs: effective visualization*, *Bull. seism. Soc. Am.*, **109**(1), 99–109.
Klotz, O. & Galitzin, P.B., 1915. Earthquake of February 18, 1911, *Bull. seism. Soc. Am.*, **5**, 206-213.

Koulakov, I. & Sobolev, S.V., 2006. A tomographic image of Indian lithosphere break-off beneath the Pamir-Hindukush region, *Geophys. J. Int.*, **164**(2), 425-440.

Kufner, S.K., Schurr, B., Haberland, C., Zhang, Y., Saul, J., Ischuk, A. & Oimahmadov, I., 2017. Zooming into the Hindu Kush slab break-off: a rare glimpse on the terminal stage of subduction, *Earth planet. Sci. Lett.*, **461**, 127-140.

Kufner, S.K., Schurr, B., Sippl, C., Yuan, X., Murodkulov, S., Schneider, F., Mechie, J. & Tilmann, F., 2016. Deep India meets deep Asia: Lithospheric indentation, delamination and break-off under Pamir and Hindu Kush (Central Asia), *Earth planet. Sci. Lett.*, **435**, 171-184.

Kulikova, G., Schurr, B., Krüger, F., Brzoska, E. & Heimann, S., 2016. Source parameters of the Sarez Pamir earthquake of February 18, 1911, *Geophys. J. Int.*, **205**, 1086-1098.

Leprince, S., Barbot, S., Ayoub, F. & Avouac, J.P., 2007. Automatic and precise orthorectification, coregistration, and subpixel correlation of satellite images, application to ground deformation measurements, *IEEE Trans. Geosci. Remote Sens.*, **45**(6), 1529-1558.

Leprince, S., Berthier, E., Ayoub, F., Delacourt, C. & Avouac, J.P., 2008. Monitoring earth surface dynamics with optical imagery, *EOS, Trans. Am. geophys. Un.*, **89**(1), 1-2.

Lozos, J.C., Oglesby, D.D. & Brune, J.N., 2013. The effects of fault stepovers on ground motion, *Bull. seism. Soc. Am.*, **103**(3), 1922-1934.

Lozos, J.C., Oglesby, D.D., Brune, J.N. & Olsen, K.B., 2015. Rupture propagation and ground motion of strike-slip stepovers with intermediate fault segments, *Bull. seism. Soc. Am.*, **105**(1), 387-399.

Lozos, J.C., Oglesby, D.D., Duan, B. & Wesnousky, S.G., 2011. The effects of double fault bends on rupture propagation: a geometrical parameter study, *Bull. seism. Soc. Am.*, **101**(1), 385-398.

Macelwane, J.B., 1926. Are important earthquakes ever caused by impact?, *Bull. seism. Soc. Am.*, **16**(4), 15-18.

Metzger, S. et al., 2017. The 2015 Mw7.2 Sarez strike-slip earthquake in the pamir interior: response to the underthrusting of india's western promontory, *Tectonics*, **36**(11), 2407-2421.

Milliner, C.W., Dolan, J.F., Hollingsworth, J., Leprince, S. & Ayoub, F., 2016. Comparison of coseismic near-field and off-fault surface deformation patterns of the 1992 Mw7.3 Landers and 1999 Mw7.1 Hector Mine earthquakes: implications for controls on the distribution of surface strain, *Geophys. Res. Lett.*, **43**(19), 1015-10124.

Milliner, C.W., Dolan, J.F., Hollingsworth, J., Leprince, S., Ayoub, F. & Sammis, C.G., 2015. Quantifying near-field and off-fault deformation patterns of the 1992 Mw7.3 Landers earthquake, *Geochem. Geophys. Geosyst.*, **16**(5), 1577-1598.

Mohadjer, S., Alan Ehlers, T., Bendick, R., Stübner, K. & Strube, T., 2016. A quaternary fault database for central Asia, *Nat. Hazards Earth Syst. Sci.*, **16**(2), 529-542.

Mohadjer, S. et al., 2010. Partitioning of India-Eurasia convergence in the Pamir-Hindu Kush from GPS measurements, *Geophys. Res. Lett.*, **37**(4), 1-6.

Negredo, A.M., Replumaz, A., Villaseñor, A. & Guillot, S., 2007. Modeling the evolution of continental subduction processes in the Pamir-Hindu Kush region, *Earth planet. Sci. Lett.*, **259**(1-2), 212-225.

Oglesby, D.D., 2008. Rupture termination and jump on parallel offset faults, *Bull. seism. Soc. Am.*, **98**(1), 440-447.

Okada, Y., 1985. Surface deformation due to shear and tensile faults in a half-space, *Bull. seism. Soc. Am.*, **75**(4), 1135-1154.

Oldham, R.D., 1923. The Pamir earthquake of 18th February, 1911, *Quart. J. Geol. Soc.*, **79**(1923), 237-245.

Ong, S.Q.M., Barbot, S. & Hubbard, J., 2019. Physics-based scenario of earthquake cycles on the Ventura Thrust System, California: the effect of variable friction and fault geometry, *Pure appl. Geophys.*, **176**, 3993-4007.

Oskin, M.E. et al., 2012. Near-field deformation from the El Mayor-Cucapah earthquake revealed by differential LIDAR, *Science*, **335**, 702-705.

Paul, F., Winsvold, S.H., Kääb, A., Nagler, T. & Schwaizer, G., 2016. Glacier remote sensing using Sentinel-2. Part II: mapping glacier extents and

surface facies, and comparison to Landsat 8, *Remote Sens.*, **8**(7), 575, doi:10.3390/rs8070575.

Paul, F. et al., 2015. The glaciers climate change initiative: methods for creating glacier area, elevation change and velocity products, *Remote Sens. Environ.*, **162**, 408-426.

Pegler, G. & Das, S., 1998. An enhanced image of the Pamir-Hindu Kush seismic zone from relocated earthquake hypocentres, *Geophys. J. Int.*, **134**(2), 573-595.

Philibosian, B. et al., 2014. Rupture and variable coupling behavior of the Mentawai segment of the Sunda megathrust during the supercycle culmination of 1797 to 1833, *J. geophys. Res.*, **119**(9), 7258-7287.

Pozzi, J.P. & Feinberg, H., 1991. Paleomagnetism in the Tajikistan: continental shortening of European margin in the Pamirs during Indian Eurasian collision, *Earth planet. Sci. Lett.*, **103**(1-4), 365-378.

Press, W.H., Teukolsky, S.A., Vetterling, W.T. & Flannery, B.P., 1992. *Numerical Recipes in C: The Art of Scientific Computing*, 2nd edn, Cambridge Univ. Press.

Robinson, A.C., 2015. Mesozoic tectonics of the Gondwanan terranes of the Pamir plateau, *J. Asian Earth Sci.*, **102**, 170-179.

Robinson, A.C., Ducea, M. & Lapen, T.J., 2012. Detrital zircon and isotopic constraints on the crustal architecture and tectonic evolution of the northeastern Pamir, *Tectonics*, **31**(2), 1-16.

Rutte, D., Ratschbacher, L., Schneider, S., Stübner, K., Stearns, M.A., Gulzar, M.A. & Hacker, B.R., 2017b. Building the Pamir-Tibetan Plateau-Crustal stacking, extensional collapse, and lateral extrusion in the Central Pamir: 1. Geometry and kinematics, *Tectonics*, **36**(3), 342-384.

Rutte, D. et al., 2017a. Building the Pamir-Tibetan Plateau-Crustal stacking, extensional collapse, and lateral extrusion in the Central Pamir: 2. Timing and rates, *Tectonics*, **36**(3), 385-419.

Sangha, S., Peltzer, G., Zhang, A., Meng, L., Liang, C., Lundgren, P. & Fielding, E., 2017. Fault geometry of 2015, Mw7.2 Murghab, Tajikistan earthquake controls rupture propagation: insights from InSAR and seismological data, *Earth planet. Sci. Lett.*, **462**, 132-141.

Schneider, F.M. et al., 2013. Seismic imaging of subducting continental lower crust beneath the Pamir, *Earth planet. Sci. Lett.*, **375**, 101-112.

Schurr, B., Ratschbacher, L., Sippl, C., Gloaguen, R., Yuan, X. & Mechie, J., 2014. Seismotectonics of the Pamir, *Tectonics*, **33**(8), 1501-1518.

Schuster, R.L., 2005. Usoi Landslide Dam and Lake Sarez, Pamir Mountains, Tajikistan, *Environ. Eng. Geosci.*, **10**(2), 151-168.

Schwab, M. et al., 2004. Assembly of the Pamirs: age and origin of magmatic belts from the southern Tien Shan to the southern Pamirs and their relation to Tibet, *Tectonics*, **23**(4), doi:10.1029/2003TC001583.

Schwartz, D.P., 2018. Review: past and future fault rupture lengths in seismic source characterization—the long and short of it, *Bull. seism. Soc. Am.*, **108**(5A), 2493-2520.

Schwartz, D.P., Haeussler, P.J., Seitz, G.G. & Dawson, T.E., 2012. Why the 2002 Denali fault rupture propagated onto the Totschunda fault: implications for fault branching and seismic hazards, *J. geophys. Res.*, **117**(B11), B11304.

Searle, M.P. & Hacker, B.R., 2018. Structural and metamorphic evolution of the Karakoram and Pamir following India - Kohistan - Asia collision, in *Himalayan Tectonics: A Modern Synthesis*, Vol. **483**, eds Treloar, P.J. & Searle, M.P., Geological Society London Special Publications.

Sharp, R., 1982. The Imperial Valley, California, earthquake of October 15, 1979, *USGS Professional Paper 1254*, p. 451.

Sharp, R.V., Budding, K.E., Boatwright, J., Ader, M.J., Bonilla, M.G. & Clark, M.M., 1989. Surface faulting along the Superstition Hills fault zone and nearby faults associated with the earthquakes of 24 November 1987, *Bull. seism. Soc. Am.*, **79**(2), 252-281.

Shelef, E. & Oskin, M.E., 2010. Deformation processes adjacent to active faults: examples from eastern California, *J. geophys. Res.*, **115**(B5), B05308.

Sibson, R., Ghisetti, F. & Ristau, J., 2011. Stress control of an evolving strike-slip fault system during the 2010-2011 Canterbury, New Zealand, earthquake sequence, *Seismol. Res. Lett.*, **82**(6), 824-832.

Sibson, R.H., 1985. Stopping of earthquake ruptures at dilational fault jogs, *Nature*, **316**(6025), 248-251.
Sippl, C. et al., 2013. Geometry of the Pamir-Hindu Kush intermediate-depth earthquake zone from local seismic data, *J. geophys. Res.*, **118**(4), 1438–1457.

Stein, R.S., Barka, A.A. & Dietrich, J.T., 1997. Progressive failure on the North Anatolian fault since 1939 by earthquake stress triggering, *Geophys. J. Int.*, **128**(3), 594–604.

Stirling, M.W., Goded, T., Berryman, K.R. & Litchfield, N., 2013. Selection of earthquake scaling relationships for seismic-hazard analysis, *Bull. seism. Soc. Am.*, **103**(6), 2993–3011.

Strecker, M.R., Frisch, W., Hamburger, M.W., Ratschbacher, L., Semiletkin, S., Zamoruyev, A. & Sturchio, N., 1995. Quaternary deformation in the eastern Pamirs, Tadzhikistan and Kyrgyzstan, *Tectonics*, **14**(5), 1061–1079.

Strom, A., 2014. Sarez Lake problem: ensuring long-term safety, in *Landslide Science for a Safer Geoenvironment*, Vol. 1, pp. 1–486, eds Sassa, K., Canuti, P. & Yin, Y., Springer.

Strom, A. & Abdrakhmatov, K., 2018. *Rockslides and Rock Avalanches of Central Asia: Distribution, Morphology, and Internal Structure*, Elsevier.

Styron, R., Taylor, M.H. & Okoronkwo, K., 2010. Database of active structures from the Indo-Asian Collision, *Eos*, **91**(20), 181–182.

Stübner, K. et al., 2013. The giant Shakhdara migmatitic gneiss dome, Pamir, India-Asia collision zone: 2. Timing of dome formation, *Tectonics*, **32**(5), 1404–1431.

Taylor, M.H. & Yin, A., 2009. Active structures of the Himalayan-Tibetan orogen and their relationships to earthquake distribution, contemporary strain field, and Cenozoic volcanism, *Geosphere*, **5**(3), 199–214.

U.S. Department of the Interior, 2016. *Landsat 8 Data Users Handbook*, NASA, U.S. Geological Survey.

U.S. Geological Survey, Earthquake Hazards Program, 2017. Advanced National Seismic System (ANSS) Comprehensive Catalog of Earthquake Events and Products: Various, https://doi.org/10.5066/F7MS3QZH, https://www.sciencebase.gov/catalog/item/52eab950e4b0444d1ce67917.

Vallage, A., Klinger, Y., Grandin, R., Bhat, H.S. & Pierrot-Deseseilligny, M., 2015. Inelastic surface deformation during the 2013 Mw 7.7 Balochistan, Pakistan, earthquake, *Geology*, **43**(12), 1079–1082.

Verma, R.K. & Chandra Sekhar, C., 1985. Seismotectonics and focal mechanisms of earthquake from Pamir-Hindukush regions, *Tectonophysics*, **112**, 297–324.

Walker, R.T., Bergman, E.A., Szeliga, W. & Fielding, E.J., 2011. Insights into the 1968-1997 Dasht-e-Bayaz and Zirkuh earthquake sequences, eastern Iran, from calibrated relocations, InSAR and high-resolution satellite imagery, *Geophys. J. Int.*, **187**(3), 1577–1603.

Wegmüller, U. & Werner, C., 1997. Gamma SAR processor and interferometry software, in *Proceedings of the Third ERS Symposium on Space at the Service of our Environment*, Vol. 414, p. 1687, ESA Special Publication.

Wells, D. & Coppersmith, K.J., 1994. New empirical relationships among magnitude, rupture length, rupture width, rupture area, and surface displacement, *Bull. seism. Soc. Am.*, **84**(4), 974–1002.

Wesnousky, S.G., 2006. Predicting the endpoints of earthquake ruptures, *Nature*, **444**(7117), 358–360.

Winsvold, S.H., Kääb, A. & Nuth, C., 2016. Regional glacier mapping using optical satellite data time series, *App. Earth Observ. Remote Sens.*, **9**(8), 867–890.

Wright, T.J., Parsons, B.E., Jackson, J.A., Haynes, M., Fielding, E.J., England, P.C. & Clarke, P.J., 1999. Source parameters of the 1 October 1995 Dinar (Turkey) earthquake from SAR interferometry and seismic bodywave modelling, *Earth planet. Sci. Lett.*, **172**(1-2), 23–37.

Yin, A., 2010. Cenozoic tectonic evolution of Asia: a preliminary synthesis, *Tectonophysics*, **488**(1–4), 293–325.

Zachariasen, J. & Sieh, K.E., 1995. The transfer of slip between 2 en echelon strike-slip faults - a case-study from the 1992 Landers Earthquake, Southern California, *J. geophys. Res.*, **100**(B8), 15281–15301.

Zhang, J., Shan, X. & Huang, X., 2011. Seismotectonics in the Pamir: an oblique transpressional shear and south-directed deep-subduction model, *Geosci. Front.*, **2**(1), 1–15.

Zinke, R., Hollingsworth, J. & Dolan, J.F., 2014. Surface slip and off-fault deformation patterns in the 2013 Mw 7.7 Balochistan, Pakistan earthquake:

implications for controls on the distribution of near-surface coseismic slip, *Geochem. Geophys. Geosyst.*, **15**, 1–17.

Zubovich, A., Schöne, T., Metzger, S., Mosienko, O.I., Mukhamediev, S., Sharshebaev, A. & Zech, C., 2016. Tectonic interaction between the Pamir and Tien Shan observed by GPS, *Tectonics*, **35**(2), 283–292.

Zubovich, A.V. et al., 2010. GPS velocity field for the Tien Shan and surrounding regions, *Tectonics*, **29**(6), 1–23.

## SUPPORTING INFORMATION

Supplementary data are available at GJI online.

**Figure S1** Correlation results from individual year-spanning scene pairs. Each pair of panels ('northward' and 'eastward') shows the respective orthogonal components of the displacement fields for a single pair of pre- and post-earthquake Landsat-8 panchromatic images captured on nearly the same Julian day of each year. Labels on the left indicate the Julian day of the post-event image & the Julian day of the pre-event image in the year before. The simplified trace of the observed SKF rupture is shown on each in red. International borders are faint black lines. Note (a) the differing spatial distribution of no-data among different scene combinations, resulting usually from cloud cover; (b) the differing gradients and other long-wavelength artefacts in each pair, indicating rotations, undulations and higher order errors in geographic registration within the images. Taking the median of all these results helps remove these non-systematic errors.

**Figure S2** Swath profiles used to characterize surface offsets along the trace of the SKF. Each panel shows a 1-km-wide swath profile through the fault-parallel (i.e. positive in the 034° azimuth) displacement field, oriented perpendicular to the simplified trace of the SKF rupture, as indicated in Fig. 8 in the main text. Points represent individual pixel values, coloured by their position across the swath. Thin black curve represents median value along swath at pixel resolution. The grey bar represents the zone within which image correlation windows include pixels on both sides of the fault, depressing the absolute value of displacement and thus contaminating the discrete offsets. We exclude that portion from our calculations of offset. Additional colored lines represent various values of offset measured using different approaches, as labeled in Fig. 7 in the main text. We derive a near-field slip distribution by comparing displacements within 1 km on either side of the fault (but excluding the zone <250 m from the fault). We quantify fault offset following a variety of common methods: differencing the intercepts of linear fits on either side (red bar), differencing maximum and minimum values of either side (blue bar), fitting an arctangent function (black dashed line) and differencing the median values within 1 km (orange horizontal bars). Our preferred slip distribution, presented in main text Fig. 8 and the accompanying supplemental spreadsheet represents the difference between median values of the swath within the 250–1000 m range on either side of the fault (represented as horizontal orange bars in this figure), for reasons explained in Fig. S3.

**Figure S3** Comparison of offsets derived in different manners. We measured fault offset from the Landsat-8 horizontal displacement field using four different common approaches in order to test and/or illustrate how these differently derived values represent fault slip. The graphs in this figure compare the values derived by each method against each other. The four methods are: (1) differencing median displacement on either side of the fault, (2) differencing the intercepts of lines fit to either side of the fault, (3) differencing extreme values of displacement on either side of the fault and (4) fitting an arctangent function to both sides of the fault. The relatively common approach of fitting lines to the profile on either side of the
fault and projecting them to the trace to measure their offset proves scattered and unreliable, as line fitting is sensitive to excursions of the displacement field within the window we choose, and results in a large number of negative (i.e. apparently right-lateral) offsets. We thus find this approach is unsuitable for uniformly determining slip measurements along the fault. Taking extreme values similarly enhances outliers that don't always appear representative of the overall offset observed, and of course they skew the entire set of offset measurements higher. The method we prefer is differencing the median values of displacement within the specified range (250–1000 m) on either side of the fault. These values should best reflect near-field displacement while reducing the noise from spuriously high pixels. Similar results are obtained when fitting an arctangent function to both sides of the fault. These latter methods result in only 4 profiles producing negative offsets (swaths at km 56–59; see Fig. S2) because of a local negative anomaly on the east side of the fault (visible & highlighted in Fig. 8 of the main text). By contrast, differencing linear models of either side of the fault results in more than 10 per cent of the profiles producing negative offsets, contrary to the consistent left-lateral slip sense we observe overall.

**Figure S4** Sentinel-1 wrapped interferograms of the Murghob earthquake showing observed line-of-sight (LOS) phase shift in centimetres, modelled LOS phase shift based on the inverted fault slip model and residual misfit between the two. Arrows indicate the satellite orbit direction (long) and look direction of the SAR (short) for the descending (track 5) and ascending (track 100) geometries. Black line represents simple fault used for initial uniformlipslip inversion; grey line represents 9-segment fault used for final slip-varying inversion. To derive these interferograms we use Sentinel-1 Precise Orbit Ephemerides data files (available from https://qc.sentinell.eo.esa.int/) to extract state vectors for the SAR images. We process the Single-Look-Complex (SLC) data in TOPS mode to combine bursts and use all three sub-swaths. To coregister the slave to the master image, we estimate the SLC offsets with intensity cross-correlation. We then refine the azimuth offset estimation with an enhanced spectral diversity method in the burst overlap regions to remove phase discontinuity errors at burst edges due to misalignment in the azimuth direction. To simulate the interferometric phase due to topography we use the SRTM 30 m digital elevation model. We generate the differential interferogram pair, multilooking with 18 looks in range and 3 looks in azimuth (~42 m) to improve coherence coverage (Figs 3 and S4). In order to aid unwrapping we smooth the phase using an iterative adaptive spectral filtering method (Goldstein & Werner, GRL 1998). Next

we unwrap the interferometric phase to give relative line-of-sight displacements using minimum cost flow and triangulation. Finally, we geocode the unwrapped phase from radar to geodetic coordinates, before downsampling the data to 100 m (Figs 3 and S4) and converting to UTM zone 43N.

**Figure S5** Input InSAR data, synthetic model and residuals from initial inversion for fault geometry. Fault dip and rake are solved for using the smallest number of planar faults with uniform slip that fairly represent observed surface discontinuity, coarsely follow the preexisting trace of the SKF, and are in rough agreement with previous source models for this earthquake. Track 100a is the ascending track of Sentinel-1 radar data; track 005d is the descending track of Sentinel-1 radar data; 'track 888d' represents the horizontal deformation detected from Landsat-8 image correlation.

**Figure S6** Input InSAR data, synthetic model, and residuals from inversion for spatially varying slip along a 9-segment fault model compatible with the full resolution discontinuity in Landsat offsets. Fault section dips and rakes are based on those determined in unform-slip/ simple fault trace inversion. This figure shows the inversion using InSAR data only.

**Figure S7** Same as Fig. S6 but with Landsat horizontal displacement field and InSAR LOS displacement fields given equal (50 per cent) weight.

**Figure S8** Same as Figs S6 and S7 but for optimum weighting between Landsat and InSAR displacement data (see Fig. S9 for optimization of weighting).

**Figure S9** RMS misfit across entire inversion model, used for optimization of data weighting between Landsat displacement field and InSAR LOS displacements.

**Figure S10** Upper panel: Distribution of 2015 slip along the SKF from the inversion of the InSAR surface displacements and Landsat-8 offsets. Black crosses denote slip patches with magnitudes less than two times the uncertainty (indicated in the lower panel). Grey lines denote the boundaries between planar segments (3-D geometric relationships shown in main text Fig. 6). Pink star indicates USGS hypocentre (projected onto fault plane) with uncertainty ellipse. Lower panel: estimates of slip uncertainty from a Monte Carlo approach of calculating the standard deciation of 100 slip inversions using perturbed data sets based upon the far-field noise characterized in the InSAR and Landsat-8 data.

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.