than normal azimuthal anisotropy. Additionally, whether such anisotropy is related to ancient lithospheric deformation or reflects current plate motion remain enigmatic. To address these questions, a number of additional efforts are required, among which the most important ones are to establish a uniform set of shear wave splitting measurements for the entire continent and to estimate the depth of the source(s) of the observed anisotropy. The most recent continental scale shear wave splitting (SWS) study (Heintz & Kennett, 2005) was conducted nearly two decades ago, since then a large amount of data has been recorded by existing and newly established stations. Another critical step that has not been taken by previous SWS studies is providing some constraints on the depth of the source area of the observed anisotropy, as the low vertical resolution is a well-known weakness of the SWS technique resulted from the steep incidence of the teleseismic shear waves.

In this study, we first produce a uniform database of individual (rather than station-averaged) teleseismic shear wave splitting measurements for the Australian continent using all the publicly available data. We then use the individual splitting parameters to estimate the depth of the source of the observed anisotropy (K. H. Liu & Gao, 2011), a practice that has not been applied to the Australian continent. Finally, by combining our measurements with previous surface wave inversion results (Debayle & Kennett, 2000;Fishwick et al., 2008), we interpret most of the observations in terms of the relative strength of two anisotropic layers, which have approximately orthogonal fast orientations and reside in the lithosphere and the transition zone between the lithosphere and asthenosphere, respectively.

The Australian continent is composed of Precambrian cratons in the west and Phanerozoic orogen in the east (Figure 1). The three main cratonic units in the west and central mainland include the West, North, and South Australian cratons. The numerous fragments were assembled into these three cratons by ∼1830 Ma (Myers et al., 1996). The two Archean cratons, Pilbara and Yilgarn, collided along the Capricorn Orogen and formed the West Australian Craton during the early Proterozoic (Myers et al., 1996). The proto-Gawler and proto-Curnamona cratons were incorporated together along the Kimban Orogen to form the South Australian Craton between ∼1845 and 1700 Ma (Myers et al., 1996). The North Australian Craton was constructed by accretion of numerous crustal blocks, which include the Kimberly, Pine Creeks, Lucas, and Altjawarra cratons, with sutures marked by several orogens between ∼1950 and 1700 Ma (e.g., King Leopold, Halls Creek, and Tennant Creek) (Betts et al., 2002;Myers et al., 1996). These three major cratons collided at ∼1300 Ma to form a part of the Rodinia supercontinent (Myers et al., 1996). At ∼830 and 750 Ma, the Centralian Superbasin developed at the junction of the three cratons, following the breakup of Rodinia (Myers et al., 1996;Walter et al., 1995). After the Rodinia breakup, the Centralian Superbasin was broken into several smaller basins and intracratonic basins were formed at ∼750-540 Ma (Myers et al., 1996).

Several orogenic belts in the eastern Australian continent, including the Thompson, Lachlan, New England, and Delamerian orogens, developed along the east margin and formed the Tasmanides, which has a NE-SW dominant structure trend (Bello et al., 2019;Betts et al., 2002). The boundary between the Precambrian cratons and the Phanerozoic orogenic belts is referred to as the Tasman Line (Direen & Crawford, 2003). The origin of the intraplate volcanism in the eastern Australia, the Cosgrove track, is still debated (Davies et al., 2015;Davies & Rawlinson, 2014;Montelli et al., 2006;Rawlinson et al., 2017). Most previous tomography studies show the presence of positive seismic velocity anomalies in the depth range from 50 to 300 km beneath the area west of 140°E (Fichtner et al., 2009;Fishwick et al., 2005;Fishwick & Reading, 2008;Kennett et al., 2013;Rawlinson & Fishwick, 2012). The depth of the lithosphere-asthenosphere boundary (LAB) is greater than 200 km beneath the cratons and reduces to 100 km beneath the eastern coast (Conrad & Lithgow-Bertelloni, 2006;Pasyanos et al., 2014). Beneath the northeastern portion of the Australian continent, the high seismic wave speed anomaly and the undulation of mantle transition zone discontinuities are interpreted to be associated with the remnant of the Pacific plate subducted beneath the Melanesian Arc between 45 and 25 Ma (Ba et al., 2020;Hall & Spakman, 2002;Simmons et al., 2012).

