Uplift and fault slip during the 2016 Kaikōura Earthquake and Late Quaternary, Kaikōura Peninsula, New Zealand

ABSTRACT The Kaikōura Earthquake uplifted Kaikōura Peninsula by ≤∼1 m. Uplift in 2016 mainly resulted from slip on an offshore thrust fault (OSTF), modelled to splay from the plate-interface, and was further influenced by slip on two newly identified faults (Armers Beach Fault, ABF; Te Taumanu Fault, TTF) mapped onshore from differential lidar (D-lidar). Forward dislocation modelling indicates that 2016 peninsula uplift can be reproduced by mean slip of ∼2.3 m on the OSTF and 0.25–0.5 m on the ABF and TTF. The variable co-seismic uplift recorded during the 2016 earthquake differs from the near-uniform (1.2 ± 0.2°) northwest tilting of MIS5c (96 ± 5 ka) and MIS5e (123 ± 5 ka) marine terraces; these ages are constrained by Optically Stimulated Luminescence (OSL) dating and correlation to sea-level curves. Tilting of Late Quaternary marine terraces can be primarily reproduced by slip rates of ∼0.8–2.7 mm/yr on the OSTF and 0.3–0.6 mm/yr on the ABF. Slip on the TTF is not required to produce tilting of the marine terraces, suggesting that it may have ruptured less frequently than the OSTF and ABF in the Late Quaternary. The OSTF links 2016 ruptures north and south of Kaikōura, with the earthquake rupturing an interconnected network of faults.


Introduction
The 2016 Kaikōura Earthquake produced uplift along at least 100 km of the northeastern coastline of New Zealand's South Island (Figure 1). Coastal uplift reached up to ∼6 m and can be partly accounted for by slip on the Hundalee, Papatea and Kekerengu-Needles faults, which ruptured the ground surface and seabed in 2016 Kearse et al. 2018;Langridge et al. 2018;Williams et al. 2018;Mouslopoulou et al. 2019;Zinke et al. 2019;Howell et al. 2020 this issue). Elastic dislocation and tsunami modelling, together with aftershock relocations, collectively indicate that some 2016 uplift also requires slip on a shallow westward-dipping (<40°) thrust fault(s) beneath the coastline and immediately offshore (Bai et al. 2017;Clark et al. 2017;Hollingsworth et al. 2017;Power et al. 2017;Mouslopoulou et al. 2019;Chamberlain et al. 2021) (Figure 1). Slip on the thrust caused coseismic uplift of the Kaikōura Peninsula by up to ∼1 m and is a key focus of this study.
Kaikōura Peninsula was uplifted and strongly shaken during the Kaikōura Earthquake Kaiser et al. 2017;Mouslopoulou et al. 2019;Reid et al. 2020). The peninsula trends across the northeast structural grain of the margin and is located at least ∼15 km from onshore fault rupture of the ground surface in 2016 with >1 m displacement (e.g. Hundalee, Papatea and Kekerengu faults). The peninsula offers a rare opportunity to study uplift in 2016 across multiple faults that did not rupture the ground surface and were either not represented in surfacetrace maps of the 2016 rupture or located with significant uncertainty (e.g. Litchfield et al. 2018;Zinke et al. 2019;Howell et al. 2020). Uplift of the peninsula in 2016 is most commonly thought to be due to slip on a northwest dipping thrust that approaches the seabed along the shelf break, ∼10 km southeast of the peninsula Mouslopoulou et al. 2019;Duffy 2020;Howell and Clark, in press) (Figure 1 cross section). Mouslopoulou et al. (2019) indicate that this thrust fault forms part of the accretionary prism complex in the southern Hikurangi margin and splays from the plate-interface transferring subduction-related slip onto upper-plate faults (Figure 1). This transfer of displacement is part of the transition from the subduction system to the strike-slip faults of the Marlborough Fault System (Figure 1) (Holt and Haines 1995;Barnes and Audru 1999;Wallace et al. 2012;Mouslopoulou et al. 2019). Kaikōura Peninsula preserves a flight of Late Pleistocene to Holocene marine terraces that record uplift rates of up to ∼1.5 mm/yr and support the notion that earthquake-related slip on faults beneath the peninsula has been occurring for >100 kyr (Bull 1984;McFadgen 1987;Ota et al. 1996;Duffy 2020;Howell and Clark, in press) (Figures 2 and  3). The combination of 2016 uplift and well-preserved marine terraces provide an opportunity to compare deformation over a range of timescales and to answer the question 'was vertical displacement during the Kaikōura Earthquake characteristic of the long-term uplift on Kaikōura Peninsula?'.
In this paper, we use a lidar-derived digital elevation model (DEM), differencing of pre-and post-earthquake lidar (hereafter referred to as Dlidar), Optically Stimulated Luminescence (OSL) dating of Late Pleistocene marine terraces and elastic dislocation forward modelling to estimate the geometry, slip and slip rates for faults that ruptured in the 2016 earthquake and during the Late Quaternary. We use uplift estimates derived from D-lidar to locate the faults that ruptured in 2016 and assess the extent to which these faults experienced repeated slip in the Late Quaternary. D-lidar reveals the presence of two previously unmapped faults with a component of reverse displacement and mainly northeast strikes, parallel to the subduction margin. Here, we refer to these new faults as the Armers Beach Fault (ABF) and Te Taumanu Fault (TTF). Additionally, we use lidar to map and reconcile previously mapped marine terraces, chart the spatial distribution of paleo-sea cliffs on  Beavan et al. (2002). B, Geology from Rattenbury et al. (2006): Basement (mostly Torlesse Terrane)blue; Late Cretaceous-Cenozoic stratadark yellow; Quaternary depositslight yellow. Active faults (black lines from Langridge et al. 2016) and November 14, 2016 Mw 7.8 Kaikōura Earthquake fault ruptures (red lines with key fault names shown, modified from Mouslopoulou et al. 2019). The epicentral location (red star) and the moment tensor of the 2016 Kaikōura Earthquake are from Chamberlain et al. (2021) and GeoNet, respectively. Onshore topographic basemap is from LINZ and bathymetry is from NIWA. Abbreviations are as follows: MFS, Marlborough Fault System; NCD, North Canterbury Domain; PKF, Point Kean Fault; OSTF, Offshore Splay Thrust Fault; TTF, Te Taumanu Fault; ABF, Armers Beach Fault; KC, Kaikοura Canyon. The localities of Figures 2 (F2) and 3 (F3), that focus on the Kaikōura Peninsula, are indicated by the rectangles labelled F2 and F3, respectively. C, Regional NW-SE cross section showing upper plate faults in the region surrounding the study area and their geometric relations to the plate interface (modified from Duffy 2020). The location of the plate interface is from Williams et al. (2013). The section shows no vertical exaggeration. the peninsula and mainland, and measure displacement across a newly recognised fault that crosses the peninsula. Recognition of four stranded sea-cliff and marine-terrace couplets, five new OSL dates from marine-terrace cover beds and marine-terrace correlation to the high resolution stacked marine isotope sea-level curve of Lisiecki and Raymo (2005), collectively provide increased confidence in our revised correlations. The available data support the view of Duffy (2020) that the thrust fault underlying the Kaikōura Peninsula repeatedly ruptured during earthquakes in the Late Quaternary. However, the 2016 profile of uplift along the peninsula does not appear to match the tilting towards the northwest of the Late Pleistocene marine terraces. We attribute these differences to rupture of a newly identified fault northwest of the peninsula (Figures 2 & 4) that slipped in 2016 and appears to have ruptured infrequently in the Late Quaternary. Our observations indicate that fault rupture during the Kaikōura Earthquake was more complex than previously thought and that, at least locally on Kaikōura Peninsula, this complexity does not appear to have been completely replicated in Late Quaternary earthquakes.

