Tectonic inheritance during plate boundary evolution in southern California constrained from seismic anisotropy

,


Introduction
Southern California hosts one of the world's best-instrumented and most thoroughly studied transform plate boundaries, yet questions as to how surface deformation transitions to convective flow and deformation at depth remain.Present-day transform motion is accommodated on the San Andreas Fault (SAF) (Fig. 1) and the major faults strands to the west, such as the San Jacinto Fault and offshore faults of the Continental Borderland; as well as to the east within the Eastern California shear zone (Fig. 1) (e.g.Bennett et al., 1996;Bourne et al., 1998;Meade & Hager, 2005;Becker et al., 2005;Chuang & Johnson, 2011;Y. Zeng & Shen, 2016, 2017;E. H. Hearn, 2019).The SAF itself forms a large-scale restraining bend in southern California such that many active faults are rangebounding thrusts and/or oblique-slip faults with potentially complex geometries and intersections at depth (e.g.Yule & Sieh, 2003;Matti et al., 1993Matti et al., , 1992)).The tectonic history of the region encompasses long-lived Mesozoic subduction, followed by extensional and transform regimes (Atwater, 1970).Subduction was accompanied by vigorous arc  (Plesch et al., 2007).SN, Sierra Nevada; WSNF, western Sierra Nevada foothills; SAF, San Andreas Fault; KCF, Kern Canyon Fault; WWF, White Wolf Fault; CR, Coast Ranges; WTR, Western Transverse Ranges; GF, Garlock Fault; ECSZ, Eastern California shear zone; SBM, San Bernardino Mountains, which form part of the ETR, Eastern Transverse Ranges; SJF, San Jacinto Fault; PR, Peninsular Ranges.
magmatism and emplacement of batholiths that now compose the Sierra Nevada and Peninsular Ranges (Fig. 1).These appear to behave as rigid blocks within the present-day kinematic reference frame, but nonetheless preserve Mesozoic magmatic fabrics and fossil ductile shear zones that may influence interpretations of geophysical images.Some subregions underwent significant rotation in addition to translation following the Mesozoic, notably an up to 110 • rotation of the Western Transverse Ranges block (Nicholson et al., 1994;Atwater & Stock, 1998;McQuarrie & Wernicke, 2005).
Given the complex tectonic history of the region, it remains a challenge to link different types of geophysical datasets and images to aspects of the present-day and/or longterm deformation regimes.In particular, there are several open questions, including: a) Do crustal and lithospheric mantle fabric reflect the present-day stress state?What is the timescale for resetting lithospheric rock fabric?b) Is strain localized on currently active faults and shear zones, and if so, to what distance and depth?Does strain localization cross the brittle-plastic transition?Is lithospheric fabric reset below this transition?Do fault communicate laterally below the brittle-plastic transition?c) How do asthenospheric processes relate to surface stress?If deep processes drive near-surface stress state, is lithospheric fabric reset to match both in response?
Here we attempt to address some of these questions through synthesizing a wide range of geophysical datasets that provide information about fabric and stress state in southern California, from the surface through the asthenosphere.There have been a number of studies comparing some of the relevant data sets to address these issues with each other already, in particular SKS splitting and focal mechanisms to GNSS velocities to infer shallow stress (Polet & Kanamori, 2002;Chamberlain et al., 2014;Yang & Hauksson, 2013) as well as surface wave results (Kosarian et al., 2011), but a comprehensive analysis across all accessible depth ranges as presented here seems to be missing.In particular, we seek to gain resolution in the lithosphere via inclusion of additional data from receiver functions and P and P n tomography.We begin by providing an overview of the geophysical data sets that can be used as stress and strain markers at different crust and mantle depths (Fig. 2) and previous results from the literature.We then describe our methods and discuss the results of cross-comparisons in light of the tectonic questions listed above.
2 Overview of geophysical datasets that can be used as stress and strain markers in southern California Geodetic constraints.Present-day surface deformation on decadal scales can be characterized by geodetic constraints (e.g.Haines & Holt, 1993;Flesch et al., 2000;Kreemer et al., 2000;Wei et al., 2010), and gradients of GNSS velocities, for example, yield estimates of the horizontal strain rate tensor at the surface.Such inferences often differ greatly in amplitude, in large part due to choices of interpolation and if and how fault localization is included (e.g.Haines et al., 2015;Sandwell & Wessel, 2016).Strain rates are also affected by co-and postseismic effects, leading to possible regional temporal bias -4-manuscript submitted to Geochemistry, Geophysics, Geosystems (e.g.Hetland & Hager, 2006;Chuang & Johnson, 2011;E. H. Hearn et al., 2013).Despite the uncertainties, however, the broad, ∼ 50 km length scale patterns of geodetically inferred surface strain rates for our study region are fairly consistent among models (e.g.Sandwell et al., 2009;Kreemer et al., 2014) and we use the SCEC Community Geodetic Model V1 for geodetic strain-rates (Sandwell et al., 2016) in the following.This model consists of an average of a range of GPS-based estimates, but preserves spatial variations on scales similar to those imaged from seismicity, discussed next.
Focal mechanisms or moment tensors.At seismogenic depths, recent and ongoing deformation can be inferred from focal mechanisms or moment tensors.These provide direct constraints on co-seismic strain only, and compressional stress orientations can in theory lie anywhere within the compressional quadrant of a single mechanism (McKenzie, 1969).There are a range of approaches to convert a set of mechanisms to a deformation tensor; for example relatively straightforward Kostrov (1974) summation for strain(-rates), or normalized stress-tensor inversions that often assume the direction of slip on the rupture plane is parallel to the orientation of the maximum shear stress (e.g.Michael, 1984).
We know that stress inversions are affected by structural heterogeneity (e.g.Gephart & Forsyth, 1984;Pollard et al., 1993;Townend & Zoback, 2004;Hardebeck & Michael, 2004) and different approaches yield different estimates on length scales 10 km (e.g.Abolfathian et al., 2020).However, Kostrov summed estimates of strain broadly agree with smoothed stress inversions on longer wavelengths (Hardebeck & Michael, 2006;Becker et al., 2005).We therefore use the stress inversion of Yang and Hauksson (2013), which is of the smoothed Hardebeck and Michael (2006) type and has been adopted as a working model for crustal stress by SCEC.In the comprehensive comparison, we also use our own Kostrov summation of focal mechanisms (Yang et al., 2012), as used by Yang and Hauksson (2013) with a temporally more limited catalog, and simple binning with 0.25 • boxes.Given the comparison with seismic anisotropy below, we here focus on the horizontal, orientational quantities that can be derived from strain (rate) or stress tensors, i.e. the major compressional or extensional axes of the principal component system, and assume temporal stationarity.
Local event shear-wave splitting.A more indirect measure of the shallow crustal stress is given by splitting of shear waves from local earthquakes, arriving on steep paths to stations above.The shear wave component polarized parallel to cracks will propagate faster than its orthogonal counterpart; since stress will control crack closure and formation, this leads one to expect fast polarization orientations that match the compressive stress orientation (e.g.Crampin & Chastin, 2003).The most recent complete analysis of local splitting in the study region is by Z. Li and Peng (2017) who note that this expected relationship breaks down over significant portions of the region, including in the vicinity of major strike-slip faults.Comparison of local shear wave splitting with stress inferences from boreholes indicates that complexities are expected, for example because of shallow fabric rather than crack-stress control (Boness & Zoback, 2006).There are other, more direct indicators of crustal stress from boreholes, but those are not available widely enough, nor are they sampling on the appropriate scales for our study (Persaud et al., 2020;Luttrell & Hardebeck, 2021).
Receiver Functions.Receiver function arrivals show strong variations in arrival amplitude and polarity with backazimuth in our study region, and arrivals with a 360 • dependence with backazimuth (first azimuthal harmonic) can be interpreted as conversions from contrasts in crustal seismic anisotropy with a plunging symmetry axis, such as dipping foliation (Porter et al., 2011;Schulte-Pelkum, Ross, et al., 2020;Brownlee et al., 2017).Unlike splitting methods such as local event shear wave, SKS, or receiver function Moho P -to-S conversion (P s) splitting analyses, this method is not a cumulative measurement of shear wave anisotropy and instead is sensitive to changes in P anisotropy (Levin & Park, 1998;Bianchi et al., 2010;Park & Levin, 2016).It therefore allows resolving the depth of anisotropic contrasts and is sensitive to shear zones with thicknesses -5-manuscript submitted to Geochemistry, Geophysics, Geosystems 2 km (Schulte-Pelkum & Mahan, 2014b, 2014a;Liu & Park, 2017).Azimuthally varying receiver function conversions have higher amplitudes and therefore better sensitivity for plunging axis symmetry, e.g.dipping foliation, compared to a horizontal symmetry axes or purely azimuthal anisotropy (Schulte-Pelkum & Mahan, 2014b;Park & Levin, 2016;Brownlee et al., 2017).Studies modeling observed receiver function arrivals in the southern and central California region find plunging symmetry axes in the crust, including at middle and lower crustal depths, with foliations paralleling major faults (Audet, 2015) or perpendicular to Farallon convergence (Porter et al., 2011).A simpler method that determines foliation strike without waveform modeling (Schulte-Pelkum & Mahan, 2014b, 2014a) shows surface fault-parallel strikes with a pervasive lithospheric fabric that was interpreted as tectonic inheritance from previous compressional and extensional episodes (Schulte-Pelkum, Ross, et al., 2020).
Regional phase P n and teleseismic P wave tomography The regional seismic phase P n propagates horizontally as a headwave in the uppermost mantle immediately beneath the Moho.Azimuthal variations in its propagation speed are used to measure horizontal axis seismic anisotropy (T.M. Hearn, 1996;G. P. Smith & Ekström, 1999;Buehler & Shearer, 2010, 2012, 2014, 2017), reflecting lithospheric mantle fabric.Studies either solve for isotropic and anisotropic structure using tomography (T.M. Hearn, 1996;Buehler & Shearer, 2010, 2014, 2017) or analyze localized azimuthal variations using station subsets (G.P. Smith & Ekström, 1999;Buehler & Shearer, 2012).The regional fast axes are broadly subparallel to the strike of the SAF, albeit with geographic variations.Joint tomography of local and teleseismic P arrival times can also be used to infer seismic anisotropy (e.g.Bokelmann, 2002), though it is subject to trade-offs between isotropic and anisotropic structure.Consistent with the P n results, a recent, regional anisotropic P wave tomography model shows broadly SAF-parallel fast orientations at lithospheric depths (Yu & Zhao, 2018).
Surface wave tomography.Surface wave tomographic studies in the region resolve isotropic heterogeneity as well as radial anisotropy, with anomalous regions of negative radial anisotropy bordered by the SAF (K.Wang et al., 2020).Continent-wide azimuthal anisotropy studies using surface waves show fast orientations somewhat similar to SKS in the upper mantle and variable anisotropy in the crust (Lin et al., 2011).Regional-scale studies find evidence of azimuthal anisotropy in phase velocities in the area (Qiu et al., 2019), though not in some more localized studies such as for the Parkfield area (X.Zeng & Thurber, 2019).A surface wave inversion for 3-D azimuthal anisotropy at regional resolution is not available to date.
SKS splitting.In addition to local event shear wave splitting, there are a multitude of studies that target splitting of much lower-frequency teleseismic core phases such as SKS (Savage & Silver, 1993).Similar to the case of local event splitting, the teleseismic splitting time delay and fast axis orientation are integrated nonlinearly (Silver & Savage, 1994;Silver & Long, 2011;Kaviani et al., 2011) over the entire raypath, in this case from the core-mantle boundary to the station.The signal is typically interpreted as being dominated by the asthenospheric mantle, with some contribution from the lithosphere and SAF shear (e.g.Savage et al., 2004) and small to negligible influence of the crust (Silver, 1996).While the stress field has nevertheless been invoked to explain SKS splitting in our study region (Polet & Kanamori, 2002), a more common explanation is large-scale mantle flow and asthenospheric convection leading to lattice preferred orientation olivine fabrics (e.g.Savage & Sheehan, 2000;Silver & Holt, 2002;Becker et al., 2006;Bonnin et al., 2010;Ramsay et al., 2016;Zhou et al., 2018).Alignment with shear as expected from absolute plate motions does not match the observed orientations well (Silver & Holt, 2002;Bonnin et al., 2010), but modeling of plate motion and density-driven mantle flow, without any deep shear localization at the plate boundary, captures SKS patterns on scales of 200 km (Becker et al., 2006).Nonetheless, small-scale variations in teleseismic splitting across the SAF are observed ( Özalaybey & Savage, 1995;Savage et al., 2004;Jiang et al., 2018), and the role of the SAF in affecting SKS splitting remains debated (e.g.Savage et al., 2004;Bonnin et al., 2010).Depth-dependent splitting studies show SAF-parallel fast orientations in the lithosphere on the northern part of the study area (north of the Garlock fault; Fig. 1), although solutions may be ambiguous (Monteiller & Chevrot, 2011;Hartog & Schwartz, 2001).
In the following, we synthesize these geodetic, geological, and seismological deformation and stress and strain markers by quantitative comparisons of their orientations across the region.