One of the most common-used seismic techniques to study seismic azimuthal anisotropy is SWS analysis (Hess, 1964;Silver & Chan, 1991). A polarized shear wave traveling through an anisotropic medium splits into two quasi-shear waves with orthogonal polarization orientations and propagate with different speeds. The refracted phases at the core-mantle boundary (PKS, SKKS, and SKS, collectively XKS) are commonly used to determine the anisotropy of the crust and upper mantle beneath the stations (e.g., S. Gao et al., 1994;Heintz & Kennett, 2005;Silver & Chan, 1991). Two XKS splitting parameters, including the polarization orientation of the fast wave (ϕ or fast orientation) and the travel time difference between the fast and slow waves (δt or splitting time), are utilized to quantify the orientation and strength of the seismic azimuthal anisotropy, respectively. The resulting measurements of the XKS technique have a good lateral resolution and poor vertical resolution (Silver & Chan, 1991). Under some conditions including the presence of spatially varying simple anisotropy, the optimal depth of seismic anisotropy can be estimated using the spatial coherency of the observed splitting parameters (S. S. Gao & Liu, 2012;K. H. Liu & Gao, 2011).

Previous laboratory experiments and field investigations suggest that the main cause of seismic azimuthal anisotropy in the mantle is the lattice preferred orientation (LPO) of anisotropic mineral, mainly olivine in the upper (Betts et al., 2002) and previous shear wave splitting measurements. The straight lines represent previous shear wave splitting measurements (light green for Vinnik et al. (1992); blue for Barruol & Hoffmann (1999); black for Heintz & Kennett (2005); red for Bello et al. (2019); orange for Chen et al. (2021); and dark green for Eakin et al. (2021)). The circles represent seismic stations with data available at Incorporated Research Institutions for Seismology Data Management Center and AusPass with non-null measurements and green circles represent stations from which no measurements observed by previous studies. The white and gray arrows show the directions of the absolute plate motion based on the NNR-MORVEL56 (Argus et al., 2011) and HS3-NUVEL-1A (Gripp & Gordon, 2002), respectively. The thick purple solid lines represent craton boundaries and the dashed black line represents the Tasman Line (Betts et al., 2002;Direen & Crawford, 2003). mantle (Silver, 1996;Zhang & Karato, 1995). In most areas, more than one mantle deformation process can lead to LPO. For the Australian continent, the most likely processes include (a) simple shear originated from the relative movement between the lithosphere and the underlain asthenosphere which causes the fast orientation to be parallel to the APM direction, (b) horizontal lithospheric compression, which leads to the LPO to be normal to the direction of the maximum shortening, and (c) regional mantle flow system which results in flow-parallel fast orientations.

Most of the previous SWS studies in the Australian continent (Figure 1) show a complex pattern of the fast orientations (Barruol & Hoffmann, 1999;Bello et al., 2019;Chen et al., 2021;Eakin et al., 2021;Heintz & Kennett, 2005, 2006;Vinnik et al., 1992). Heintz and Kennett (2005) measured SWS at 190 stations covering the whole continent. The difference between the APM direction calculated based on the NUVEL1A-HS3 model (Gripp & Gordon, 2002) and the fast orientations observed at most of the stations was attributed to seismic anisotropy induced by the asthenospheric flow system around the lithospheric roots. Measurements in the Phanerozoic eastern Australia are consistent with the dominant strike of the orogenic belts, indicating lithospheric anisotropy. Using data recorded at 20 stations deployed on both sides of the Tasman Line, Heintz and Kennett (2006) reported the absence of observable azimuthal anisotropy. They proposed that the lack of anisotropy is caused by a two-layered anisotropy system with fast orientations that are perpendicular to each other. The two layered anisotropic structure with nearly orthogonal fast orientations is consistent with previous surface wave inversion studies, which proposed that at 100 km depth, the fast orientation is almost east-west, and at 200 km depth, it changes to north-south (Debayle, 1999;Debayle et al., 2005;Debayle & Kennett, 2000;Fishwick et al., 2008;Fishwick & Reading, 2008). Recently, Bello et al. (2019) obtained shear wave splitting measurements at 32 stations in the southeastern part of the Australian continent and the Bass Strait. The observed large splitting times and complex fast orientations were interpreted as the combined effects of frozen-in lithospheric anisotropy and present-day asthenospheric flow around lithospheric roots. Chen et al. (2021) reported shear wave splitting measurements at 4 stations beneath the Yilgarn Craton in west Australia and the fast orientations of the station averaged measurements are almost orthogonal to the APM direction calculated based on the NUVEL1A-HS3 model (Gripp & Gordon, 2002). They attributed the observations to the complex mantle flow resulted from LAB topography variation. Eakin et al. (2021) presented shear wave splitting measurements obtained with SKS and PKS phases from 25 stations in the BILBY array in central Australia and concluded that seismic anisotropy is mostly from lithospheric fabric as opposed to asthenospheric flow.