Kaikōura Earthquake and regional tectonic setting
The 2016 M w 7.8 Kaikōura Earthquake nucleated on The Humps Fault, near the township of Waiau, at a depth of Figure 2. A, Hillshade of the Kaikōura Peninsula area illustrating its geomorphology and geographic names referred to in the paper. The location of the Armers Bay Fault (annotated ABF) is indicated by the red dashed line. Locations, with annotated lab numbers for OSL ages are shown as red filled circles. State Highway 1 (SH1) mapped as a white line. B, Hillshade model from A draped with coloured, labelled benches and sea cliffs. OSL sample locations are labelled with ages. Note that marine terrace MtB (coloured orange) lies above two sea cliffs on northwest and southeast sides of the ABF, and is vertically displaced (down to the SE) by the fault.
13 ± 2 km and propagated to the northeast for about 180 km ( Figure 1) (Cesca et al. 2017;Kaiser et al. 2017;Ando and Kaneko 2018;Nicol et al. 2018;Lanza et al. 2019;Chamberlain et al. 2021). The earthquake focal mechanism is oblique, resulting from a combination of right lateral and reverse slip (GeoNet;www.geonet.org.nz/ earthquake/technical/2016p858000, last accessed October 30, 2021: USGS; www.earthquake.usgs.gov/ earthquakes/eventpage/us1000778i, last accessed October 30, 2021; also see Kaiser et al. 2017) (Figure 1 inset), while the faults that ruptured the ground surface during the earthquake are also dominated by right-lateral and reverse slip (Kearse et al. 2018;Litchfield et al. 2018;Nicol et al. 2018;Williams et al. 2018; Zinke    Howell et al. 2020). The earthquake ruptured ≥17 faults with >1 km length at the ground surface, producing one of the most complex historical ruptures globally Litchfield et al. 2018;Nicol et al. 2018;Mouslopoulou et al. 2019;Zinke et al. 2019;Hamling 2020;Howell et al. 2020) (Figure 1). Although the earthquake initiated in North Canterbury where fault displacement rates are low (<2.5 mm/yr), the maximum fault slip of ∼12 m at the ground surface and the greatest energy release during the earthquake occurred north of Kaikōura, in the region of the Kekerengu and Papatea faults ( Figure 1) (Cesca et al. 2017;Kaiser et al. 2017;Kearse et al. 2018;Howell et al. 2020). The plate-interface underlies the 2016 surface ruptures at depths of ∼20-35 km (cf Williams et al. 2013) and its role in the earthquake continues to be debated. It is possible that the plate-interface did not accommodate significant (e.g. >∼2 m) slip during the earthquake (e.g. Hamling et al. 2017), although many consider that the interface accommodated some slip (e.g. ∼2-10 m) during (Bai et al. 2017;Duputel and Rivera 2017;Furlong and Herman 2017;Hollingsworth et al. 2017;Wang 2017;Mouslopoulou et al. 2019) and after the earthquake (Wallace et al. 2018;Mouslopoulou et al. 2019). Despite the divergence of views, there is general agreement that the majority of seismic moment released during the earthquake was generated by slip on upper-plate faults and that slip occurred on at least one blind thrust beneath coastal Marlborough and offshore. Clark et al. (2017) and Mouslopoulou et al. (2019) modelled slip on a thrust fault that dips ∼35°northwest and nominally intersects the sea-bed southeast of Kaikōura Peninsula. This fault is referred to as the 'Point Kean Fault' in some studies (e.g. Clark et al. 2017;Chamberlain et al. 2021), which is marked by a clear seabed scarp that is interpreted to have ruptured in 2016 (Litchfield et al. 2018). The mapped offshore Point Kean fault is however downthrown to the northwest (Litchfield et al. 2018) and cannot be responsible for uplift of Kaikōura Peninsula. Therefore, we believe that the Point Kean fault should not be used to define the location or name of the thrust(s) underlying Kaikōura Peninsula. Instead, Mouslopoulou et al. (2019) used seismic-reflection profiles, tsunami modelling, aftershock relocations and sea-bed bathymetry to approximately locate and model slip on a northwest-dipping thrust that intersects the plateinterface at c. 30 km depth ( Figure 1). Here, we follow Mouslopoulou et al. (2019) and refer to the main thrust that resulted in uplift of the Kaikōura Peninsula as the Offshore Splay Thrust Fault (OSTF) and in our model it projects to the seabed ∼3 km east of the Point Kean fault.
Kaikōura Peninsula is located south of the Hope Fault and the Marlborough Fault System (MFS) ( Figure 1). The Cretaceous to Cenozoic bedrock exposed on Kaikōura Peninsula has been folded into a series of kilometre-scale anticlines and synclines that trend northeast, parallel to the faults in the accretionary prism complex offshore (Campbell et al. 2005;Rattenbury et al. 2006). On the peninsula, folding of the middle Miocene Waima Siltstone Formation produces bed dips of up to ∼80°and records NW-SE contraction consistent with the direction of plate convergence (Campbell et al. 2005;Rattenbury et al. 2006). The folded strata have been abraded by wave action to form a flight of marine benches (here referred to as marine terraces) of Late Pleistocene and Holocene age (McFadgen 1987;Ota et al. 1996;Campbell et al. 2005;Duffy 2020; in press) ( Figure 2). These terraces were approximately horizontal when formed and are presently up to ∼105 m above mean sea-level, with each bedrock bench locally overlain by up to ∼6 m of beach deposits and loess cover beds (McFadgen 1987;Ota et al. 1996;Campbell et al. 2005). The marine terraces were most likely cut during sealevel high stands and subsequently raised above sea-level (and the influence of wave action) by tectonic uplift. The highest two Late Pleistocene terraces, which are the focus of this study, are tilted towards the northwest at 1.2 ± 0.2° (Suggate 1965;Ota et al. 1996;Duffy 2020) (Figures 3 and 5B). The tilting and uplift recorded by the Late Pleistocene marine terraces on Kaikōura Peninsula are proposed to record slip on a reverse fault southeast of the peninsula (Ota et al. 1996;Barrell 2015;Duffy 2020). In this paper we further explore the relationships between cumulative slip on the offshore thrust fault (i.e. OSTF) and tilting of these marine terraces.

