(2019). The decade-long Machaze-Zinave aftershock sequence in the slowly straining Mozambique Rift. Geophysical (1),

be considered if a large earthquake were to occur in the southern East African Rift system today.


I N T RO D U C T I O N A N D B A C KG RO U N D
The continental East African Rift System (EARS) is the surface expression of the separation of the Nubian and Somalian tectonic plates. It can be followed southwards from Ethiopia to Malawi, south of which extension rates are low (<2 mm yr −1 ; Stamps et al. 2008;Saria et al. 2014), and the focus of extension is less distinct (Chorowicz 2005;Bird et al. 2006). In the south of Malawi two branches of the EARS have been hypothesized based on the distribution of seismicity and faulting. These two branches are (i) the Luangwa Rift, trending northeast-southwest through Zimbabwe and Botswana (e.g. Scholz et al. 1976;Modisi 2000;Kinabo et al. 2007), and (ii) the Mozambique Rift, which continues southwards (e.g. Fairhead & Henderson 1977;Fonseca et al. 2014;Domingues et al. 2016;Stamps et al. 2018;Fig. 1a). Faulting in the Mozambique Rift may therefore represent the southernmost and least mature portion of the continental EARS. A band of seismicity, whose trend mirrors that of the Neoproterozoic Mozambique Mobile Belt and Karoo age tectonic structures, occurs from the Cainozoic Urema Graben in the north to the Machaze district (administrative region between 22.2 • and 21.4 • S, bounded in the south by the Save River) of southern Mozambique (Forster 1975;Fonseca et al. 2014;Fig. 1a). In the south of the Machaze district, there is little moderate or high magnitude seismicity (>M w 3) recorded in the USGS catalogue, possibly suggesting a change in strain rate or crustal rheology. Understanding how and where the EARS extends through this region is important for understanding the development of continental rifts.
In rift environments, pre-existing structures exhibit a control on the distribution of faulting (e.g. Milani & Davison 1988;Versfelt & Rosendahl 1989;Kinabo et al. 2008;Muirhead & Kattenhorn 2018). This is because the frictional strength of faults, that is, the differential stress required for slip on pre-existing faults, can be lower than the stress required for fault formation (Scholz 2002). Faults therefore represent planes of weakness and may be reactivated in stress conditions that are not orientated optimally for failure. This phenomenon has implications for incipient rifting, as large pre-existing structures can control the stress distribution and geometry of a rift (e.g. Versfelt & Rosendahl 1989; The Machaze-Zinave aftershock sequence  (Fonseca et al. 2014;Reeves et al. 2016). (b) The red, green and blue boxes show the cropped extents of Sentinel-1 tracks T79, T72 and T174 that are used in this study, respectively. The USGS focal mechanisms for the Machaze (M) and Zinave (Z) earthquakes are shown, with the earthquakes that occurred since the Machaze event [coloured by time, see sub-panel (c) (circles show USGS solutions, squares are events observed by Fonseca et al. (2014)]. The black line shows the location of the Machaze fault from Copley et al. (2012). (c) Magnitude-time distribution of earthquakes in the Machaze-Zinave region. The dashed lines denote the observation period of Fonseca et al. (2014) 2008). In some regions, however, pre-existing structures are absent or not suitably orientated for reactivation. In these cases, extension will not be fully accommodated through fault reactivation and new structures may form.
Coseismic slip induces stresses on the surrounding rock Stein et al. 1994) that can promote or inhibit failure, either immediately or after a period of time (Tuttle et al. 2002;Nostro et al. 2005). Indeed, fault interaction plays an important role in influencing earthquake recurrence and fault growth in extensional settings (e.g. Nicol et al. 2010;Hodge et al. 2018b). In the 2009 Karonga earthquake sequence, Malawi, each new earthquake occurred within a region of positive static stress change from the previous event. This process of segmented fault rupture allowed for the lateral transition of slip on faults over ∼40 km, while being confined to the upper ∼10 km (Biggs et al. 2010;Fagereng 2013).
By studying incipient rifts we can observe rift development processes in action and aim to understand the growth and interaction of faults (Scholz et al. 1993;Walsh et al. 2002;Hodge et al. 2018b). Understanding how faults within a system interact is important in determining the induced or reduced hazard an earthquake might represent for its vicinity (e.g. Lin & Stein 2004;Wedmore et al. 2017). This is especially true in regions with low strain rates, which can produce long aftershock sequences (Stein & Liu 2009), elevating the hazard for up to hundreds of years after the main shock.