Seismic data recorded at 289 broadband stations are used in this study with the recording period from 1980 to February 2018. Among the stations, 156 are from the Incorporated Research Institutions for Seismology (IRIS) Data Management Center (DMC) and 133 are from the Australian Passive Seismic Server (AusPass) (Figure 1). The detailed data processing procedure to obtain the splitting measurements is described in K. H. Liu and Gao (2013), which is based on the minimization of the transverse energy technique (Silver & Chan, 1991). To balance the data quality and quantity, we utilized earthquakes with a cutoff magnitude of 5.6 for events shallower than 100 km and 5.5 for the deeper events. To ensure that all the usable XKS signals are included, the epicentral distance range for data requesting is set as 120°-180°, 95°-180°, and 84°-180° for PKS, SKKS, and SKS, respectively (Figure 2). Seismograms were windowed 5 s before and 20 s after the predicted XKS arrival times calculated based on the IASP91 model and were bandpass filtered using a four-pole two-pass filter with corner frequencies of 0.04 and 0.5 Hz.

The individual splitting measurements were manually verified and the data processing parameters, including the starting and end times of the XKS window and the frequencies of the band-pass filtering were adjusted if necessary. The ranking of each measurement (A: excellent, B: good, C: unusable, and N: null) was determined by visually checking and considering the following factors: (a) high S/N XKS arrivals on both the original radial and transverse components and significantly reduced energy on the corrected transverse component, (b) the goodness of waveform fitting and similarity between the fast and slow components after correction, (c) the presence of an elliptical particle motion pattern on the original and the linearity of the corrected particle motions, and (d) the significance and uniqueness of the minimum on the contour map of the corrected transverse component (K. H. Liu et al., 2008). Station averaged splitting times were computed as simple means, and station averaged fast orientations were calculated using directional statistics (Mardia & Jupp, 2000) and are termed as circular means.

BA ET AL.

10.1029/2022JB025866 5 of 17

A total of 511 well-defined (Quality A or B) individual measurements were observed at 116 stations (Table S1), including 436 from 74 IRIS stations and 75 from 42 AusPass stations (Figure 3). The station averaged measurements are shown in Figure 4 and Table S2, and the difference between the fast orientations and the APM direction (Gripp & Gordon, 2002) is shown in Figure 5. The circular mean of the fast orientations over the 511 individual measurements in the whole area is 48.5 ± 49.8°, and the simple mean of the splitting times is 0.85 ± 0.33 s, which is smaller than the global average of 1.0 for continents (Silver, 1996). The standard deviations (at 95% confidence level) for the fast orientation measurements range from 1.0° to 22.5° with a mean value of 9.8°, and those for the splitting time measurements range from 0.03 to 2.38 s with a mean value of 0.24 s. In Figure 6, the mean splitting times in overlapping 1° × 1° block with a moving step of 0.1° are presented.