Data and methods
Magnitude of uplift, uplift rates and tilting on Kaikōura Peninsula for the 2016 earthquake and the Late Quaternary have been determined using two lidar DEMs, one collected pre-earthquake and a second post-earthquake. The two lidar models have been used to derive differential lidar (D-lidar) which represents the net displacements of the ground surface in the 2016 earthquake. Dating of marine terrace cover beds by Optically Stimulated Luminescence (OSL) was used to support the assigned ages of these terraces, while fault slip and slip rates were derived using elastic dislocation modelling.

Lidar DEM
A lidar-generated DEM collected following the earthquake in December 2016 has been used to remap Late Pleistocene marine terraces originally mapped by Ota et al. (1996) and Duffy (2020). These lidar data were processed to produce a metre-scale resolution hillshade model (AAM NZ 2016; see description in following paragraph) ( Figure 2A and Figure 3). The lidar DEM has estimated absolute horizontal and vertical accuracies of ±0.5 m and ±0.1 m, respectively, however, local relative uncertainties are likely smaller than 0.1 m. As the lidar was collected after the earthquake, the altitudes of marine terraces include up to 1 m of earthquake uplift, which is only a small proportion of the marine terrace altitudes of ∼40-105 m amsl and does not significantly impact the geometries of these terraces. These terraces formed as sub-horizontal benches dipping at low angles (<1°) away from the peninsula axis and are separated vertically by paleo sea-cliffs ( Figure 2B). The terrace and paleo sea-cliffs are locally dissected by small gullies and are discontinuous across the ABF (Figure 2A and Figure  4). Terrace-sea cliff couplets have been used to correlate the marine terraces across the ABF. Underpinning their use as a correlation tool is the assumption that the same couplets formed on each side of the fault. Some of our terrace correlations differ from previous interpretations and enable fault vertical displacements post terrace formation to be estimated on one terrace (MIS5c, named MtB in this paper).

Differential lidar (D-lidar)
D-lidar is our primary tool for measuring uplift across, and immediately inland of, Kaikōura Peninsula during the 2016 earthquake ( Figure 4A). The D-lidar model was produced by differencing the altitude of 1 × 1 m pixels using pre-earthquake lidar (captured for Environment Canterbury Regional Council by Aerial Surveys in 2012 with a vertical accuracy of <+/−0.1 m at 95% confidence level), and a post-earthquake lidar (AAM New Zealand, Nov 20 2016 to January 2017, for Environment Canterbury Regional Council, that has a vertical accuracy of <+/−0.10 m at 95% confidence level); the 2012 lidar has been reprocessed to the same geodetic datum as the 2016 lidar (NZVD2016). Due to the impact of gravity induced slope failures, growth of vegetation and horizontal tectonic displacement of about 2 m in the study area (Zinke et al. 2019;Howell et al. 2020), the best Dlidar uplift measurements are from unvegetated flat surfaces. In particular, the most reliable data are from the centre lines of roads and grassed flat pastures. Where topography is sub-horizontal, uplift estimates from the D-lidar typically have 2σ uncertainties of ±0.06-0.18 m at individual grid nodes, however, Figure 5. Slip models derived from forward modelling for the A, lidar vertical displacements due to the 2016 Kaikōura Earthquake and B, long-term uplift measured on the MtB marine bench cut around the peninsula. Normalised slip models for the OSTF are shown in both cases along depth (grey lines). The best-fit slip model with the minimum residual between data and predicted data is shown in red. The medium panels show, in each case, the comparison between the observed (black line and dots) and the predicted (red line and dots) displacements. While the lower panels, show the slip distribution at depth of the best-fit slip models along each fault (left panel) and the predicted vertical displacement pattern at the ground surface (right panel). Profile A indicates the location of data presented in Figure 7.
averaging vertical-displacement values for 20 or more adjacent nodes typically reduces the uncertainty by at least a factor of four. For the purposes of this study we assume mean 2σ uncertainty on the averaged 2016 uplift estimates of ±0.05 m. The uncertainties are sufficiently small to permit identification of 'steps' in the D-lidar model of ∼0.2-0.5 m, which we attribute to sub-surface slip on the ABF and TTF faults during the 2016 earthquake ( Figure 1 map and cross section & Figure 4). These faults have been mapped onshore using uplift gradients identified on D-lidar profiles as a result of the 2016 earthquake (e.g. Figure 4B-F). The uplift measured from D-lidar is consistent with estimates of 2016 coastal uplift from intertidal marine biota and the Kaikōura Tide Gauge (see Figures 1 and 4A for location) Mouslopoulou et al. 2019;Reid et al. 2020), supporting the assertion that uplift recorded by the D-lidar occurred during the earthquake.