Methods
We select published oriented data sets and compilations.Pairs of data sets are compared as shown in Fig. 3 as an example.Orientations of collocated pairs of data with amplitude and strike information are plotted against each other in map view as in Fig. 3a.
A second view for each comparison pair shows the angular deviation on a local footprint as background color (Fig. 3b) to highlight geographical variations in angular agreement or disagreement.Statistics of the angle difference over the entire area in which data from both methods are available are calculated and plotted in an inset (example Fig. 3b).For this study, we assume all data with azimuthal information to be orientational, or axial, rather than vectorial (e.g., N-S, rather than N or S).Most data used in this study are indeed axial (fast polarization or propagation orientations, stress and strain axes, shortening or extensional axes).Fault orientations are reduced to unsigned strike, ignoring dip direction.Azimuthal receiver function arrival strikes are signed and offer dip direction information if additional assumptions on the layers bounding a contrast are made, such as sign of a velocity contrast, whether anisotropy is stronger above or below an interface, or whether anisotropy has a fast or slow symmetry axis (Schulte-Pelkum, Ross, et al., 2020).In this study, we use unsigned strikes for receiver functions to compare to other axial quantities, but display the signed orientation in Fig. 6c; the short arrows point towards downdip if the signal is from the top of an anisotropic layer with slow axis symmetry (such as schists; (Brownlee et al., 2017) or from a slow-over-fast dipping interface between isotropic layers.
To test whether specific datasets are consistent with the kinematic regime of San Andreas transform right-lateral shear, we plot them relative to each other as well as relative to two possible reference orientations coinciding with the instantaneous and finite strain ellipses, assuming simple shear strain geometry.The right-lateral motion along the Pacific-North America plate boundary is oriented NW-SE, with the Pacific plate moving NW relative to North America.Instantaneous shear strain in a medium deforming under distributed viscous simple shear deformation is expected to show compression and extension at 45 • to the shear plane, in this case, N-S compression and E-W extension.
During progressive simple shear, however, the lengthening orientation rotates into the plane of shear (e.g.McKenzie & Jackson, 1983), so datasets that potentially track high strain fabrics such as ductile foliations or gouge fabrics along faults may show alignment with the finite strain ellipse (SAF-parallel, NW-SE).
We also compare some datasets to the inferred orientation of Farallon paleofoliations.Farallon convergence was oriented NE relative to the North American margin, and we assume paleofoliations formed perpendicular to Farallon convergence orientations, consistent with the general orientation of major Mesozoic-to-early-Tertiary terrane boundaries, faults, and ductile foliations exposed at the surface presently (Ernst, 1970;Hamilton, 1969;Dickinson, 1970;Jacobson et al., 1996).In most of the study area, these paleo-

Mantle/asthenospheric depths
In Figure 3, we compare fast axis orientations from SKS splitting (i.e.inferred orientation of the polarization plane of the fast propagating pulse) from a compilation (update of Becker et al., 2012) to the uppermost mantle layer solution of an isotropic and anisotropic Rayleigh wave inversion using ambient noise as well as earthquake phase dispersion (Lin et al., 2010).Sensitivity of SKS extends from the core-mantle boundary to the surface but is assumed to be dominated by the asthenospheric mantle.The surface wave inversion has depth constraints down to ∼ 100 km (Lin et al., 2010).This can be compared to an estimated LAB depth ranging mostly from 60-80 km (half-width of distribution) in the area inferred from Sp receiver functions (Lekic et al., 2011;Levander & Miller, 2012).
Both sets of observed fast axes are rotated counterclockwise from the SAF (Fig. 4 and Supplementary Information).The misorientation of SKS relative to the average SAF strike outside of the Big Bend is −52 (Fig. 3b).The P tomography study by Yu and Zhao (2018) shows only low-amplitude azimuthal anisotropy at depth nodes of 125 and 200 km.Cross-comparisons plots for all datasets and models can be accessed through Fig. 4 which also summarizes the median, standard deviation, as well as the degree of unimodality of the angular misfit distributions.

Lithospheric mantle
Depth-specific constraints for the lithospheric mantle are extracted from P n tomography, localized P n azimuthal velocity variations, and P tomography.While some anisotropic receiver function signals are from structures below the Moho, most are from crustal depths, and receiver function inferred strikes will be discussed in the next section.inset and strike bar fill color in Fig. 6c).The sign of the azimuthal polarity flips contains information about foliation dip at a contrast or interface dip, but the dip sense is ambiguous.The phase of the A 1 amplitude signal is shown as short one sided arrows in Fig. 6c.
-9-manuscript submitted to Geochemistry, Geophysics, Geosystems   (Plesch et al., 2007) in Fig. 6c and d.The angular deviation histograms (insets in Fig. 6b and d) are peaked at 0 • and −2 • deviation, respectively.A rotation of receiver function strikes from SAF-parallel to E-W in the Western Transverse Ranges is matched in both cases.The comparison to Farallon paleofabric shows a higher peak at zero deviation, in part due to station distribution such as the dense Salton Sea experiment (Klemperer, 2011;Barak et al., 2015) showing consistent strikes in the notably presently internally undeformed Peninsular Ranges block near the international border.
The comparison to local fault strikes shows reduced scatter such as the improved matches in the northern border of the Western Transverse Ranges and other localized strikes such as the White Wolf Fault (Fig. 1).A comparison of receiver function strikes to average SAF strike shows large deviations in the Western Transverse Ranges and in the Eastern California shear zone (Fig. 4).Despite the concentration of receiver function strikes above 20 km depth, there is little alignment of the strikes to local event shear wave splitting fast axes (Fig. 4), suggesting that microcracks are not a dominant mechanism for generating anisotropy detected by A 1 arrivals in receiver functions.
As noted in Schulte-Pelkum, Ross, et al. ( 2020), vertical foliation that may be expected from subvertical shear along transform faults would produce a 180 • -periodic arrival with backazimuth (A 2 ) rather than the A 1 signal considered here.In our study area, the A 2 signal amplitude is smaller than that of than A 1 (Schulte-Pelkum, Ross, et al., 2020).If we assume that Farallon paleofoliations have gentle to intermediate dips and that foliation reset by SAF-age transform motion should be subvertical and possibly concentrated along the SAF, then the strong and geographically distributed A 1 signal appears to favor paleofabric.