Besides Quality A and B measurements, 175 null measurements, which are characterized by strong XKS arrivals on the original radial but no observable XKS arrivals on the original transverse components, were observed at 81 stations in the study area, including 25 stations without well-defined measurements (Figure 7). Most of these null measurements have an event back azimuth that is parallel or perpendicular to the fast orientation of SWS measurements and thus do not suggest the absence of significant azimuthal anisotropy (K. H. Liu & Gao, 2013;Silver & Chan, 1991;Yu et al., 2018). Based on the tectonic setting and the characteristics of the observed splitting parameters, we divide the study area into five regions (Figure 3) to facilitate discussion. Such a division does not imply that the observed anisotropy is entirely dependent on lithospheric domains. As a matter of fact, as demonstrated below, the shear wave splitting measurement pattern for a given region is not always clearly defined by the tectonics boundaries, indicating that lithospheric fabric is not the only origin for the observed seismic anisotropy. The mean fast orientation of each region is shown by the rose diagram in Figure 4, and the average splitting times for each region can be found in Table S3. Information about the original and corrected radial and transverse components, the particle motion patterns, as well as the contour plot of remaining energy on the corrected transverse component for each of the 511 Quality A and B and 175 null measurements can be found in the Data Availability Statement section.

Region A, which includes the Archean-Proterozoic West Australian Craton, has a lithospheric thickness of about 200 km which is among the thickest among all the five regions (Figure 4; Pasyanos et al., 2014;Pirajno & Bagas, 2008). Only 17 pairs of station-averaged shear wave splitting parameters were reported by previous studies in this area with different splitting times (Figure 1; Chen et al., 2021;Heintz & Kennett, 2005;Vinnik et al., 1992). In this study, we observed a total of 91 individual measurements from 30 stations in this region. The mean splitting time is 0.86 ± 0.35 s and the circular mean of the fast orientations is 67.8° ± 30.9° (Figure 3). Based on the lithospheric thickness and observed fast orientations, Region A was further divided into two sub-regions. In the area with thickened lithosphere north of 30°S (Region A1) (Figure 4), the observed fast orientations are mostly NE-SW, which is similar to the direction of the APM under the no-net-rotation frame (Argus et al., 2011). However, in the southern part of the craton (Region A2), the fast orientations are mostly E-W which is different from the APM (Argus et al., 2011), and are not consistent with measurements of the previous studies (Heintz & Kennett, 2005;Vinnik et al., 1992) but are consistent with the recent study of Chen et al. (2021). Azimuthal variations in the measurements observed at stations CP11_7J and BLDU_AU may suggest the existence of complex anisotropy in this region.

Region B is the North Australian Craton, including the Proterozoic orogen and basins (Pirajno & Bagas, 2008). A total of 175 well-defined individual measurements were obtained at 19 stations located in this region. The mean fast orientation is 9.75° ± 38.64°, which is parallel to the APM direction in the fixed hotspot frame (Gripp & Gordon, 2002) and consistent with previous SWS studies (Figure 1; Eakin et al., 2021;Heintz & Kennett, 2005). Previous surface wave tomography studies also observed mostly N-S fast orientations at ∼200 km depth (Fishwick & Reading, 2008;Kennett et al., 2004;Simons et al., 2002). In the central part of this region, a few E-W fast orientations are also observed in this study, which are parallel to the south boundary of the North Australian Craton. The mean splitting time is 0.80 ± 0.29 s which is the smallest among the five regions. Five PKS measurements observed at Station FITZ_AU show NW-SE fast orientations, and nine SKKS and SKS measurements show E-W fast orientations, suggesting that a two-layered anisotropy structure could exist in the western part of the North Australian Craton.

Region C is the South Australian Craton. There are 82 high quality (Quality A or B) measurements observed at 18 stations. In comparison, only six station-averaged measurements were reported by previous shear wave splitting studies in this region (Figure 1; Eakin et al., 2021;Heintz & Kennett, 2005), among which two stations have individual measurements. The fast orientations in this region are mostly E-W with a mean value of 68.47° ± 19.39°, which are consistent with results from previous tomography anisotropy studies in the lower upper mantle (e.g., Fishwick & Reading, 2008;Kennett et al., 2004;Simons et al., 2002) and are approximately perpendicular to the APM direction in the fixed hotspot frame (Gripp & Gordon, 2002). The mean splitting time in this region is 0.95 ± 0.27 s, which is comparable to the averaged value of 1.0 s for global continents (Silver, 1996). The null measurements observed in the southeastern portion of this region are most from event with a back azimuth that is perpendicular or parallel to the fast orientations. They provide independent verification of the Quality A or B measurements and are not suggestive of the presence of mantle isotropy in this region (Figure 7).