Optically Stimulated Luminescence (OSL) dating
We have primarily used OSL dating and correlation to the sea-level curve of Lisiecki and Raymo (2005) to estimate the ages of Late Pleistocene marine terraces. When the age of the marine terrace is known (or can be estimated), sea-level curves provide a means of estimating the altitude (relative to present sealevel) at which the marine bench formed, the absolute uplift and the time-averaged uplift rates since their abandonment. To estimate uplift, we have used the detailed stacked sea-level curve of Lisiecki and Raymo (2005) because their 1 kyr oxygen isotope data resolution is precise for the last ∼130 ka and, when used in conjunction with Antonioli et al. (2004) and Rabineau et al. (2006), provides uncertainties for temporal changes in absolute sea-level. The broad shape of sea-level curves is stable between different publications and curve selection does not impact our first-order estimates of uplift and uplift rates.
The Late Pleistocene terraces have estimated ages between ∼56 ka and ∼130 ka primarily inferred from correlation to sea-level high stands (Suggate 1965;Bull 1984;Ota et al. 1996;Campbell et al. 2005; Duffy 2020) ( Figure 5). To provide additional age constraints for the marine terraces we have OSL-dated five loess samples from cover beds resting on the bedrock benches beneath the terrace surfaces (Table 1 and Figure 5). The dated samples were all collected from within 30 cm of the bedrock bench, and represent minimum ages for the cessation of the wave abrasion processes that formed each bench. In eastern New Zealand, because loess is most often deposited during cold-climate time intervals (Milne and Smalley 1979;Berger et al. 2001), at best the loess samples will provide an age for the cool stadial period that immediately post-dates the warm-climate period during which the underlying marine terrace was cut. OSL dates are from bulk samples comprising multiple grains of feldspar using the Victoria University of Wellington laboratory (see Table 1 for further details). OSL ages are quoted with 1σ uncertainties and generally increase in age with greater height above sea-level. OSL sample WLL1252 (33.1 ± 1.9 ka, Table 1) is a notable exception, as it was sampled from the highest terrace and produced the youngest OSL age. This age highlights the potential for OSL ages to be reset (i.e. due to Late Pleistocene exposure to sunlight) and/or for loess accumulation rates to vary in space and time. The OSL age for sample WLL1252 is not considered further in this paper.