Seismogenic crust and surface
As noted, local shear wave splitting is thought to show fast polarization axes parallel to maximum compressive stress due to compression-parallel orientation of microcracks.We compare local S splitting with compressive stress from a focal mechanism inversion (Yang & Hauksson, 2013) in Fig. 7a and b.As identified by Z. Li and Peng (2017), there is significant misorientation to the maximum horizontal compressive stress axis along the SAF and other major transform fault strands as well as in additional regions such as the Peninsular and Coast Ranges.In most of these areas of misalignment, local S fast axes are closer to Farallon paleofabric (Fig. 6c).
In contrast, a comparison of the same focal mechanism-based maximum compressive stress axes with GNSS-derived compressive strain rates (Sandwell et al., 2016) shows a close match, with a strongly peaked distribution centered on −6 • (Fig. 7c,d).Geodetic inferred deformation and focal mechanism-derived stress orientations agree well with a regional, roughly N-S compressional/E-W extensional state (Fig. 7c,d).Extensional horizontal strain based on surface geodesy as well as extensional stress from focal mechanisms are roughly E-W, ∼45 • counterclockwise from average SAF strike (Fig. 8a,b).
This is broadly what is expected if the whole region were in a stress and strain-rate state with the major compressive axes oriented N-S and assuming the crust deforms as a viscous continuum.When inferred from focal mechanism inversions, that average azimuthal value is close to 8 • (Yang & Hauksson, 2013).Surface instantaneous strain as well as focal mechanisms therefore are as expected for the present transform motion, an alignment that is borne out on smaller scales for most faults in the area besides those recently affected by major earthquakes (Becker et al., 2005).
Local shear wave splitting results sample the same depth range as focal mechanisms, but as discussed in Section 2, they have the potential to be parallel to either the compressive stress orientation (in the case of microcracks) and/or the direction of finite strain (in the case of fault gouge fabrics).In some parts of the study area, fast local S splitting polarizations match the focal-mechanism-based stress and GPS-derived strain (∼N-S fast, parallel to maximum compressive horizontal stress).However, there are significant areas of misorientation along, for example, the SAF and San Jacinto fault and in internally undeformed blocks such as the Peninsular Ranges (Fig. 7a,b), as discussed by Z. Li and Peng (2017).For these areas of mismatch to the stress field, we hypothesize that shallow splitting observations are dominated by paleofabric (Coastal Ranges, Peninsular Ranges; Fig. 6a) and/or present-day fault-parallel (Fig. 6c) fabric.Laboratory ultrasonic studies suggest that microcracks in the uppermost crust above crack closure pressure depths may orient themselves along preexisting mineral fabric (Rasolofosaon et al., 2000).We thus propose that shallow splitting of local event shear waves is partly controlled by current stress state, but also partly by rock fabric, some of the latter inherited since it is found to dominate splitting in presently internally non-deforming blocks.
Proceeding to the entire crust past seismogenic depths and into the creeping deformation regime, we find that receiver function A 1 arrivals show a pervasive fabric that deviates from the near-surface E-W extension/N-S compression orientation (Fig. 8c), with foliation strikes paralleling present-day surface fault traces as well as paleofabric strikes (Fig. 6).It is unlikely that the fabric imaged by receiver functions is strictly due to presentday faulting.Since the inferred fabric is distributed widely geographically (no clear dependence with distance from faults) and in depth, Schulte-Pelkum, Ross, et al. (2020) previously interpreted the receiver function signal as inherited fabric.This conclusion is bolstered by the match between receiver function strikes and paleofabric (Fig. 6).
The Moho depth in the study region mostly varies between 25-32 km, with a few local Moho roots reaching depths of 35-40 km under the western foothills of the Sierra Nevada, San Bernardino mountains, and Peninsular Ranges (Zhu & Kanamori, 2000;Yan & Clayton, 2007;Miller et al., 2014;Ozakin & Ben-Zion, 2015).Receiver function A 1 station maxima are concentrated at crustal depths, but reach into lithospheric mantle depths (Fig. 6a) with similar strikes to those at shallower depths (Schulte-Pelkum, Ross,  et al., 2020).Changes in seismically imaged anisotropic symmetry can be due to changes in mineralogy and do not necessarily imply a change in deformation regime (Bernard & Behr, 2017;Brownlee et al., 2017), and the converse also holds.Nevertheless, the similarity between receiver function and P n inferred fast strikes (Fig. 5) suggests similar fabrics in the deep crust and uppermost lithospheric mantle.
Like receiver functions strikes, P n fast orientations deviate from shallow E-W/N-S orientations (Fig. 8d).The area showing SAF-parallel P n fast axes is broad and shows no clear dependence on distance to the SAF.No rotation of P n fast axes are seen in the compressive Big Bend in the Eastern Transverse Ranges.In contrast, the Western Transverse Ranges show P n strikes consistent with compressive Farallon-age paleofabric that was subsequently rotated.This contrast between the western and eastern part of the Big Bend as well as the lack of localization of P n anisotropy relative to major fault strands supports inherited fabric rather than present day deformation as a source of lithospheric mantle fabric.Teleseismic P wave derived anisotropy at 75 km depth supports this interpretation (Yu & Zhao, 2018) (Fig. 4).Tomographic inversions solving for isotropic and anisotropic structure are in principle subject to trade-off effects between the two types of anomalies.However, the trade-off tests presented in the tomographic studies (Buehler & Shearer, 2010, 2012, 2014, 2017;Yu & Zhao, 2018) and the agreement of tomographic P n results with those using azimuthal variations between station pairs (G.P. Smith & Ekström, 1999;Buehler & Shearer, 2012) lend some confidence to the results and our interpretation.While isotropic P n velocity models may be impacted by local Moho topography, it is less likely that the fast axes are distorted significantly.
We thus interpret the crustal and lithospheric mantle deformation indicators as a pervasive lithospheric fabric that may largely reflect formation processes (accretion and intrusions during long-lived Farallon subduction) rather than recent deformation processes.Wholesale rotation of fabric at lithospheric mantle depths along with their block surfaces such as in the Western Transverse Ranges and a lack of effects of recent processes such as compression in the Eastern Transverse Ranges on deep crustal and lithospheric mantle fabric suggests that deep fabric is preserved rather than reset.
At deeper levels in the mantle, surface wave tomography and teleseismic splitting show a close correspondence and are not parallel to either paleofabric or SAF strike orientation.SKS fast axes deviate more from SAF strike (Fig. 8f) than surface wave fast azimuths at periods sensitive to the mantle (Fig. 8e).It is likely that SKS is sensitive to greater depths (asthenosphere) than the surface wave mantle solution (Lin et al., 2010), and the latter may be averaging over lithospheric as well as deeper anisotropy.Although the roughly E-W orientation for SKS matches surface extension (Fig. 8a), modeling suggests it is due to asthenospheric broader-scale mantle circulation (Becker et al., 2006).
Whether and how this mantle circulation, mostly due to sinking of the formerly subducting slab into the mantle, may control the present-day near-surface stress state remains an open question (Humphreys & Coblentz, 2007;Ghosh et al., 2013;Chamberlain et al., 2014;Kosarian et al., 2011).Short-wavelength perturbations to SKS likely reflect influence from shallower, lithospheric structure (e.g.Savage et al., 2004;Bonnin et al., 2010;Jiang et al., 2018), but may also be associated to smaller length scale convection than what was considered by Becker et al. (2006) and/or reorientation of olivine fabrics (e.g.Zhou et al., 2018;W. Wang & Becker, 2019).
We show histograms of the angular misalignment as in Figs.3-7 for all quantities discussed above as well as additional axial data sets in Fig. 9.This comparison is further expanded and summarized in Fig. 4 whose quantification of the statistical moments can, however, only incompletely capture the character of alignment where, for example, a block-wise match of orientations can leads to a smeared out distribution (Fig. 9).Each comparison pair appears twice in Fig. 9, and we show median and standard deviation as well as area of coverage in addition to the angle difference histogram.The quantities compared here are the same as ones discussed above.The comparisons show the same  trends of agreement between surface strain and focal mechanism-based sets with deeper mantle estimates from SKS and surface waves on the one hand, and between lithospheric estimates (receiver function studies, P n , lithospheric depths from P tomography) on the other hand.We conclude that even if mantle strain is transmitted through the lithosphere to the surface, lithospheric fabric itself is preserved from prior tectonics and independent of current lithospheric deformation.
Our interpretation is summarized in near the SAF and in the Coast Ranges and Peninsular Ranges blocks.In the areas of misalignment with stress, the shallow splitting results may be interpreted as fault-parallel (SAF), but they may also show paleofabric, in particular in areas with little present-day internal deformation (Peninsular Ranges block).