Region D includes the area between the three cratonic units. The eight stations located in Region D resulted in 40 measurements, with seven stations situated between the North and the South Australian cratons and one station (WRKA_AU) situated between the West and North Australian cratons. The circular mean of the fast orientations for Region D is 68.37° ± 9.74°. The seven stations between the North and South Australian cratons show E-W station averaged fast orientations, which is mostly parallel to the strike of the boundaries between the two cratons and is consistent with the results of Eakin et al. (2021), while station WRKA_AU has a N-S station averaged fast direction (Figure 4). The abrupt change of the fast orientation between Regions D and B can be observed on the boundary of the North Australian Craton. The splitting times in Region D are also larger than those in Region B, with the mean value of 0.94 ± 0.12 s.

Region E includes the area to the east of the Tasman Line, consisting of units formed in the Phanerozoic, and the lithosphere in this area is significantly thinner than western and central Australia (Betts et al., 2002;Direen (Gripp & Gordon, 2002). The grey arrows represent the APM direction. The triangles represent stations at which well-defined measurements were observed. & Crawford, 2003;Pasyanos et al., 2014). There are 123 well-defined individual measurements recorded at 41 stations. The circular mean of the fast orientations is 42.47° ± 56.61°, which is mostly parallel to the structural trend of the orogenic belts in this region. Our results in this region are in generally agreement with those from previous shear wave splitting studies in fast orientations but show smaller splitting times (Bello et al., 2019;Heintz & Kennett, 2005). The fast orientations are mostly NW-SE and E-W on the southeast coast and islands, which are different with the NE-SW fast orientations obtained in the rest of the sub-area (Figure 7). The mean splitting time of this region is 0.95 ± 0.27 s and is not consistent with those from a previous study which observed larger splitting times of 0.66-2.7 s at approximately the same stations from SKS and SKKS phase in the southeast part of Region E (Figure 1; Bello et al., 2019). Null measurements at several stations are observed, which cover a large azimuthal range and thus may suggest the absence of seismic anisotropy in the southeast coastal area (Figure 7).

In the following we attribute the observed splitting to the LPO of anisotropic mantle minerals especially that of olivine, although it has long been recognized that mathematically it is possible that anisotropy can be mimicked by isotropic heterogeneities (Backus, 1962;Fichtner et al., 2013). Specifically, it has been proposed that isotropic heterogeneities such as lithospheric dykes (e.g., S. Gao et al., 1997;S. S. Gao et al., 2008), preferably oriented cracks (Crampin & Chastin, 2003), and oriented melt pockets (e.g., Bastow et al., 2010) can lead to azimuthal anisotropy. However, several observations are not consistent with the notion that isotropic heterogeneities are the dominant cause of the observed XKS splitting in the study area which is composed of some of the oldest continental blocks on Earth. First, oriented cracks are mostly the main contributing mechanism for the brittle upper crust and therefore they can only contribute to a small fraction of the observed ∼1 s splitting at stable continental areas. Second, oriented melt pockets are commonly found in tectonically active and extended areas such as rifts and are rarely found in stable continental areas. Third, while lithospheric dykes cannot be completely ruled out, their contribution is small in stable cratonic areas where the cooling of the magmatic dykes led to great reduction in the velocity contrast between the dykes and the surrounding lithospheric material. Therefore, while an isotropic heterogeneity origin of the observed splitting cannot be completely ruled out, it is clear that the observed splitting is largely the result of LPO of anisotropic mantle minerals.