Fault-slip modelling
The magnitude of fault slip and slip rates at depth beneath Kaikōura Peninsula in 2016 and during the Late Quaternary were estimated using forward elastic-dislocation modelling (cf Okada 1985) applied to the OSTF, ABF and TTF. Elastic-dislocation modelling has been widely used to reproduce permanent uplift of Late Quaternary marine terraces (e.g. Anderson and Menking 1994; Jara-Munoz et al. 2019). Our models assume that all uplift resulted from co-seismic fault slip and do not incorporate viscoelastic or plastic deformation manifest as fault creep. Viscoelastic or plastic deformation could only account for up to 30% of the co-seismic moment release (Melbourne et al. 2002), and its exclusion does not impact our first-order conclusions.
For each forward model, uplift (vertical displacement) at the ground surface was predicted using a range of fault slip values. The predicted uplift was then compared to the observed uplift and models that produced results within the uncertainty bounds of the data were selected. We, therefore, produced a set of slip-models with displacement patterns within the uncertainty bounds of the data that are considered permissible, with the best-fit slip-model having the lowest average residuals to the data. To estimate surface deformation and slip at depth, we calculated Green functions on 0.25 × 0.25 km rectangular fault patches along the plate interface, OSTF, ABF and TTF faults, using the SDM2011 software (available at https://www.gfz-potsdam.de/en/section/physics-of-ea rthquakes-and-volcanoes/infrastructure/tool-develop ment-lab/) (Wang et al., 2012).
Preliminary modelling suggests that slip on the plate-interface (Williams et al. 2013, interface model) and on upper-plate faults located >5 km from the peninsula, have negligible impact on uplift of the peninsula ( Figure S1 supplementary material). For this reason, we have only included in our modelling the sections of the OSTF, ABF and TTF which are located within 5 km from the peninsula. Using the preferred geometries and slip sense characteristics of the modelled faults, we have evaluated the optimal combination of slip on the OSTF, ABF and TTF faults to produce the observed uplift profile along the axis of the peninsula for the 2016 earthquake ( Figure 5A) and the Late Pleistocene marine terraces ( Figure 5B).
For the results presented in Figure 5 all faults are assumed to be rectangular with planar dips and pure dip-slip motion (rake angle: 90°) (Figures 5 and S1). Geometries and depths for all modelled faults are presented in Table 2. The location and geometry of the OSTF is adopted from Mouslopoulou et al. (2019) ( Table 2). The OSTF has a modelled dip of 35°and splays from the plate-interface at ∼30 km depth; the proposed intersection of the OSTF and the plate interface is ∼30 km northwest of the study area ( Figure  1C). The approximate location, size and dip of the OSTF is constrained using a combination of aftershock relocations, seismic-reflection profiles, sea-bed bathymetry and tsunami modelling (Mouslopoulou et al. 2019 and references therein). Given the available data, it is possible that the OSTF is listric (e.g. Duffy 2020), however, our sensitivity testing suggests that the modelled fault slip typically varies by <30% for a range of planar and listric fault geometries. The locations and strike of the ABF and TTF faults were determined from the D-lidar data and assumed to have slipped from 500 m below the ground surface to depths of 7 and 4 km, respectively. Dips of 60°a nd 55°were assumed for the ABF and TTF faults, respectively. A ± 10°change in the dip angle of these onshore faults results in small (<20%) differences in the modelled slip that does not impact the first-order results of this study.
Uniform slip models for the OSTF were not able to reproduce the observed uplift gradients along the peninsula ( Figure S2). Thus, in order to introduce some slip variability along the OSTF with depth, we have divided the fault plane into four up-dip sections ( Figure 5; 0-4.5 km, 4.5-8.5 km, 8.5-15.5 km and 15.5-30 km; Figure S1 supplementary material), each assigned with uniform slip. The smaller onshore faults (ABF and TTF) were assumed to slip uniformly. Forward calculations were performed in order to select potential slip models that reproduce the surface peninsula uplift within the uncertainty intervals of ±10 m and ±0.1 m for the uplifted marine terraces and the 2016 coseismic uplift, respectively. For the long-term uplift, total slip of 25-800 m and 0-100 m was tested for the OSTF and ABF, respectively, using 25 m steps in slip between each candidate model. TTF was excluded from modelling of the tilt of Late Pleistocene terraces as its sense of slip was not compatible with their observed tilts. For modelling of the 2016 coseismic displacements, we tested the slip values of 0.25-6 m every 0.25 m for the OSTF and 0.25-1 m every 0.25 m for the ABF and TTF. For both sets of models (longterm and co-seismic), slip values were tested equally along all four pre-defined depths on the OSTF. A set of slip-models (with a range of mean-slip for each fault section) was obtained replicating the observed uplift patterns, within the range of uncertainties. 3.93 ± 0.18 0 4.02 ± 0.22 38.5 ± 1.7 42.3 ± 2.4 ‡ All samples analyzed at the Victoria University of Wellington OSL laboratory with measurements taken using blue luminescence of fine-grained feldspar produced during infrared stimulation. †All ages for multi-grain samples using the single Aliquot Regeneration method, reported with 1 sigma uncertainties. *The dose rate of WLL1253 is the average dose rate of the remaining four samples in the

Kaikōura Earthquake uplift
In the region of the Kaikōura Peninsula, 2016 coastal uplift ranges from zero near Peketā ∼8 km southwest of the peninsula to about a metre on the peninsula itself ( Figure 4A). Inspection of the uplift colour ramp in Figure 4A shows an uplift low, west and northwest of the peninsula. A series of profiles (see Figure 4D-F) show that up to 0.5 m decrease in uplift north and west of the peninsula, primarily occurs across a 0.5-1 km wide zone that can be mapped across alluvial plains and passes offshore north and south of the peninsula ( Figure 4A). The mapped north-westward down step in uplift, defines a monocline which ranges in trend from north-northeast to east-northeast. Here, we propose that monoclinal warping of the alluvial terraces formed due to slip at depth on a southeast dipping reverse fault ( Figure 5A cross section), which we refer to as the TTF. The fault is assumed to project to the ground surface along the base of the monoclinal limb. Field and lidar mapping immediately following the 2016 earthquake did not reveal fracturing or rupture of the ground surface along the TTF and this structure is mapped for the first time here. We infer that slip on the TTF during the 2016 earthquake was blind reverse slip, although it is also possible that shear-strains are entirely accommodated by monoclinal folding with little or no fault-slip at depth. Southeast of the TTF and monocline structure, 2016 uplift varies from about 0.7-0.9 m for ∼6 km, with the highest values typically in the immediate hanging wall (northwest wall) of the ABF ( Figure  4B profile 4 & 5A left graph labelled predicted data). Southeast of the ABF, uplift decreases abruptly to ∼0.6 m with this value approximately remaining constant on the southeast side of the fault to the south-eastern end of the peninsula. Inspection of a series of profiles trending parallel ( Figure 4B) and perpendicular ( Figure 4C) to the axis of the peninsula show up to 0.3 m vertical displacement on the ABF distributed across several hundred metres normal to fault strike. This displacement can be observed across the width of the peninsula and coincides with a ∼45 m change in the height of topography which is highest in the hanging wall (northwest side) of the fault. The trend of the ABF, as mapped from northwest-southeast D-lidar uplift profiles, is northeast and subparallel to the strike of bedrock bedding (see Campbell et al. 2005& Rattenbury et al. 2006. No surface rupture (e.g. cracking, fissuring and fault displacement) was mapped along the ABF in 2016 and we propose that the differential displacement recorded by the D-lidar reflects blind slip on the ABF at depth in 2016.

Correlation, age and uplift rates of marine terraces
Late Pleistocene marine terrace heights and spatial distributions of Ota et al. (1996) and Duffy (2020) are similar, while the mean ages they assigned to these terraces differ by up to 19 ka. Our interpretations have common elements to both previously published studies of the marine-terraces, however, through mapping of stranded marine-cliffs and the recognition of the ABF, our terrace correlations between the northwest and southeast ends of the peninsula may be more reliable than previous publications. The revised correlations are further supported by existing dating (Ota et al. 1996), cover-bed stratigraphy (Ota et al. 1996) and new OSL dating (Table 1 and Figure 6).
The location of the ABF in Figures 2 and 4 is based on our interpretation of the D-lidar data. The fault coincides with two deeply incised gullies and a northeast-trending step in the topography up to 45 m in height (Figures 2 and 3). The northwest side of the step in topography is higher than the southeast and the up-step direction is in the same sense as the uplift and modelled 2016 displacement on the ABF. For the interpretations of Ota et al. (1996) and Duffy (2020) the highest three marine terraces (I/1, II/2 and III/3) were not correlated across the topographic step. Here, we propose that the poor correlation of the older marine terraces across the topographic step is due in part to slip on the ABF (either at the surface and/or at depth, as occurred in 2016), which raised the northwest side of the fault relative to the southeast side. We propose that marine terraces II/2 and III/3, previously interpreted to be separate terraces cut during different MISs, are the same terrace that has been vertically displaced by 23 ± 5 m on the ABF at the ground surface ( Figure 3). In support of this argument, we have mapped marine-terrace and sea-cliff couplets immediately adjacent to the fault. Counting back from the Holocene marine-terrace and the associated Holocene sea-cliff on each side of the fault, we find that previously mapped marine terraces II/2 and III/3 occupy the same stratigraphic position in the flight of terraces, supporting our new proposed correlation. In addition, re-examination of stratigraphic auger holes from Ota et al. (1996) on their terraces II and III show that the cover bed thickness is 3-3.5 m on both terraces. Finally, the cover beds on both terraces are comparable and comprise ∼0.5 m of pebble to fine cobble gravel (clasts 1-4 cm diameter) at the base of the sequence and resting on the bedrock strath, which supports the view that these two terrace surfaces could have formed synchronously.
Using the new correlation across the ABF, we map three main marine terraces (from lowest to highest, MtA, MtB, MtC; and an additional subsidiary terrace, MtAa) on the peninsula itself ( Figure 2B and Figure 3). MtC correlates with terrace I and 1 of Ota et al. (1996) and Duffy (2020), respectively, and is the highest terrace of the peninsula sequence. An amino acid date of 110 ± 20 kyr from Ota et al. (1996) remains the best available direct constraint for the MtC terrace age; our OSL date from this terrace is the youngest of those sampled and is interpreted to date loess that accumulated (or was exposed to sunlight) long after the terrace was abandoned by the sea. Correlating MtC with MIS5e and using the sea-level curve of Lisiecki and Raymo (2005) we estimate the age of this terrace to be 123 ± 5 ka which, within the uncertainties, is consistent with the 110 ± 20 ka age from Ota et al. (1996) and the 123 ± 5 ka estimate of Duffy (2020). We have two OSL dates from the second-highest terrace, here referred to as MtB (Table 1). These dates produced an average age of 87.7 ± 4.8 ka for the two MtB terrace samples (see Table 1, Figure 6), which is here proposed to have formed during MIS5c (96 ± 5 ka, Table 3). OSL samples from cover beds on the MtA and MtAa terraces returned ages of 38.5 ± 1.7 ka and 42.3 ± 2.4 ka, respectively. These dates are minimum values and the ages of the MtA and MtAa terraces are relatively poorly constrained. Here, we infer that the MtA terrace correlates with MIS5a with an age of 83 ± 5 ka, while MtAa is typically only ∼5 m lower in altitude than the MtA and may also be from MIS5a. It is also possible that terrace MtAa correlates with MIS3 with an age of 55 ± 5 ka. Given the uncertainty in their ages, the MtA and MtAa marine terraces have not been used in this paper to determine the rates of uplift.
Using the altitudes and correlation to the sea-level curve (Lisiecki and Raymo 2005), we have estimated the ages and maximum uplift rates for the MtC and MtB marine terraces. Maximum uplift rates for both terraces were measured in the hangingwall of the ABF where, based on our revised age estimates, values of about 1 mm/yr are recorded (Table 3). These rates are comparable to the 1-1.3 mm/yr and 1-1.5 mm/yr estimated by Ota et al. (1996) and Duffy (2020) for Late Pleistocene terraces, respectively. Our Late Quaternary uplift rates are also comparable to the ∼0.9 mm/yr calculated here for the Holocene (using 6 m amsl for the altitude of the base of the Holocene sea cliff and an age of 6.5 kyr for its abandonment), and the ∼1.5 mm/yr uplift rate calculated from Holocene marine terraces by McFadgen, (1987). In addition, we have calculated vertical separation rates (i.e. either fault displacement at the ground surface or folding at the ground surface with fault slip at depth) of 0.18-0.31 mm/yr for the ABF measured across the MtB terrace (23 ± 5 m over 96 ± 5 ka, Table 3).  Lisiecki and Raymo (2005) with absolute sea level relative levels averaged from sources summarised in Antonioli et al. (2004) and Rabineau et al. (2006). Grey polygon on the sea-level curve shows the 2σ uncertainties on the altitude of sea level relative to the present. Details of OSL ages (red circles) are given in Table 1. Altitudes of sea-level formation and their uncertainties together with terrace ages are tabulated in Table 3. See Figure 2 for map and legend of the marine terraces.

Fault slip and slip rates
Uplift on Kaikōura Peninsula in the 2016 earthquake and during the Late Quaternary are widely considered to reflect slip beneath the peninsula on a north-westward dipping reverse or thrust fault (e.g. Ota et al. 1996;Barrell 2015;Mouslopoulou et al. 2019;Duffy 2020). Our fault-slip models constructed to account for the observed uplift include a combination of slip on the OSTF, ABF and TTF ( Figure  6). They provide a first-order indication of the slip and slip rate on these faults beneath Kaikōura Peninsula. Forward modelling of the D-lidar data resulted in 279 slip-models within the uncertainty interval of ±0.1 m ( Figure 5A). The OSTF best-fit slip model for 2016 uplift of Kaikōura Peninsula ranges from a maximum slip of 6 m in the upper ∼5 km of the crust to 2-2.5 m below Kaikōura Peninsula, at ∼5-15 km depths km, however, many of the models also have maximum slip at depths of 5-15 km. Mean-slip for the OSTF is 2.3 m (with a range of mean slip for all models of 0.9-2.5 m), while for the ABF and TTF mean slip is 0.5 and 0.25 m (mean-slip range: 0.25-0.5 m for the ABF and TTF), respectively. Because the modelled OSTF fault extends to the subduction interface (at 30 km depth) and has a low dip (35°), it is the primary contributor to the empirically measured long-wavelength co-seismic 2016 uplift of about 0.4-0.7 m, which extends from ∼8 km south to 10 km north of Kaikōura Peninsula ( Figure 5A). Results are also consistent (within 10 cm) with co-seismic horizontal displacements recorded 5 days after the earthquake by a campaign GPS station located on the peninsula (ID: 1163; Hamling et al. 2017).
Late Quaternary uplift of the peninsula decreases towards the northwest, resulting in tilting of the marine terraces towards the northwest. Our interpretation of the northwest tilt of the marine terraces is that it was potentially produced by slip on the OSTF and ABF faults. A total of 170 slip models were obtained from forward modelling of the Late Pleistocene marine terraces ( Figure 5B). The best-fit slip-model (residual: 1.8 m) indicates that long-term uplift of the peninsula can be accommodated by mean slip of ∼190 m over the OSTF slip surface (mean-slip range: 80-230 m), with the ABF contributing about 50 m of slip at depths of 0.5 to ∼8 km (mean-slip range: 25-50 m). The peninsula tilting is primarily (∼88%) accounted for by slip on the OSTF, with the remainder (∼12%) due to slip on the ABF. The range of mean-slip of 80-230 m on the OSTF modelled here, is lower than the 364-403 m (mean 386 m) proposed in the geometric listric-fault model of Duffy (2020) for the MtC terrace. We have further tested these by developing our own slip models with listric fault geometries for the OSTF. In these listric models dip angle was decreased linearly from 60°to 30°from the surface to 25 km depth. Our listric fault slip models produced a mean-slip range of ∼150-270 m (best mean-slip model: 250 m) ( Figure S3), and are broadly consistent with our planar-fault models and the Duffy (2020) listric model. Moreover, sensitivity analysis of uniform slip-models for different geometries of the OSTF ( Figures S1 and S2), indicate that complexity in slip variability or geometry at shallow depths (<10 km), is required to reproduce the tilt, in agreement with Duffy (2020).
Given the uplift of the peninsula during the 2016 Kaikōura Earthquake, we infer that the observed Late Quaternary uplift primarily resulted from repeated earthquakes. Our models have implications for slip rates and earthquake recurrence intervals on the OSTF and ABF. If mean slip on the OSTF for all models ranges between 80 and 230 m since formation of the MtB terrace at 86-106 ka (2σ uncertainties), then the fault-slip rate would be ∼0.8-2.7 mm/yr. This slip rate is within the range of 0.8-5.4 mm/yr and the 1-3 mm/yr proposed by Duffy (2020) and Barrell (2015), respectively. If the 2016 modelled mean slip of 2.3 m (+0.2 −1.4 m) and our calculated slip rate of 0.8-2.7 mm/yr are both representative of the Late Quaternary earthquake activity on the OSTF, then the earthquake recurrence interval would be ∼340-3300 years. The average recurrence interval for the best-fit model with slip rates of 1.8-2.2 mm/yr is ∼900-2000 years. Similarly, if mean slip on the ABF for all models ranges between 30-50 m at depth since formation of the MtB terrace at 86-106 ka (2σ uncertainties), then the fault-slip rate would be ∼0.3-0.6 mm/yr and the recurrence interval ∼500-1300 years. Average recurrence intervals for the OSTF and ABF are typically longer than values of ∼400-800 years for the Late Holocene (McFadgen 1987;Duffy 2020;Howell and Clark, in press). This disparity in recurrence intervals could arise if the OSTF and ABF typically do not rupture together and/or may reflect elevated Late Holocene rates of earthquakes that co-rupture the OSTF and ABF. 6.5 ± 0.5 6 ± 0.5 0 6 ± 0.5 -0.9 ± 0.1 -

Discussion
The available data allows us to compare Late Quaternary and 2016 uplift to determine whether deformation during the Kaikōura Earthquake was representative of the long-term average deformation of Kaikōura Peninsula. On Figure 7 uplift for the MtB marine terrace (blue filled squares) and the 2016 uplift (grey squares) is normalised to the maximum values for each time interval which, in both cases, occurs at ∼29.5 km distance. The uplift values for the 2016 and MtB profiles approximately overlap at distances southeast of 32.2 km, with 2016 uplift (and its associated 2σ uncertainty bounds) higher than Late Pleistocene values northwest of this point (Figure 7). The difference in uplift profiles increases towards the northwest because northwest tilting of the Late Pleistocene (MtC and MtB) marine terraces is significantly larger than the first-order north-westward decrease in uplift in 2016. The relatively low 2016 uplift gradient along Kaikōura Peninsula between the ABF and TTF faults ( Figure 4B and Figure 7) could be due to a number of factors. First, a greater proportion of the total slip may have occurred on deeper parts of the OSTF in 2016 compared to the Late Quaternary. This possibility is partly borne out by the normalised slip profiles in Figure 5A & B (top graphs) which show more slip at >10 km depth in 2016 compared to the Late Quaternary. However, the depth distribution of slip is model dependent and may require additional work to conclusively demonstrate this point. Second, the low uplift gradients between the ABF and TTF in 2016 (compared to the Late Quaternary) could arise because during 2016 these faults ruptured in conjunction with the OSTF to produce an uplifted block, with such synchronous rupture on the ABF and TTF faults being rare in the Late Quaternary. Independent of exactly how the slip was distributed on the causal faults, it is clear that slip during the Kaikōura Earthquake in the Kaikōura Peninsula study area did not replicate the general long-term pattern of faulting in this region (Figure 7). This difference suggests that 2016 fault slip in the Kaikōura Peninsula study area was more complex and included more faults than was typically the case in the Late Quaternary. The observed complexity of fault rupture during large magnitude earthquakes is dependent on the number of faults that rupture and the resolution of the technique used for fault-rupture mapping. The resolution of fault mapping for historical earthquakes depends on a range of factors, including the depth of faulting, the type and extent of vegetation, the date of the earthquake and the number of techniques available for resolving the fault slip (Duffy et al. 2013;Nicol et al. 2016;Howell et al. 2020). None of the three faults mapped and modelled in this paper were observed to rupture the ground surface in the Kaikōura Earthquake (cf. Litchfield et al. 2018). The reasons for this are two-fold. First, the OSTF does not appear to correlate with the Point Kean Fault and there is presently no evidence that it ruptured the ground surface or the sea-bed in 2016 (further seabed imagery may be necessary to confirm the lack of seabed rupture). The apparent absence of a seabed scarp might partly be due to some fault slip being accommodated by folding in the near surface (e.g. <5 km depth). Despite the absence of a seabed scarp, relocated aftershocks (Mouslopoulou et al. 2019;Chamberlain et al. 2021), and modelling using coastal uplift and tsunami data Power et al. 2017;Gusman et al. 2018;Mouslopoulou et al. 2019; this study), provide strong evidence for rupture of a thrust fault beneath Kaikōura Peninsula in 2016. However, the rupture dimensions and slip magnitude remain poorly resolved compared to faults that ruptured the ground Figure 7. Normalised displacements along profile A (see Figure 5 for location) for the 2016 Kaikōura earthquake as derived from the lidar data (grey dots) and for Late Quaternary uplift as determined from the MtB marine terrace (blue line and filled squares). Grey polygon enclosing the lidar data defines the estimated 2σ uncertainties on the normalised uplift, as estimated from the Dlidar dataset. Dashed blue line shows the inferred normalised uplift for the MtB marine terrace northwest of State Highway 1. The location of the dashed line is based on interpretation of Figure 3, and is inferred as studies have not been conducted to verify the marine origin or the age of the surface.
surface. Second, the ABF and TTF did not rupture the ground surface and have vertical displacements too small (and/or spatially distributed) to be routinely resolved by field studies and remote mapping analyses (e.g. Nicol et al. 2018;Zinke et al. 2019;Howell et al. 2020). For the Kaikōura Earthquake, we suggest that most faults with >1 m surface displacements were at least partly resolved by post-earthquake field mapping, however, it was not always possible to unambiguously map faults that did not rupture the ground surface or have surface displacements of <1 m. Our data highlight the likelihood that the complex network of faults that were identified as rupturing in 2016 (e.g. Litchfield et al. 2018), represent a minimum number of faults and that the actual ruptures may be more numerous and complex than those already mapped by field observations (e.g. Litchfield et al. 2018;Nicol et al. 2018;Williams et al. 2018) or by remote techniques (e.g. Zinke et al. 2019;Howell et al. 2020). It is however worth noting that the unmapped faults likely represent a small proportion (e.g. <10%) of the total seismic moment generated in the earthquake and, in this regard, may not be a significant component of the earthquake.
The Kaikōura Peninsula is located between the northern ruptures (in the Marlborough Fault System) and the southern ruptures (in the North Canterbury Tectonic Domain) of the Kaikōura Earthquake, in an area where few faults were mapped to have ruptured the ground surface in 2016 (Litchfield et al. 2018). The paucity of faults at the surface in the Kaikōura area may partly reflect the accommodation of slip at depth on the blind OSTF Gusman et al. 2018;Mouslopoulou et al. 2019). We have modelled slip on three faults (ABF, TTF & OSTF) that intersect at depth, while the OSTF is modelled to extend to the subduction megathrust ( Figure 5). South of the peninsula the ABF and TTF are projected  Chamberlain et al. (2021). Faults in the 3D model B have been truncated at 20 km depth which is the approximate lower limit of the majority of the aftershocks used to define the fault planes (Mouslopoulou et al., 2019;Chamberlain et al. 2021) . In the cross section C the fault planes are projected below 20 km depth to the plate interface. In both the 3D model and the cross section the dip of the OSTF is 35°, consistent with the dip used in the forward modelling.
to intersect in map-view close to the northern wall of the Kaikōura Canyon, while 3D models of the OSTF (or its modelled equivalent) show it extending to the northern edge of the canyon (e.g. Mouslopoulou et al. 2019) (Figure 1B). At their closest map-view point the OSTF and Hundalee Fault surfaces may be separated by as little as 6 km along strike at the seabed. Given the proximity of the faults, that the maximum displacement of ∼4 m on the Hundalee Fault occurs at the coastline (Williams et al. 2018), and that the strike and dip of the OSTF and Hundalee Fault are comparable (compare Mouslopoulou et al. 2019, slip vectors calculated by Zinke et al. 2019), we propose that the OSTF and Hundalee Fault form part of a single larger fault (Figure 8). The Hundalee-OSTF structure could have a strike length of 70-100 km and account for at least 25% of the total seismic moment for the Kaikōura Earthquake (Mouslopoulou et al. 2019; this study). The Hundalee-OSTF structure extends north of Kaikōura Peninsula where it has been modelled to link at depth with the Papatea and Hope faults (Mouslopoulou et al. 2019) (Figure 8). Therefore, we propose that the OSTF hard links the southern and northern rupture domains and, as such, the faults that ruptured in 2016 are part of a fully connected network (Gusman et al. 2018;Mouslopoulou et al. 2019;Chamberlain et al. 2021) (Figure 8). In this hard-linked fault model the Papatea Fault connects, and transfers displacement between, the OSTF and the Kekerengu Fault, which may help explain why the Papatea Fault has such a high average displacement (∼6.4 m; Langridge et al. 2018) for its length (∼19 km) compared to the international literature (e.g. Wesnousky 2008). In addition to representing a key linking structure, the Hundalee-OSTF fault forms the southernmost thrust in the subduction complex, and is the only Hikurangi subduction-related thrust fault that can be mapped both offshore and onshore. Given its potential tectonic significance at the southern termination of the subduction system and onshore-offshore location, this fault warrants further study.