Potential origins and implications of preserved paleofabric
Several lines of evidence support our interpretation that most levels of the crust and lithospheric mantle preserve inherited fabrics as opposed to fabrics that were reset by SAF system deformation.Firstly, although lithospheric orientations are broadly SAFparallel as well as Farallon paleofabric-parallel, the rotation of these strikes in the Western Transverse Ranges only matches the expected rotations of Farallon paleofabrics and is inconsistent with SAF-related shear.Secondly, the strength of the receiver function, If fabric is due to accretionary tectonics and schist underplating during Farallon convergence and subsequent extension, they may be expected to show shallow to intermediate dips (Chapman, 2017;Jacobson et al., 1996).In contrast, most strike-slip faults taking up the present-day transform motion are thought to have intermediate to vertical dips, and foliation from transform shear would be expected to be subvertical.Past receiver function studies in the area that attempted to constrain foliation dip angle (or--19-manuscript submitted to Geochemistry, Geophysics, Geosystems thogonal to symmetry axis plunge) through waveform modeling inferred a wide range of dips (Porter et al., 2011;Audet, 2015).
While regional-scale surface wave inversions for crustal azimuthal anisotropy are not yet available (Qiu et al., 2019), crustal radial anisotropy shows anomalous negative values (vertically polarized shear wave speeds higher than horizontally polarized shear wave speeds) in the mid-to lower crust, mostly south and west of the SAF (K.Wang et al., 2020).There are several possible explanations for this observation.Steeply dipping schists are one, although K. Wang et al. (2020) argue that the isotropic average velocities in the mid-to lower crust are too high for schists.Another possibility are mafic intrusive dikes.A third possibility is preservation of steep fabrics associated with intrabatholithic strike-slip faults activated during oblique subduction in the Late Cretaceous, as described for several surface exposures in the Sierra Nevada (Busby-Spera & Saleeby, 1990;Kistler, 1993).Although they seem localized in surface exposures, the receiver function strikes in the southern Sierra Nevada show broad consistency at lower crustal depths (Fig. 6c) and are parallel to the Kern Canyon fault (Fig. 1), which itself reactivated a Cretaceous ductile shear zone (Nadin & Saleeby, 2010;Busby-Spera & Saleeby, 1990).
Cretaceous transpressional strike-slip shear zones are also exposed in the eastern San Gabriel mountains within the Transverse Ranges (May & Walker, 1989).A fourth option for generation of negative radial anisotropy is suggested by K. Wang et al. (2020) based on xenolith work by (Bernard & Behr, 2017) showing that feldspars can align with a foliationperpendicular fast axis.In this case, the lower crust would have high isotropic velocities as well as radial anisotropy.This case would also match the receiver function observations if the foliation is gently dipping, since strikes inferred by the receiver function method are strikes of the plane perpendicular to the symmetry axis, regardless of whether that axis is fast or slow (Schulte-Pelkum & Mahan, 2014b; Brownlee et al., 2017).
Whether foliation is steeply or shallowly dipping, the presence of rotated strikes in the rotated Western Transverse Ranges block compared to dominant NW strikes in the unrotated Central and Eastern Transverse Ranges supports a paleofabric interpretation.
P n sampling the mantle just below the Moho shows fast axes that are broadly parallel to SAF and paleo-subduction strike.Similar results are given for 75 km depth in anisotropic teleseismic P tomography (Yu & Zhao, 2018); the neighboring solution grid depths are at 25 km and 125 km depth, and the 75 km depth slice shows the strongest mantle anisotropy.For P n , the prominent NW-SE strike is seen across a broad zone of 100 km width through the Big Bend region, showing no obvious change in azimuth along with the SAF.In contrast, P n fast axes are rotated in the Western Transverse Ranges block.This picture is consistent with previously proposed bottom accretion of a stack of Farallon lithospheric mantle slices (Luffi et al., 2009), and/or with strong shearing of North American-affinity mantle in response to flat-slab subduction, as suggested by several studies on western US mantle xenoliths (Z.-X. A. Li et al., 2008;Behr & Smith, 2016;D. Smith et al., 2004;Usui et al., 2003).Similar to receiver functions, the P n orientations appear to support a paleofabric interpretation over SAF motion-induced alignment.
The Transverse Ranges are related to the restraining bend in the SAF, and a number of related mantle anomalies have been proposed previously.A roughly contiguous isotropic Transverse Ranges fast upper mantle anomaly is imaged in body wave tomographic studies (Humphreys & Dueker, 1994;Schmandt & Humphreys, 2010;Yu & Zhao, 2018) and was interpreted as downwelling lithosphere, while surface wave tomography show a fast mantle root limited to the Western (K.Wang et al., 2018Wang et al., , 2020) ) or West--20-manuscript submitted to Geochemistry, Geophysics, Geosystems ern and Eastern (Barak et al., 2015) Transverse Ranges.A study that attempted to invert for depth dependent SKS splitting attributed a local path of E-W fast axes in the Big Bend area to local shortening (Monteiller & Chevrot, 2011).The downwelling interpretation would be contradicted by the strong azimuthal anisotropy seen in anisotropic P n (Buehler & Shearer, 2017) and P (Yu & Zhao, 2018) tomography.It appears difficult to reconcile all tomographic results from the Transverse Ranges area, but the uppermost mantle is also notoriously difficult to resolve.Irrespective of what the isotropic structure may be, the anisotropic results in conjunction with the crustal receiver function results appear to favor a paleofabric interpretation.
While it may be tempting to interpret the correspondence between E-W surface extension and E-W SKS splitting fast axes in a pure shear or instantaneous strain framework within SAF-parallel shear, mantle flow modeling studies (Steinberger, 2000;Becker et al., 2006;Zhou et al., 2018;W. Wang & Becker, 2019) show that slab sinking to the east of the study area is likely responsible for aligning asthenospheric fabric sampled by teleseismic splitting.