Splitting parameters recorded at densely spaced seismic stations in the central part (Region B and D) and southeast coast (southeastern part of Region E) of the Australian continent provide ideal circumstances to estimate the depth of the center of the anisotropy layer by applying the spatial coherency approach (S. S. Gao & Liu, 2012;K. H. Liu & Gao, 2011). A detailed description of the procedure together with a FORTRAN program can be found in S. S. Gao and Liu (2012). The approach is based on the Fresnel-zone principle and can be considered as an extension of the intersecting Fresnel-zone approach proposed by Alsina and Snieder (1995). It is based on the realization that when the individual splitting measurements are projected to the "real" depth, measurements from different event-station pairs reach the maximum spatial coherency. In practice, the procedure searches for the anisotropy depth by calculating the spatial variation factors, assuming that the depth ranges from 0 to 400 km for computing the location of ray piercing points, with a 5 km interval. The spatial variation factors are the weighted sums of the standard deviation of the splitting parameters at ray piercing points. The optimal depth corresponding to the minimum variation factor can be observed on the resulting curve, which is 220 and 100 km beneath the central and southeast Australian continent, respectively (Figure 8).

Beneath the central Australian continent (Regions B and D), compared with the lithospheric thickness (Figure 4) and previous tomographic studies (Conrad et al., 2006;Fichtner et al., 2009;Pasyanos et al., 2014), the estimated optimal depth of 220 km suggests that the observed anisotropy is mostly from the lithosphere-asthenosphere transition zone (Debayle & Ricard, 2013), as proposed elsewhere such as the eastern United States and the eastern Himalayan Syntaxis (L. Liu et al., 2019;Yang et al., 2017). This interpretation is different from that of a recent shear wave splitting study in the same area by Eakin et al. (2021), who did not report measurements with N-S orientations at three of the stations and consequently proposed a lithospheric origin for the mostly E-W oriented anisotropy. In the southeastern part of Region E, although the estimated resulting depth is shallower than that beneath central Australia, it is deeper than the lithosphere-asthenosphere boundary proposed by the LITHO1.0 model (Figure 4; Pasyanos et al., 2014), suggesting that the anisotropy beneath this area is from the upper asthenosphere. In other regions of the study area, the spatial coherency approach cannot be used mostly due to the large distance between the seismic stations and limited azimuthal coverage of the seismic events.

To provide additional constraints on the depth distribution of the observed seismic anisotropy, we next compare our splitting measurements with results from global and regional scale surface wave tomography studies (e.g., Debayle & Kennett, 2000;Fishwick et al., 2008;Fishwick & Reading, 2008;Kennett et al., 2004;Simons et al., 2002). Debayle and Kennett (2000) presented a 3-D model for the SV wave heterogeneity and seismic azimuthal anisotropy. They revealed a drastic change of the fast orientations between the upper 150 km and a layer below 150 km, and proposed that the anisotropy in the upper 150 km is related to deformation frozen in the lithosphere and that in the lower layer is due to the northward motion of the Australian Plate. Simons et al. (2002) proposed an azimuthal anisotropic shear wave speed model, suggesting a very complex pattern of anisotropy in the top 150 km and a smoother pattern at larger depths. Their measurements show that the fast orientations in the deeper layer are not corresponding to the plate motion direction everywhere. The measurements of Kennett et al. ( 2004) agree with those from a previous study (Debayle & Kennett, 2000) and show mostly E-W fast orientation at ∼100 km and N-S at 200 km depths (Figure 9). Two surface wave tomography studies (Fishwick et al., 2008;Fishwick & Reading, 2008) provided azimuthal anisotropy at the depths of 75, 150, 250, and 100, 150, 200 km, respectively. Their measurements were in general agreement with previous studies that suggest two layers of anisotropy with perpendicular fast orientations. Besides the regional surface wave tomography constraints on the azimuthal anisotropy beneath the Australian continent, several global seismic azimuthal anisotropy studies similarly indicated that the upper layer is dominated by E-W and the lower layer is mostly characterized by N-S fast orientations (Debayle et al., 2005;Debayle & Ricard, 2013;Long & Becker, 2010;Montagner, 2007;Schaeffer et al., 2016), although the boundary between the two layers varies among different studies.