The Mozambique rift
On 2006 February 22, an M w 7.0 normal faulting earthquake occurred in the Machaze region of Mozambique, in the southern East African Rift (Fenton & Bommer 2006;Yang & Chen 2008;Copley et al. 2012). This normal-faulting earthquake was one of the largest earthquakes in continental Africa for a century and demonstrates that the region is subject to extensional stresses. On 2016 September 22, 10.5 yr after and <10 km from the Machaze earthquake, an M w 5.6 normal faulting earthquake occurred in the Zinave National Park (USGS; http://earthquake.usgs.gov).
The Machaze-Zinave region is at or near the intersection of three major inherited structures that obliquely intersect ( Fig. 1a): (i) the ∼E-W Limpopo Belt, hosting the Okavango and Limpopo Dyke Swarms (∼180 Ma and 728 ± 3-1683 ± 18 Ma; Le Gall et al. 2002;Jourdan et al. 2006), (ii) the ∼NE-SW trending Urema Graben north of Machaze (Steinbruch 2010) (iii) and the ∼NNW-SSE Mazenga Graben to the south (Fig. 1a), which runs parallel to the eastern limit of the Kaapvaal Craton (Domingues et al. 2016). The Kaapvaal and Zimbabwe Cratons are south-west and north-west of the Machaze region, respectively. This region of the Mozambique coastal planes is thought to have a 5-10 km thick sedimentary sequence of post-Jurassic age overlying a Precambrian crystalline basement (Gwavava et al. 1992;Salman & Abdula 1995).
We chose to investigate the Zinave earthquake using Interferometric Synthetic Aperture Radar (InSAR) and seismology, together with a reanalysis of InSAR for the Machaze event, to better understand the tectonics of southern Mozambique. The region is poorly instrumented and difficult to access, and so satellite geodesy provides an ideal tool to investigate faulting. In this work we investigate the spatial and temporal relationships between the 2006 Machaze and 2016 Zinave earthquakes. The temporal evolution of seismicity following the Machaze earthquake is analysed, and we conclude that the Zinave earthquake is part of an aftershock sequence.
We process ENVISAT interferograms to investigate the surface deformation of the Machaze event and three overlapping tracks of Sentinel-1 data for the Zinave earthquake. We then invert the observations of surface deformation for earthquake parameters using uniform and distributed slip models. For the Zinave earthquake, we compare the geodetic models to seismological estimates of the source parameters from an inversion of body wave waveforms. The Coulomb stress transfer following the Machaze earthquake is then calculated to investigate whether the Zinave fault was brought closer to failure following the Machaze main shock. We discuss the role that the subsurface frictional properties have in the assessment of hazard, and the style of slip that could occur in this region. We relate the occurrence of these earthquakes to fault growth, the role pre-existing structures have on controlling incipient rifting and the evolution of the East African Rift.

S E I S M I C I T Y I N T H E M O Z A M B I Q U E R I F T
The Machaze earthquake occurred on an ∼N-S striking fault that dipped unusually steeply (70 • -75 • ) to the west (Fenton & Bommer 2006;Copley et al. 2012;Attanayake & Fonseca 2016). The fault plane is thought to be a reactivated pre-existing structure within the crystalline basement (Yang & Chen 2008;Copley et al. 2012). The earthquake surface rupture was >15 km in length with up to 2 m of vertical offset (Fenton & Bommer 2006). Models based on teleseismic waveforms and ENVISAT interferograms show that coseismic slip extended down to ∼25 km depth with diminishing amplitude between 10 km and the surface (Copley et al. 2012). This reduction of slip is interpreted to be a result of the juxtaposition of velocity weakening (crystalline basement) and velocity strengthening (sedimentary) lithologies (Copley et al. 2012). The shallow slip deficit was later partially recovered post-seismically through afterslip in the upper 10 km (Copley et al. 2012). The control of pre-existing structures and upper-crustal rheology on slip during the earthquake cycle has been similarly observed elsewhere in the EARS (e.g. Versfelt & Rosendahl 1989;Ring 1994;Kolawole et al. 2018).
From 1990 to 2006 there was no moderate magnitude (e.g. >M w 4.5) seismicity observed in the Machaze region, but between the Machaze event and 2018, 134 magnitude 3 or greater earthquakes were recorded (USGS catalogue; Figs 1b and c). Between 2011 and 2013 there was a deployment of 30 seismometers throughout southern Mozambique (Fonseca et al. 2014). This deployment observed persistent seismicity with 143 earthquakes greater than a magnitude 0.9 in the Machaze region, indicating that the area remained seismically active (Fig. 1).

Aftershock sequences
Static and dynamic stress changes from large earthquakes cause aftershocks, which can pose a significant seismic hazard in their own right (Wiemer et al. 2002;Parsons et al. 2008;Stein & Liu 2009). We test whether the Zinave earthquake belongs to an aftershock sequence following the Machaze earthquake or is an independent event. The empirical Omori law (Omori 1894) states that the rate of aftershocks is inversely proportional to the time since the main shock. The modified version of this relationship (such that p is not equivalent to 1) is given by eq. (1) (Utsu 1961;Dieterich 1994). (1) In the modified Omori law K and c are constants that describe the productivity and the time delay of the sequence (a function of catalogue completeness), respectively. The parameter p is related to the physical heterogeneity of the aftershock region and can be time dependent (e.g. Kisslinger & Jones 1991;Helmstetter & Shaw 2006). In our analysis, we fix p to be equal to 1 for simplicity and do not consider the effects of aftershocks of aftershocks. The Omori law fits the temporal evolution of the recorded seismicity since the Machaze earthquake well with K = 13.4 and c = 0.54 (Fig. 1e). This observation suggests that the Zinave earthquake is an aftershock that occurred 10.5 yr after the Machaze earthquake. The high c value for the sequence indicates that the seismic record is likely missing events early in the sequence, which introduces uncertainty in other model parameters and supports our use of p = 1.

Aftershock duration in low strain rate regions
The duration of an aftershock sequence is observed to be inversely proportional to the strain rate such that regions with low strain rates have long aftershock sequences (Stein & Liu 2009). This inverse relationship is consistent with a rate-and-state description of fault friction, where the aftershock duration (t a ) is a function of the normal stress (σ n ) and the rate of shear stressing across a fault. The rate of shear stress can be approximated using the relative velocity across a fault with a simple geometry to give eq. (2) (Savage & Burford 1973), where A is a constitutive parameter (Dieterich 1994), w is the fault width, μ is the rigidity and v is the relative velocity across the fault.
Using typical values, Stein & Liu (2009) determine an expected aftershock duration (in years) of 314/v (for v in mm yr −1 ). In southern Mozambique, where v 2 mm yr −1 , we may therefore expect aftershocks, and thus an significant seismic hazard, to last up to ∼150 yr (Goda 2012;Hodge et al. 2015). The 2016 September 22 M w 5.6 Zinave earthquake occurred at 20:06:11 UTC. We use observations of teleseismic body waves and the surface deformation from InSAR to investigate the fault geometry and earthquake kinematics. We invert the InSAR data for coseismic uniform and distributed slip models.

Body wave inversion
We jointly inverted P and SH waveforms to obtain the strike, dip, rake, centroid depth and source time function of the Zinave event.
Seismograms recorded at epicentral distances of 30 • -80 • were bandpass filtered at 15-100 s. These measures remove complexities due to lithospheric reverberations, interactions with the core and short-wavelength complexity in the source process or velocity structure. We can then model the earthquake as a point source, using the MT5 program of Zwick et al. (1994) (based on the algorithm of McCaffrey et al. 1991 andAbers 1988). This procedure is now routine and its application in Africa is described in detail in Craig et al. (2011). P waveforms were weighted 2:1 against the SH waveforms to account for their lower amplitudes, and all seismograms were weighted according to azimuthal coverage. The velocity structure at the source was specified as a simple half-space velocity model, with V p = 3.9 km s −1 and V s = 2.2 km s −1 , which is based on the crustal parameters used in the InSAR inversions (Section 3.3). The velocity model is relatively slow to reflect the thick sedimentary deposits.
Our best-fitting model is shown in Fig. 2. Body wave modelling indicates that, if the causative fault dipped west (as implied by the InSAR analysis below), then the earthquake occurred on a 175 • striking normal fault, with a dip of 66 • , rake of −69 • , depth of 7 km and moment of 9.5 × 10 16 Nm. The seismic estimate of the strike of this earthquake is comparable to estimates of the Machaze earthquake (175 • ; Yang & Chen 2008; 172 • CMT), but the dip is ∼10 • shallower. The resolution of the model is hampered by the lack of clear waveforms at stations in the Atlantic and Indian oceans and in Antarctica. The errors in our parameters for this event are typical for this technique and are strike ±10 • , dip ±5 • , rake ±10 • and centroid depth ±4 km.

Sentinel-1 InSAR
To investigate the surface deformation associated with the Zinave earthquake we produced 79 Sentinel-1 interferograms, 34 of which are coseismic, from three tracks (Table A1). Interferograms were processed using the GAMMA software (Werner et al. 2000) within the LiCSAR facility (González et al. 2016). We used the 30 m SRTM DEM (Farr & Kobrick 2000) to remove topographic phase contributions and down-sampled interferograms to a 100 m final pixel size. Each interferogram was filtered using a power-law filter (strength 0.85; Goldstein & Werner 1998) and pixels with coherence values less than 0.8 were masked out.
In tropical regions, atmospheric water vapour can cause artefacts in interferograms, particularly in regions of high topography (e.g. Webley et al. 2004;Ebmeier et al. 2013;Parker et al. 2015). We see signs of turbulent atmospheric delays in some interferograms and test using the high-resolution European Centre for Medium-Range Weather Forecasts (HRES ECMWF) weather model through the GACOS facility to correct for them (Yu et al. 2017b,a). However, the model was unable to account for the observed delays with the corrections increasing the standard deviation in our data (Fig. A1).
Weather models are best suited to predicting the stratified component of atmospheric delays often related to topography. In this region of southern Mozambique there is a maximum of ∼100 m elevation change over ∼100 km 2 . As such, we do not expect topographycontrolled atmospheric delays to be a large source of noise. We use a pair-wise logic approach to identify atmospheric delays and acquisitions with strong delays were removed (e.g. Massonnet & Feigl 1995;Ebmeier et al. 2013). To increase the signal-to-noise ratio, we stacked coseismic interferograms from each track, ensuring each acquisition contributed equally (Table A1, Fig. A2). Each stack was then de-ramped to remove long-wavelength atmospheric, orbital or ionospheric delays.
The Sentinel-1 interferograms show an approximately northsouth trending, ∼20 km by ∼10 km region of up to ∼5 cm positive line-of-sight range change (away from the satellite), beside a smaller region of ∼1 cm negative line-of-sight range change to the east. The boundary between these two regions is relatively sharp, and the pattern is broadly consistent with that caused by an approximately north-south striking normal fault. All Sentinel-1 interferograms were highly coherent, except for coseismic ones in the river bed directly above the fault (Fig. A2). We suggest this loss of coherence is caused by liquefaction, a phenomenon that was extensive during the Machaze event (López-Querol et al. 2007). Chains of post-seismic interferograms for each track, covering up to 3 months after the earthquake, show no significant post-seismic deformation (Fig. A3).
The distribution of slip during an earthquake is never uniform (e.g. Reilinger et al. 2000;Simons et al. 2011;Sangha et al. 2017). However, the inversion for fault geometry with variable slip is nonlinear and computationally expensive. As such, we initially perform a nonlinear inversion for the fault geometry, with uniform slip, and then use this fixed geometry to linearly invert for the distribution of slip on the fault plane (e.g. Pedersen et al. 2003;Walters et al. 2009;Bie et al. 2017).

Zinave earthquake uniform slip modelling
We use the analytical solution for slip on a rectangular dislocation in an elastic half-space (Okada 1985) to model the surface displacements associated with the Zinave earthquake. We invert the observations of line-of-sight surface displacement using a Bayesian approach, incorporating the Markov chain Monte Carlo Metropolis-Hastings algorithm to explore the parameter space (Hastings 1970;Mosegaard & Tarantola 1995;González et al. 2015), and report posterior probability density functions for each model parameter. We assume all measurement errors are Gaussian. Each interferogram stack is subsampled prior to inversion, based on the interferogram variance, using the quad-tree approach (Jonsson et al. 2002). Initial conditions and bounds for the fault model parameters are given in Table A2. We use a Poisson's ratio of 0.25 and a shear modulus of 10 GPa for the crust to account for the thick sedimentary sequence.
The location of the Zinave earthquake is covered by three Sentinel-1 tracks (two ascending, tracks 72 and 174, and one descending, track 79). We perform five separate inversions using different subsets of the data to test the robustness of our solutions, and the influence of additional data sets on model parameter certainty, as each geodetic data set has different noise and sensitivity (Table A1, Fig. 3). All of our inversions are consistent with a northsouth trending, westward dipping normal fault (Tables 1 and A3). To compare the results of each inversion, we calculate the standard deviation of the 2σ confidence bounds of each model parameter, 508 R. Lloyd, J. Biggs and A. Copley Zinave, 22/9/2016 175°/66°/289°/7 km/9.544 x 10 16 Nm 5741d 0 60s assuming the distributions are Gaussian ( Table 2). Results of inversions using individual interferograms (from T174, T79 and T72) have the largest standard deviations, indicating they are least well constrained. The overall best constrained model is that from the inversion of InSAR data from all three tracks, with between a 38 and 79 per cent decrease in the standard deviation of the parameter confidence bounds ( Table 2).
The best-fitting model is for a fault that is 16 km long (95 per cent confidence: 14.8-16.7 km), with a dip of 59 • (54 • -67 • ), bottom depth 7.6 km (5.6-8.0 km) and downdip width of 5.1 km (1.4-5.7 km). All of the modelled fault parameters, 95 per cent threshold values and root-mean-square misfit values are given in Tables 1 and A3. Overall the model is able to describe the observed surface deformation well; however, there is a small residual near the southern end of the fault (Figs 4c, f and i). The spatial pattern of this residual is consistent between all of the tracks, and so is likely a feature of the rupture that is not captured by the uniform slip model.

Zinave earthquake distributed slip modelling
The model with uniform slip can reproduce the first-order patterns of surface deformation caused by the Zinave earthquake. However, systematic residuals suggest the slip was greater or shallower at the southern end of the fault. To test this, we refine our model by inverting for variable slip using the best combination of interferograms found in the previous section (stacks from ascending tracks 72 and 174, and descending track 79).
We use the plane defined in Section 3.3 (striking 168 • and dipping 59 • ), but enlarge it to extend between the surface and a depth of    10 km, with a length of 24 km (Fig. 5). This plane is then divided into 100 (10 by 10) smaller patches, each of which is 2.4 km horizontal by 1 km vertical in size (see Fig. A4 for comparison to model resolution). For a fixed fault geometry, the slip on each patch can be linearly related to the observations of surface deformation (eq. 3; Wright et al. 2004;Funning et al. 2005). Eq.
(3) relates the surface displacement observation at each pixel location (described by X and Y) (d) to the slip on each fault patch (m) via the Green's function (G) of line-of-sight displacements from 1 m of slip on each patch, using the elastic dislocation model (Okada 1985). The formulation includes ∇ 2 , which is the Laplacian smoothing operator to avoid sharp slip variations (Jonsson et al. 2002), and κ 2 which is a prescribed scalar smoothing factor (Funning et al. 2005;Biggs et al. 2006). Eq. (3) is solved using a non-negative least-squares algorithm (Bro & De Jong 1997). In all inversions κ = 5 × 10 5 , which was selected from a trade-off curve to minimize the roughness without oversmoothing (Fig. A5a). A ramp (described by a first-order polynomial with constants a and b in the x and y directions, respectively) and an static offset (c) for each interferogram stack is also solved for in the inversion to account for long wavelength noise contributions (Biggs et al. 2007).
We investigate whether the rake varies across the slip region by performing three inversions: one with variable rake (G is formulated such that we can solve for variable rake, that is, non-negative slip in orthogonal components, between −45 • and −135 • , separately for each patch), one with the rake fixed from the geodetic uniform slip distribution model (−66 • ), and one with a rake fixed to the body wave inversion result (−71 • , Fig. A6). We then compare the fit of the variable rake model considering the increase in model parameters using an F-test.  All three inversions show no significant difference in their rootmean-square misfit (Table A4). However, in the variable rake model, the rake within the main slip region varies unrealistically, between 30 • and 70 • over a distance of ∼6 km (Fig. A7). A fixed rake model is a nested variant of the variable rake model (i.e. both models are described by the same equation, but the fixed rake model has fewer degrees of freedom), and as such we can perform an F-test to test whether the inclusion of variable rake is justified. The F-statistic value for this test is 0.005 considerably lower than the 5 per cent probability threshold value of 1.32. The increase in degrees of freedom for the variable rake model is therefore not statistically justified and rejected.
The slip distribution for the Zinave earthquake using a fixed rake of −66 • , and fit to the data are shown in Fig. 5. In comparison to the uniform slip model, there is a greater area and magnitude of slip (up to 0.12 m) at the southern end of the fault compared to the north (∼0.1 m), as suggested by the residuals to the uniform slip model (Fig. 5). The region with >2 cm of slip is ∼14 km by ∼5 km with top and bottom boundaries at 4 km and 8 km depth, respectively. For this slip distribution, the seismic moment, using a shear modulus of 32 GPa, is 1.85 × 10 17 Nm.
The geodetically and seismically derived values of dip, depth and rake for this earthquake agree within the error on both methods (Fig. 3). The seismically derived seismic moment for this event, 9.5 × 10 16 Nm, is 47 per cent smaller than the geodetic estimate of the moment (1.85 × 10 17 Nm). The seismic moment of an earthquake is difficult to determine accurately with both methods, and this degree of mismatch is typical (e.g. Weston et al. 2011).

ENVISAT InSAR
We processed 2 ENVISAT scenes (06/06/2004 and 07/05/2006) to produce an interferogram spanning the 2006 Machaze earthquake using the ROI-PAC software (Rosen et al. 2004). These scenes were selected following Copley et al. (2012), but we remove topographic phase contributions using a 30 m DEM (Farr & Kobrick 2000), which was not available in 2012. The interferogram was filtered with strength 0.9 (Goldstein & Werner 1998). We resample the interferogram to 90 m pixel spacing and use the HRES ECMWF atmospheric model through the GACOS facility to correct for atmospheric delays in the interferogram (Yu et al. 2017b,a). We find that the atmospheric correction reduces the standard deviation of the interferogram from 12.1 cm (90 m pixel, no correction) to 9.6 cm (90 m pixels, atmospheric correction applied).
The use of the 30 m DEM considerably extended the connected coherent region (coherence >0.1) of the ENVISAT interferogram compared to the interferogram of Copley et al. (2012) which used a 90 m DEM, especially in the far-field and some of the signal of the fault motion. However, we were still unable to unwrap the fault within ∼8 km due to poor coherence and/or high deformation rate.

SPOT cross-correlation
Surface deformation can be measured using the cross-correlation of optical images (e.g. Leprince et al. 2007Leprince et al. , 2008. To better constrain the near-field displacement field, we make use of two SPOT5 images, acquired on 2001 August 3 and 2008 August 26 (Copley et al. 2012). We use the results of Copley et al. (2012) who cross-correlate the images using the Cosi-Corr programme (Leprince et al. 2007(Leprince et al. , 2008. Due to the non-nadir viewing geometry of the satellite, vertical displacements contribute to the apparent eastwest displacement (Copley et al. 2011). For the orientation of the Machaze earthquake and SPOT satellite viewing geometry, we calculate that the vertical displacement will destructively interfere with the east-west component. Furthermore, the subjectivity in choice of tie-points for the cross-correlation, required because of uncertain satellite orbits, reduces the reliability of the long-wavelength displacements (Copley et al. 2011(Copley et al. , 2012. The displacement discontinuities, however, are robust. The acquisition dates cover both coseismic and post-seismic (30 months) periods, so the relative contribution cannot be distinguished. Nonetheless, the SPOT data set provides crucial observations of the near-field deformation (Fig. 6).

Uniform slip modelling of the Machaze earthquake using only InSAR data
In order to compare the Machaze and Zinave earthquakes, we perform inversions of observations of the Machaze earthquake using the same methodology as above. For the Machaze earthquake, we perform three inversions of the 06/06/2004-07/05/2006 ENVISAT interferogram: one with the GACOS atmospheric correction (90 m resolution; Figs 6a-c) and two without atmospheric corrections, one full resolution (30 m pixels) (Figs A8 a-c) and one resampled to 90 m (Figs A8 d-f).
We use the same inversion methodology and crustal properties as for the Zinave earthquake. Initial conditions and inversion bounds are given in Table A2. The results of all three of our inversions show that the far-field displacements captured by the ENVISAT interferograms can be well described (2.1-2.4 cm rms) by a single 24 × 7 km steeply dipping (78 • ) fault, striking 176 • , with top depth ∼7 km below the surface (Fig. 6, 95 per cent probability values can be found in Table 3). The inversion using the interferogram with corrections for atmospheric delays is our preferred model (Table 3). The best-fitting model is a fault with 6.3 m slip at a rake of −80 • . A top depth of 7 km contrasts with observations of faulting at the surface (Fenton & Bommer 2006); however, this is unsurprising as the ENVISAT data used in our inversion is incoherent close to the fault trace, resulting in little sensitivity to shallow slip. As such, we repeat the inversion including the SPOT data (Section 4.4). Our modelling results are otherwise consistent with previous seismological and geodetic estimates of the fault geometry and slip (Table 3; Copley et al. 2012;Fonseca et al. 2014).

Distributed slip modelling of the Machaze earthquake, incorporating InSAR and SPOT data
Large faults are generally complex with multiple segments (e.g. Nabelek 1985;Fletcher et al. 2014;Hamling et al. 2016). Models of the coseismic surface deformation from InSAR observations for the Machaze earthquake suggest rupture across a change in fault strike (Copley et al. 2012). In order to ascertain the total slip on the Machaze fault, we perform a distributed slip model with the fault location identified by displacement discontinuities in post-seismic deformation and the SPOT data (Copley et al. 2012). We otherwise use the same method as in Section 3.4, but with the fault divided into 3 km by 3 km patches that extends down to 30 km dipping at 75 • .
We perform two inversions, (1) using just the ENVISAT interferogram (κ = 7 x 10 5 , Fig. A5b), and (2) a joint inversion using the ENVISAT interferogram and SPOT observations (which are resampled to 50 m; κ = 1.3 x 10 6 , Fig. A5c), with InSAR and SPOT data weighted equally. We solve for a ramp in the ENVISAT data and an offset in all data sets. The difference between these two models is that model 1 contains primarily far-field coseismic deformation, but omits the near-field displacement, while model 2 includes the near-field observations, but also includes 30 months post-seismic deformation. Previous observations indicate that post-seismic deformation is primarily shallow afterslip on the fault (Copley et al. 2012).  The slip distribution of model 1 indicates that the majority of the slip is concentrated towards the northern end of the southern fault, with some displacement on the northern segment (Fig. 7). The displacement is primarily normal in the main rupture area with rake values ∼−80 • . In this model, however, the slip in the upper cells is small, in contrast to field observations of surface offsets of up to 2 m (Fenton & Bommer 2006), probably as a result of the lack of near-field observations.
The results of the joint inversion of ENVISAT and SPOT data shows the same general pattern of slip distribution as model 1, although the slip in the uppermost cell is now 2 m (Fig. 7, resolution  7. Slip distribution of model 1 (a, ENVISAT) and model 2 (b, ENVISAT and SPOT) for the Machaze earthquake. The estimated depth of the sedimentbasement interface is shown (Salman & Abdula 1995;Gwavava et al. 1992). Profiles to the right of subpanels (a) and (b) show the integrated moment release with depth.
shown in Fig. A9). The slip on the southern end of the northern fault has also increased by upto 1.5 to 4.5 m. This additional slip could be either (1) post-seismic deformation, (2) coseismic deformation that is not detected by the ENVISAT interferogram due to its limited coherence or (3) an artefact, as a result of the SPOT data image-processing methodology. The greatest misfit to the data are residuals in the ENVISAT observations towards the southern end of the southern fault, where the interferogram is coherent closest to the fault trace (Fig. 6i). This may be because the SPOT data includes surface deformation caused by afterslip that is modelled, but not detected by the ENVISAT observations. Our preferred model of the coseismic slip distribution is that which is based upon the ENVISAT observations of the displacement with a fixed fault geometry based on the SPOT discontinuities (model 1). Uncertainties as to the geometry and timing (co-versus post-seismic) of the displacement in the SPOT observations introduces difficulties in how to interpret the result of model 2, although they do broadly agree with model 1, and the additional shallow slip in model 2 is consistent with post-seismic deformation.

Comparisons to seismology
The CMT and USGS estimate the coseismic moment of the Machaze earthquake to be 4.5 × 10 19 and 4.6 × 10 19 Nm, respectively. Using a shear modulus of 32 GPa, we calculate a seismic moment of 4.4 × 10 19 Nm for our preferred model (model 1) and 5.3 × 10 19 Nm for model 2. These are 5 per cent smaller and 15 per cent larger than the seismic estimates, respectively, but within the range that may be expected given the different velocity structures and possible inclusion of post-seismic deformation in the geodetic observations (Weston et al. 2011).
In both the geodetic distributed slip models, the majority of the moment release occurred at depths of less than 15 km (model 1: 95 per cent, model 2: 94 per cent), with a geodetic centroid of 8 km. This depth estimate is shallower than the seismological estimates of the earthquake centroid depth of 15 ± 3 km (Yang & Chen 2008), 14.8 km (Attanayake & Fonseca 2016) and 13 km (Copley et al. 2012). The discrepancy is likely a result of the low depth resolution in our data due to the lack of near-field geodetic observations.
The depths of Machaze aftershocks suggest that this region has a seismogenic thickness of 20 km or greater (Yang & Chen 2010;Craig et al. 2011). Therefore, as the Machaze earthquake ruptured down to ∼20 km, it may represent an upper bound to the size of an earthquake in this region.

Coulomb stress change
Earthquakes cause static stress changes in the surrounding rocks, bringing nearby regions closer to, or further from, failure. In theory, for failure to occur the coseismic static stress change plus any pre-existing differential stresses must be greater than a threshold Coulomb failure criterion Toda et al. 2005). The Coulomb stress change ( σ c ) is a function of the normal ( σ n ) and shear stresses ( σ t ) on a fault plane and is defined as where μ is the effective coefficient of friction that incorporates the effects of pore fluid pressure. In this formulation positive shear stress is in the direction of slip and positive normal stresses are unclamping. We use the Coulomb 3.3 software package to calculate the σ c following the Machaze earthquake (Toda et al. 2011;Lin & Stein 2004) using slip on a rectangular dislocation model in an elastic half-space (Okada 1985). We use μ = 0.4 in our calculations and confirm that the results show little sensitivity within the range μ = 0.2-0.6. To understand whether the Machaze earthquake could have brought the Zinave fault closer to failure, we calculate both the Coulomb and normal coseismic stress changes, and test both uniform and distributed slip models of the Machaze earthquake (Fig. 8).
We calculate the Coulomb stress change at 5 km depth, comparable to that of the peak slip on the Zinave fault, and use the results from the uniform slip model for the Zinave fault to specify receiver fault geometry (Fig. 4, Table 1). Our preferred solution uses the distributed slip model based on the ENVISAT and SPOT data as it contains the most complete measure of the slip that occurred during and after the earthquake and will therefore give the best estimate of the overall stress change. Although the calculated stress changes close to the Machaze fault may be unreliable given the lack of near-field observations, they are likely to be relatively robust at distances comparable to the Machaze-Zinave separation (∼10 km). We find that the stress pattern does not vary greatly over the uncertainty range of the receiver fault parameters (strike: 166 • -169 • , dip: 54 • -67 • ) given by the posterior probability density functions of our uniform slip model.
The pattern of Coulomb stress change produced by the uniform and distributed slip models of the Machaze earthquake is broadly similar. For both, at 5 km depth, the stress change is asymmetric, with the largest normal and Coulomb stress changes occurring over ∼20 km to the east (<−1 MPa) and west (>1 MPa) of the Machaze fault (Fig. 8). Lobes of positive Coulomb stress change also extend over ∼15 km to the north and south of the fault. In all models the region of the Zinave earthquake experienced a localized ∼0.2 MPa increase in Coulomb stress. This Coulomb stress change is much lower than the stress drop of the event, indicating that the fault was dominantly releasing pre-existing shear stress. Positive Coulomb stress changes are more spatially heterogeneous for the distributed slip model than the uniform slip model, with positive changes not predicted as far away from the fault to the west.
The spatial pattern of aftershocks generally follow regions of positive Coulomb stress changes . The distribution of recorded Machaze aftershocks (2006)(2007)(2008)(2009)(2010) shows clusters in the regions of positive Coulomb stress change including near the Zinave fault (Fig. 8).

Rheological implications
Coseismic slip in the Machaze earthquake was largely contained within the crystalline basement, with post-seismic deformation, primarily afterslip, within the overlying sedimentary sequence. This observation led Copley et al. (2012) to conclude that the sediments were velocity strengthening, and the basement velocity weakening. In contrast, we find slip during the Zinave earthquake occurred primarily within the sedimentary sequence.
A transition between velocity weakening and velocity strengthening rheologies can be explained by (1) laterally variable material properties, (2) temporal variations in material properties or pore pressure or (3) stress-or time-dependent rheology. Spatial variations in frictional properties relating to lithology have been invoked to explain spatial variation in steady-state fault stability elsewhere [e.g. 2014 South Napa earthquake (Floyd et al. 2016) and the 2003 Tokachi-oki earthquake (Miyazaki et al. 2004)]. In the case of the Machaze-Zinave sequence, there are limited geological constraints on the subsurface but either depth variations in the sediment basement interface or physical variations within the sedimentary deposits are plausible. Indeed, shallow aftershocks observed by Copley et al. (2012) support the presence of stick-slip regions within the sedimentary sequence. In large subduction zone earthquakes aftershocks are similarly seen to concentrate in the areas dominantly deforming by aseismic creep and are generally thought to represent stuck asperities on a sliding fault surface (e.g. Igarashi et al. 2003;Bürgmann et al. 2005;Moreno et al. 2010). Alternatively, temporal changes, particularly variations in pore fluid pressure, can alter the effective friction (μ ) and thus control fault strength (Nur & Booker 1972;Copley 2017). This is a plausible mechanism in the Machaze region, which is a flood plain with highly seasonal rainfall, underlain by a thick sedimentary sequence of high permeability sands (López-Querol et al. 2007).
A steady-state, that is, rate-dependent, friction law assumes that afterslip occurs on stable faults in velocity strengthening regimes, and that the friction coefficient depends on the slip velocity only (e.g. Perfettini & Avouac 2004;Hsu et al. 2006). Within this framework, spatial heterogeneity would be required to explain afterslip in otherwise velocity-weakening zones. However, when the full rate-and-state law is considered, the distinction between velocity weakening and velocity strengthening is ambiguous (Helmstetter & Shaw 2009). When the steady-state constraint is not applied, fault behaviour is strongly influenced by post-seismic stresses, with both afterslip and earthquake nucleation possible within a homogenous medium.

Fault interaction and rifting
The steep dip of the Machaze earthquake fault (∼75 • ) and nucleation within the Precambrian basement suggests that the earthquake represents a reactivation of an existing structure, rather than the formation of a new fault (Yang & Chen 2008). The exploitation of pre-existing structures is common in nascent, and mature, rift settings (Versfelt & Rosendahl 1989;Ring 1994;Kinabo et al. 2008;Laó-Dávila et al. 2015;Domingues et al. 2016;Hodge et al. 2018a;Muirhead & Kattenhorn 2018). Pre-existing structures represent weaknesses, which, even if not optimally orientated, will fail preferentially to the formation of new faults (Scholz 2002). In contrast to Downloaded from https://academic.oup.com/gji/article-abstract/217/1/504/5298869 by University of Bristol Library user on 07 March 2019 the Machaze earthquake, slip during the Zinave earthquake and several aftershocks (February-April 2006;Yang & Chen 2010 occurred within the post-Jurassic sedimentary layer on an ∼60 • dipping plane as predicted by Andersonian mechanics. The Machaze-Zinave sequence raises several questions about fault growth in continental rift settings. The Zinave fault could (1) act as a linking fault between horizontally offset basement structures (e.g. Hodge et al. 2018b) and/or (2) represent the upward continuation of faulting on pre-existing basement structures through un-faulted sediments.
Reactivation of pre-existing structures has implications for the determination of the stress field in extensional environments using earthquake slip vectors (Delvaux & Sperner 2003;Bird et al. 2006;Saria et al. 2014). If, throughout rifting, large earthquakes occur on pre-existing structures they will not fully represent the present day principle stress orientations. Similarly, however, slip during aftershocks may represent a combination of local stresses from the mainshock and regional tectonics or interseismic strain accumulation. McKenzie (1969) demonstrated that, in a triaxial stress regime, slip vectors from shallow events provide limited constraint on the orientation of the greatest principle stress.

C O N C L U S I O N S
The Zinave earthquake occurred on a 60 • dipping normal fault between 4 and 8 km deep. The earthquake occurred in a region of positive Coulomb stress change associated with the 2006 Machaze event, indicating that it was brought closer to failure by the preceding earthquake. The depth of the Zinave earthquake suggests that it is contained within sedimentary deposits, at a depth co-incident with the coseismic slip deficit and post-seismic afterslip following the Machaze event. The occurrence of afterslip and aftershocks at the same depth suggests either spatially and/or temporally variable frictional properties, or that the fault exhibits time-or stress-dependent rheology.
A comparison to the modified Omori law for aftershock decay indicates that the Zinave earthquake is part of a prolonged aftershock sequence following the Machaze earthquake. Long aftershock sequences should be expected following large earthquakes in low-strain regions, suggesting that the seismic hazard in the least mature portions of the East African Rift is underestimated. The Machaze-Zinave sequence demonstrates that magnitude 4-5 earthquakes following magnitude >7 events should be expected for decades, with associated seismicity lasting for up to ∼150 yr.

A C K N O W L E D G E M E N T S
RL was supported by an NERC studentship tied to the LiCS (Looking inside the Continents from Space) consortium (NE/K010956/1). JB and AC are also supported by LiCS, as well as the NERC Centre for the Observation and Modelling of Earthquakes, Volcanoes and Tectonics (COMET). JB was supported by the Natural Environment Research Council (NERC) funded RiftVolc project (NE/L013932/1, Rift volcanism: past, present and future). Atmospheric corrections were made using the Generic Atmospheric Correction Online Service for InSAR (GACOS) facility (http://ceg-research.ncl.ac.uk/v2/ gacos/). The COMET Geodetic Bayesian Inversion Software (GBIS, Version 1.0) is available here: http://comet.nerc.ac.uk/gbis/. The Sentinel-1 and ENVISAT data used in this study are available through ESA. We would also like to thank Ake Fagereng, Maximilian Werner and John Elliott for useful discussions. Zwick, P., McCaffrey, R. & Abers, G., 1994. MT5 program, IASPEI Software Library, 4.

S U P P O RT I N G I N F O R M AT I O N
Supplementary data are available at GJ I online.
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 article.

A P P E N D I X A : D I S T R I B U T E D S L I P M O D E L R E S O L U T I O N
The resolution matrix for the distributed slip model, R, is dependent on the smoothed Green's function, where G s is [Gκ 2 ∇ 2 ] T (Jonsson et al. 2002;Funning et al. 2005). The resolution matrix is given by The diagonal values of R represent the model resolution. We define the horizontal and vertical resolution length scales for each patch as the distance in each direction over which the values are greater than 1/e of the maximum value of the resolution matrix (Funning et al. 2005;Biggs et al. 2006). For perfectly resolved models R will be an identity matrix. The distributed slip model for the Zinave earthquake has 7.5 and 5 km horizontal and vertical resolutions in the peak slip region, respectively (Fig. A4). The inversion for the slip distribution of the Machaze earthquake using ENVISAT data (model 1) has a more variable resolution, but is ∼14 km in the horizontal and vertical in the peak slip region (Fig. A8).       Table A3. Fault parameters from the geodetic inversions, with 95 per cent confidence, and USGS focal mechanism for comparison. Fault locations (X and Y, UTM) in the geodetic inversions are for the middle of the bottom of the fault. Root-mean-square (rms) misfit is the joint rms if more than one data set is used in the inversion.