Conclusions
We compare orientations of stress-and deformation-associated observables such as surface deformation, stress state inverted from focal mechanisms, local event and teleseismic shear wave splitting, receiver function azimuthally varying conversions, P n and local/teleseismic P tomography, and surface wave azimuthal anisotropy, in addition to geological information such as fault strike and block rotations since Farallon subduction.
The deformation indicators separate into three classes, with a near-surface match between geodesy and focal mechanisms and some local event shear wave splits, a lithospheric depth range comprising receiver functions, P n , and local event splitting in other areas, and an asthenospheric component detected by mantle surface waves and SKS splitting that is driven by mantle circulation and can be speculated to transmit through the lithosphere to drive stress at the surface.Notably, lithospheric rock fabric and anisotropy do not appear to be reset to reflect present-day deformation and instead appear to persist since the time of accretion and intrusion during long-lived subduction.Fabric from the upper crust through the lithospheric mantle appears to have been preserved and rotated along with surface block rotations through temporal changes in deformation regimes.

Figure 1 .
Figure 1.Map of the study area with shaded topography and elements referred to in text (larger context shown in inset with main map area outline drawn).Black lines are faults from the SCEC Community Fault Model (CFM) version 5.3(Plesch et al., 2007).SN, Sierra Nevada;

mFigure 2 .
Figure 2. Perspective view of the study area with elevation and topography (surface color) and fault traces.Vertical dimension shows conceptual sketch of observables (blue) and their depth sensitivities (arrows) as well as presumed mechanisms (italic script) that may dominate azimuthally varying observations in different depth layers (color).
orientations have a similar strike to that of the SAF system.Exceptions are blocks that underwent post-Farallon rotations such as the Western Transverse Ranges (e.g.McQuarrie & Wernicke, 2005).One may expect deep-seated Farallon-convergence-related paleofoliations to have a shallower dip than transform shear fabric.However,Porter et al.