Since the APM-parallel fast orientations in the lower layer are approximately perpendicular to the E-W fast orientations in the upper layer, the splitting times observed using the SWS technique reflect the splitting time difference between the two layers, and the observed fast orientation is the same as that of the layer with the greater splitting time (Silver & Savage, 1994). To provide a more quantitative comparison between surface wave tomography and SWS measurements, in Figure 9 we overlap our splitting measurements with results from surface tomography by Kennett et al. (2004) at the depth of 125 and 200 km. At the shallower depth of 125 km, the two sets of measurements both show mostly E-W fast orientations in the continental interior and are less consistent in the peripheral areas. In contrast, at the 200 km depth, the splitting results and surface wave anisotropy observations are generally consistent with each other in the peripheral areas with thin lithosphere. These relationships can be explained by the fact that in the continental interior, the thick lithosphere possesses E-W anisotropy with a splitting time that is greater than that of the N-S asthenospheric anisotropy and thus results in E-W fast orientations from SWS analysis; in contrast, in the peripheral areas, anisotropy from lithospheric fabric is weak and thus the observed anisotropy from the SWS technique mostly reflects anisotropy from mantle flow that is either related to APM or deflection by the topography of the continental keel, as discussed below.

It is well established that simple shear caused by the relative movement of the lithosphere and asthenosphere can induce seismic azimuthal anisotropy, with the fast orientation being parallel to the shear direction (Silver, 1996;Zhang & Karato, 1995). Many SWS observational studies attributed the observed APM parallel fast orientations to the simple shear in the transition layer from the lithosphere to asthenosphere associated with the plate motion, as being proposed in North America and South Africa (e.g., Fischer et al., 2010;K. H. Liu et al., 2014;Reed et al., 2017;Yang et al., 2017).

In a fixed hotspot frame, the Australian continent is the fastest moving continent among all the major continental blocks on Earth, presently at a rate of 83 mm/yr toward the north as determined by the HS3-NUVEL-1A model (Gripp & Gordon, 2002). The fast orientations observed in the North Australian Craton (Region B) are mostly N-S (Figure 4), which are parallel to the HS3-NUVEL-1A model. Such an approximate parallelism is also observed in Region E and the northern part of Region A. These observations suggest that simple shear in the rheologically transitional layer between the lithosphere and asthenosphere is a contributing anisotropy forming mechanism beneath three (northern, western, and eastern) of the four peripheral regions of the Australian continent.

Vertically coherent deformation of the lithosphere can lead to azimuthal anisotropy with fast orientations being orthogonal to the maximum horizontal shortening direction (Silver, 1996;Silver & Chan, 1991). The largely E-W fast orientations observed in Regions C and D suggest that anisotropy associated with N-S oriented lithospheric shortening has a greater strength than that produced by the APM, resulting in a net fast orientation that is identical to the fast orientation of the former. Surface wave tomography studies show that the fast orientation for both the shallow (<150 km) and deeper (>150 km) layers are N-S oriented beneath the North Australian Craton, and the deeper layer possesses much larger anisotropy (Figure 9). If we assume that the N-S oriented anisotropy with a splitting time of ∼1 s observed at the northern tip of the continent is purely APM oriented and is representative of APM induced anisotropy across the entire continent, we can conclude that the splitting time caused by (Kennett et al., 2004). The background color indicates the SV wave velocity at 125 and 200 km of the AuSREM model (Kennett et al., 2013), and the orientation and length of the red and black bars represent the fast orientation and strength of anisotropy from the tomography study and observed station-averaged SWS measurements, respectively. lithospheric fabric is about 2 s. For a lithospheric thickness of 200 km, this splitting time corresponds to an anisotropy of ∼4%, a value that is comparable from laboratory measurements (Silver, 1996;Silver & Chan, 1991).

The hypothesis that E-W oriented lithospheric fabric causes a stronger anisotropy than that related to APM in most areas of Regions C and D can also be independently confirmed by results of surface wave tomography studies (Figure 9), which suggest that the fast orientations for the shallower and deep layers are approximately normal to each other, and the degree of anisotropy in the shallower layer is greater than that in the lower layer in Regions C and D where E-W fast orientations are observed. Additionally, the degree of anisotropy revealed by surface wave tomography in the top 150 km in the two regions is about 4%, which is consistent with the estimated value from SWS measurements.