Conclusions
We have measured uplift along the Kaikoura Peninsula during the 2016 earthquake and in the Late Quaternary from marine terraces. D-lidar analysis of uplift in 2016 resulted in the discovery of two new blind faults (Armers Beach Fault, ABF; Te Taumanu Fault, TTF) which produced local changes in uplift of 0.2-0.5 m. Modelling of the stepped 2016 uplift indicates 2.3 m mean-slip on the offshore thrust and 0.25-0.5 m mean-slip on both the ABF and TTF structures. The pattern of uplift during the 2016 earthquake appears to differ from Late Quaternary uplift of MIS5c (96 ± 5 ka) and MIS5e (123 ± 5 ka) marine terraces. The MIS5c marine terrace is vertically displaced 23 ± 5 m by the ABF, while both terraces are tilted to the northwest. Northwest tilting of the marine terraces is primarily attributed to ∼80-230 m of slip and a slip rate of ∼0.8-2.7 mm/yr on the OSTF during the Late Quaternary. Slip on the ABF contributes to tilting and displacement of the marine terraces, while slip on the TTF is not required to account for the marine terrace deformation. Normalised uplift profiles along the axis of the peninsular indicate that northwest tilting of the marine terraces is greater than the 2016 uplift gradient. We attribute the lower uplift gradient in 2016 to slip on the TTF, which appears to have ruptured too infrequently in the Late Quaternary to impact the tilting of the marine terraces. The OSTF and Hundalee Fault could form a single structure that links the 2016 ruptures north and south of Kaikōura, producing a hard-linked network of faults.