Figure 3 .
Figure 3. (a) Fast axis orientation and delay time of SKS splits (compilation of Becker et al., 2012, updated as of 2020) and Rayleigh wave upper mantle fast axis and percent anisotropy (orange bars; Lin et al., 2011).Elevation shading and faults are shown for orientation, with faults in green for surface traces and red for blind faults from the SCEC CFM5.3 as in Fig. 1.(b) Signed angular difference comparison between the two orientation sets (color coding of bars in angle sign definition as in (a); amplitude information is not used in this case).Inset shows histogram of the angular differences, with mean(median) and standard deviation indicated; N is number of measurements pairs per bin, N0 total number of pairs.

Figure 5
Figure5shows the angular alignment between P n fast orientations fromBuehler and Shearer (2017) and SAF-strike-averaged orientations outside of the Big Bend region applied to the entire study area.Agreement is close in a broad corridor (order 100-200 km width) around the SAF, with no localization close to the SAF.Notably, the P n fast axes do not change strike along with the SAF in the Big Bend area, and there is no indication of constrictional E-W strikes due to the restraining bend as proposed byMonteiller and Chevrot (2011).The area closest to the SAF trace that shows significant misalignment with average SAF strike is in the Western Transverse Ranges.To test the possibility that lithospheric mantle fabric is not controlled by presentday shear associated with the SAF system and instead persists from previous tectonic episodes, we compare P n fast axes to Farallon convergence-perpendicular strikes including block rotations since 36 Ma from McQuarrie and Wernicke (2005) with interpolations byPorter et al. (2011).These orientations are a proxy for lithospheric fabric strikes that developed under compression, accretion, and pluton emplacement during long-lived subduction and were subsequently rotated along with the surface of the respective blocks as proposed by e.g.Atwater and Stock (1998).Agreement is improved in the western part of the Western Transverse Ranges and similar to the SAF comparison case elsewhere.P anisotropy(Yu & Zhao, 2018) shows the highest degree of seismic anisotropy in the lithospheric mantle layer (75 km depth grid nodes; compare 60-80 km LAB depths fromLekic et al. (2011)), with fast azimuths fairly well aligned with SAF strike or Farallon paleofabric in the SAF corridor (Fig.4).P anisotropy at 75 km depth also shows a strong falloff of amplitudes in the Western Transverse Ranges.