Besides lithospheric fabric, deep continental roots may modulate the mantle flow system and generate seismic anisotropy that is parallel to the edge of the continent, as observed along the southern and eastern edges of the North American craton (Fouch et al., 2000;Refayee et al., 2014;Yang et al., 2017). This mantle flow system may contribute to the observed anisotropy with edge parallel fast orientations observed along the western, southern, and eastern margins of the Australian continent. In southern West Australian Craton, Chen et al. ( 2021) also suggested that the observed APM orthogonal anisotropy can be attributed to the complex mantle flow due to the variation of the depth of the LAB. In particular, in the southern part of Region E (southeast coast and the Bass Strait), the large difference between APM direction (Gripp & Gordon, 2002) and fast orientations (Figure 5), as well as results from the spatial coherency analysis (Figure 8b), which indicate an upper asthenospheric origin of the observed anisotropy, suggests that the anisotropic layer may be caused by the mantle flow system flowing around the cratonic root, a mechanism that has been suggested by Bello et al. (2019) based on SWS measurements in the Bass Strait area.

The study has several global implications. First, when SWS measurements are interpreted, it is essential to consider contributions from both the lithosphere and asthenosphere, a task that is hindered by the poor vertical resolution of the splitting measurements. To partially overcome this limitation, jointly interpreting the splitting measurements with results from anisotropic seismic tomography, lithospheric thickness measurements, and depth estimation using the spatial coherency approach can be effective, as demonstrated in this study. Second, the conventional wisdom that (a) when the fast orientation is parallel to the APM, the observed anisotropy must be from the asthenosphere, and (b) when the fast orientation is parallel to old lithospheric fabric, the observed anisotropy must be related to past orogenic events might be proven to be wrong in many continental areas because both the lithosphere and asthenosphere can possess strong anisotropy. Third, the observed anisotropy could still be related to mantle flow even if the fast orientations do not agree with the APM. This is because the APM is model dependent, and the mantle flow field can be modulated by the shape of the bottom of the lithosphere, slab subduction, and global and regional density heterogeneities.

The complex spatial distribution of seismic azimuthal anisotropy beneath the Australian continent is revealed by systematic spatial variations of the observed splitting parameters using shear wave splitting analysis. The difference of the fast orientations and APM direction and the spatial coherency analysis suggest that the anisotropy beneath Australia is caused by several factors, including the APM induced mantle flow, the keel deflected flow system around the lithosphere roots, and lithospheric fabric resulted from horizontal compression. The South Australian Craton and the orogenic belt north of it possesses strong E-W oriented lithospheric anisotropy with a lithospheric contribution of 2 s to the splitting times, while APM induced N-S oriented anisotropy dominates in the North Australian Craton. The western, southern, and eastern margins of the continent are dominated by root deflected mantle flow that produces margin parallel fast orientations. The relative strength and approximate orthogonality of the APM induced anisotropy and frozen-in lithospheric fabric provide a viable explanation for the puzzling observation of pervasive existence of weak and spatially variable azimuthal anisotropy in this fast moving continent. network codes of 1P (https://doi.org/10.7914/SN/1P_2011), 6F (https://doi.org/10.7914/SN/6F_2008), 7B (https://doi.org/10.7914/SN/7B_1993), 7D (https://doi.org/10.7914/SN/7D_1997), 7G (https://doi.org/10.7914/ SN/7G_2000), 7J (https://doi.org/10.7914/SN/7J_2005), 7K (https://doi.org/10.7914/SN/7K_2007), AU (International Miscellaneous Stations), II (https://doi.org/10.7914/SN/II), IU(https://doi.org/10.7914/SN/IU), Sy (International Miscellaneous Stations) were used for access to waveforms. Related figures for all the 511 Quality A and B and 175 null measurements are publicly accessible at https://figshare.com/articles/dataset/ AuSWS_measurements_AB_N_zip/21830430.