Fig. 5 )Figure 5 .
Fig.5); k-n) P tomography fast axes(Yu & Zhao, 2018) at depths of 25, 75, and 125 km, respectively; o) mantle depth surface wave anisotropy(Lin et al., 2011) (cf.Fig.3); p) absolute plate motion (APM) orientations in the spreading-aligned reference frame(Becker et al., 2015); and q) SKS splitting compilation fast axes (update ofBecker et al., 2012) (cf.Fig.3).Clicking on any of the symbols in the PDF version of this article provides a link to web based cross-comparison plots such as in Fig.3.(To be implemented: For now, the reader can access all the plots from,

Figure 6 .Figure 7 .
Figure 6.(a) Strikes and amplitude of largest A1 harmonic receiver function arrival at each station (in orange), compared to Farallon paleofabric (in pink) as in Fig. 5c).Inset shows depth distribution of A1 arrivals.(b) Orientation comparison for (a) as in Fig. 3b.(c) As in (a), but here compared to local fault strikes from SCEC CFM 5.3 (Plesch et al., 2007) and colored by conversion depth.Short arrows show receiver function inferred up-or downdip direction, colored by the strike of that vector to allow easier geographic subsetting.(d) Orientation comparison for (c).

Figure 8 .
Figure 8.Comparison of six data sets in a common reference system.The reference orientation is average SAF strike rotated 45 • counterclockwise to match approximate SKS and surface extension orientations.Note that the background color scale used to represent axis misorientation is different from that used in Figs.3-7; purple/blue shades now show alignment, and we do not distinguish between CW and CCW rotations, i.e. |∆α| ∈ [0; 90 • ].(a) GPS extension as inFig.7c.(b) Focal mechanism extension as in Fig. 7a,c.(c) Receiver function strikes as in Fig. 6.(d) Pn fast axis as in Fig. 5. (e) Mantle surface wave fast axis as in Fig. 3.(d) SKS fast axes as in Fig. 3.

Figure 9 .
Figure 9. Angular difference histograms (−90 . . .90 • ) as in e.g.Fig. 3b for all data sets discussed above.The top right triangle has as diamonds the median deviation, with the diamond's size scaled with the inverse of the standard deviation and the diamond centered on the median.The lower left has the same coloring, but the square size scales with the area of coverage out of the whole study area.See Fig. 4 for an abbreviated representation of an extended comparison set.

Fig. 10 .
Figure 10.Schematic illustration of dominant mechanisms for strain and seismic anisotropy in the southern California lithosphere and asthenosphere.Only surface strain and focal mechanisms are in accord with compression/extension (black thin arrows) from relative plate motion (black bold arrows).Inherited fabric influences the lithosphere up into the uppermost crust (white arrow is sketched to show approximate orientation of Farallon convergence).Asthenospheric shear reflects large scale circulation with entrainment from a slab sinking far to the east (curved/dashed arrow).

P
n , and P signals does not show any obvious dependence with distance from major transform fault strands such as the SAF or San Jacinto fault.Third, azimuthally varying receiver function arrivals from all levels of the crust show strikes that are consistent with dipping foliation.The receiver function signal is not explained by purely vertical or horizontal foliation, but may be caused by foliation dips ranging from subhorizontal (∼ 20 • ) to steep.The phase of the receiver function signal is consistent with imaging the upper interfaces of anisotropic layers with dominantly NE-dipping fast foliation planes (e.g.schists).
Mantle depth, surface-wave derived fast axes are slightly less rotated relative to the average strike of the SAF (−36 • ±18 • , peak 23%).The average angle difference between SKS -fast and mantle surface wave-fast axes is 12 • ±15 • at at a strong (23%) unimodal peak.The largest misalignment is in the Western Transverse Ranges and exceeds 30 • • (SAF is NW-SE, SKS WSW-ENE), with a stangion where both data sets provide constraints.The "area of influence" of each dataset is seen in the circular regions of Fig.3band chosen based on a rough estimate of averaging width, e.g.Fresnel zone of SKS.