Intra-oceanic island arcs are believed to have significantly contributed to the volume of continental crust since the Archean (Taylor and McLennan, 1985; Rudnick, 1995). The bulk composition of the continental crust is andesitic, whereas the island arc crust is basaltic, which is a well-known paradox. Therefore, if the continental crust is indeed formed by arc magmatism, significant modification of the island arc crust is required to transform the basaltic arc crust into the andesitic continental crust (Rudnick, 1995; Behn and Kelemen, 2006; Jagoutz and Schmidt, 2012; Jagoutz and Behn, 2013; Jagoutz and Kelemen, 2015). The Kohistan-Ladakh arc (KLA), sandwiched between the Karakoram in the north and the western Himalayan syntaxis in the south (Figure 1), is recognized as a Mesozoic intra-oceanic island arc, preserving a complete arc section (Tahirkheli, 1979; Bard, 1983; Burg, 2011). Thus, knowledge of the nature of the crust beneath the KLA and its correlation with other modern island arcs and orogenic belts is critical to understanding the process of the formation of a juvenile continental crust from island arcs (Miller and Christensen, 1994; Jagoutz and Kelemen, 2015).

The crustal structure beneath the KLA is poorly investigated as compared to the rest of the Himalayas and Tibet, probably due to its inaccessible and rugged terrain. Over the last few decades, several attempts have been made to understand the gross composition and seismic velocity structure, primarily using thermodynamic modelling and laboratory measurements of exposed rocks from the Kohistan arc (Chroston and Simmons, 1989; Miller and Christensen, 1994; Kono et al., 2009; Almqvist et al., 2013; Jagoutz and Behn, 2013). These studies reveal that the Kohistan arc represents a typical island arc crust characterized by a thin upper/mid crust overlying a very thick mafic high-velocity layer (HVL) with Vp >7 km/s beyond 40 km depth. This unusually thick mafic lower crust is interpreted as a density-unstable component of an arc crust that undergoes episodic modifications, probably via foundering, to create a typical continental crust of andesitic composition (Jagoutz and Behn, 2013; Jagoutz and Kelemen, 2015). However, the evidence of the arc-like preserved crust beneath the KLA has not yet been established using seismological observations, particularly because of the non-availability of broadband seismic stations, especially in Kohistan. The majority of the geophysical studies are confined to the eastern end of the Ladakh arc, which provides evidence of crust with low seismic velocities and partial melts (Rai et al., 2006, 2009; Caldwell et al., 2009). These results indicate strong crustal variability beneath the KLA that requires a high-resolution imaging of the region.

In this study, I compute a shear wave velocity (Vs) model of the crust beneath the KLA with a lateral resolution of ∼30 km using a trans-dimensional Bayesian analysis (Bodin et al., 2012a; Hawkins and Sambridge, 2015) of a combination of ambient noise and earthquake data. Using the model, this study provides evidence of an arc-like preserved crust with a thick, high-velocity lower crust beneath the Kohistan and a modified arc crust with comparatively low velocities beneath the Ladakh arc that resembles the Karakoram and southern Tibet.

The Kohistan-Ladakh terrain represents a fossil Cretaceous-Tertiary island arc complex that presumably formed over a north-dipping intra-oceanic subduction zone near the equator and, subsequently, was obducted between the colliding Indian and Eurasian plates at about 50 Ma (Tahirkheli, 1979; Bard, 1983; Jagoutz and Schmidt, 2012). The KLA is bounded by two suture zones (Figure 1). In the south, the Indus Suture Zone (ISZ), also referred to locally as the Main Mantle Thrust (MMT), separates the KLA from the western Himalaya. The Shyok Suture Zone (SSZ), or Main Karakoram Thrust (MKT), bounds the northern margin of the KLA and separates the Karakoram batholith in the north. The presence of these two suture zones, bounding the KLA, is interpreted as a result of two, almost parallel, north-dipping subduction zones (Supplementary Figure S1), one forming the KLA and another forming the Karakoram batholith, during the Early Cretaceous (Burg, 2011; Bouilhol et al., 2013; Jagoutz et al., 2015). The Nanga Parbat, sandwiched between the Kohistan arc in the west and the Ladakh arc in the east, is a fault-bounded, antiformal massif with the youngest exposed migmatites and granites (as young as 0.7 Ma) and the fastest exhumation rates up to ∼13 mm/year (Zeitler et al., 1993; Crowley et al., 2009; Govin et al., 2020).

The Kohistan arc is broadly divided, from south to north, into three segments (Jagoutz and Schmidt, 2012). These include: 1) Southern Plutonic Complex: southernmost part of Kohistan, which exposes mafic to ultramafic arc lower crust formed during 120–85 Ma (Figure 1), 2) Chilas Complex: mafic to ultramafic intrusions representing mid-crustal arc rocks created by the arc rifting at around 85 Ma, and 3) Gilgit Complex: dominantly calc-alkaline batholith with inclusions of volcanic and sedimentary rocks, which exposes mid/upper arc crust and records a magmatic history from 120–50 Ma. Despite similarities with the Kohistan arc, the Ladakh batholith is characterized by a volumetrically large volcanic and sedimentary unit with an increasing continental signature, which shows differences in tectonic processes eastward from the Kohistan (Rolland et al., 2000; Burg, 2011). The northeastern end of the Ladakh arc is bounded by the right-lateral Karakoram Fault (KKF).

This study is built upon the previous work of Kumar et al. (2022) by reanalysing the data with a focus on the KLA. The dataset comes from ∼530 broadband seismic stations located in the Pamir-Hindu Kush, western Tibet, western Himalaya, and adjoining areas. Fundamental mode Rayleigh wave group velocity measurements from period 5–60 s are obtained using analysis of ambient noise cross-correlations (Bensen et al., 2007; Schimmel et al., 2011) and earthquake waveforms recorded over the seismic stations (Supplementary Figure S2). The spatial resolution of group velocity dispersion is estimated by synthetic checkerboard tests (e.g., Rawlinson and Spakman, 2016). For this, input models with alternating positive and negative velocity anomalies of three different sizes (0.2°× 0.2°, 0.3°×0.3°, and 0.5°× 0.5°) are used. The recovered checkerboard models are presented in Supplementary Figures S3–S5, which show that the input models with sizes 0.3° × 0.3° and 0.5° × 0.5° are reasonably recovered, while the 0.2°× 0.2°-sized model suffers from lateral smearing in the region encompassing the KLA and adjoining area.

Assuming an optimum lateral resolution of 0.3°× 0.3° (∼30 km), as discussed above, the study area is divided at a grid interval of 0.25°× 0.25° and a two-step inversion scheme is adopted to compute the Vs model. First, a 2-D tomographic inversion is performed using the Bayesian trans-dimensional tree approach of Hawkins and Sambridge (2015), which provides group velocity maps in the period range of 5–60 s (Supplementary Figure S6). Subsequently, the group velocity dispersion at individual grid-nodes is inverted for the Vs-depth relation using a 1-D trans-dimensional Bayesian inversion of Bodin et al. (2012b). Finally, a 3-D shear wave velocity model is prepared by interpolating individual 1-D models. Detailed information about the data analysis and inversion scheme used in this study is described in Kumar et al. (2022), and a summary is provided in Supplementary Text S1.

Supplementary Figure S6 presents the group velocity maps at some representative periods (5–60 s) in terms of a perturbation from the regional mean. At shorter periods (5–20 s), which represent the upper crust, the group velocities within the KLA remain higher compared to sedimentary basins, e.g., Peshawar Basin, Shivalik sediments in the foothills of the western Himalaya, and the Tarim Basin, situated north of the KLA. At the intermediate period (30 s), sensitive to the middle crust, low-velocity zones (LVZs) appear beneath the Nanga Parbat and eastern Ladakh, which may indicate the presence of partial melts and/or aqueous fluids (e.g., Caldwell et al., 2009). For longer periods, i.e., more than 40 s, LVZs are observed primarily beneath the Nanga Parbat, the Karakoram, and eastern Ladakh. In contrast, the Kohistan arc, particularly its southern part, and western Ladakh (west of 76°E) show higher group velocities. The LVZs at these longer periods are possibly indicative of large crustal thickness and/or deep crustal melts beneath eastern Ladakh, i.e., the Baltoro batholith (Caldwell et al., 2009). To constrain the depth variation of shear wave velocities, the dispersion data at multiple periods is inverted for the Vs-depth relation and discussed in the following sections.

The seismic structure of island arcs provides along-strike variability in the arc crust, which provides evidence of arc modifications leading to the formation of modern continental crust [reviewed in Rudnick (1995) and Calvert (2011)]. Modern island arcs are typically characterized by a high-velocity lower crust, indicating the predominant mafic nature of the arc crust as compared to a more felsic continental crust (Calvert, 2011). Chroston and Simmons (1989) and Miller and Christensen (1994) reconstructed the seismic structure of the Kohistan arc before its obduction onto the Indian plate using geothermobarometry and thermodynamic modelling of exposed rock samples. Their studies demonstrated that Kohistan is representative of a typical arc-like crust with a thin upper and middle crust overlying a thick HVL (Vp ∼ 6.8–7.8 km/s) in the lower crust. Beyond 40 km depth, these studies predict strong inflection in Vp (>7.0 km/s), Vs (>4.0 km/s), and density (>3.3 g/cc). The crustal thickness was estimated to be around 50 km. Jagoutz and Behn (2013), using a large sample of exposed rocks from Kohistan, reconstructed a seismic structure showing a negatively buoyant lower crust at depths exceeding 40 km. Further evidence of the densified arc crust beneath the Kohistan comes from gravity measurements by Verma and Prasad (1987) showing a locally high Bouguer anomaly (Supplementary Figure S17B) in southern Kohistan, which was interpreted as the result of high-density shallow crustal intrusions from the deep crust/upper mantle. Based on these studies, it is evident that the Kohistan arc represents a well-preserved paleo-arc crust. A surface wave study by Hanna and Weeraratne (2013) provided an average 1-D velocity of Nanga Parbat and the KLA, but could not image the lateral variability in the interior of the KLA because of low ray coverage away from the Nanga Parbat Massif. On the other hand, several studies in the Ladakh arc show evidence of the contrasting structure of the crust, dominated by low seismic velocities, possibly due to melting (Caldwell et al., 2009; Rai et al., 2009).

The crustal velocity model presented in this study, which has a lateral resolution of ∼30 km, provides a unique opportunity to investigate the crustal variability within the KLA. Based on the contrasting velocity structure beneath the KLA, discussed in the earlier section, the area can be divided into major blocks, and their average velocities can be compared to previously computed models. In Figure 3, average 1-D Vs models are provided for different regions of the KLA and compared with the laboratory-derived model of Jagoutz and Behn (2013). The Kohistan arc is divided into two blocks: 1) southern Kohistan (Chilas + Southern Plutonic Complex), which exposes the mafic/ultra-mafic lower crustal rocks (blue dots in Figure 3A), and 2) northern Kohistan, which includes the Gilgit complex, dominated by the calk-alkaline batholith and sediments (green dots in Figure 3A). The Ladakh batholith is presented in red dots. These colored dots in Figure 3A represent grid-nodes of the velocity model within the KLA. The density-depth relation of the Kohistan arc, before its obduction on the Indian plate (Jagoutz and Behn, 2013), is shown in Figure 3B, which shows the density-unstable lower crust at depths >40 km. The 1-D average Vs beneath the Kohistan arc (Figure 3C) shows a uniform upper/middle crust of Vs ∼3.5–3.6 km/s, and a thick HVL (>20 km) in the lower crust (Vs > 4.0 km/s) at depths >40 km in the southern Kohistan, consistent with the density-unstable lower crust. The lower crustal velocities in northern Kohistan are significantly lower as compared to southern Kohistan, indicating northward thinning and/or the absence of the anomalous high-density lower crust. Note that the observation of the anomalously thick densified arc lower crust beneath southern Kohistan matches well with the locally observed high Bouguer anomaly (Verma and Prasad, 1987), as discussed earlier. Additionally, the distribution of exposed rock properties on the surface of the Kohistan arc shows a similar pattern. In southern Kohistan, mafic/ultramafic arc lower crustal rocks are exposed, whereas the upper crustal rocks, e.g., calc-alkaline batholith, are exposed in northern Kohistan (Figure 1). This correlation between the deep structure and the exposed rock properties in the Kohistan arc possibly indicates the southward thrusting of the arc’s lower crustal rocks over a north-dipping Indian crust during the India-Asia collision, which resulted in the obduction of mafic/ultramafic crust in the southern Kohistan. Using the Vp/Vs range of 1.70–1.83 for the Kohistan arc (Miller and Christensen, 1994) and the P-wave model of Jagoutz and Behn (2013), the possible range of shear wave velocities for the pre-collisional arc crust is shown as red horizontal bars in Figure 3C. This study’s average crustal velocity values in southern Kohistan are remarkably similar to those obtained by laboratory measurements. At depths between 20 and 30 km, the model of Jagoutz and Behn (2013) shows higher velocities, possibly due to the presence of an ultramafic horizon within a layered mafic intrusion of the order of 1 km thickness, which is difficult to resolve in seismic experiments (see Miller and Christensen, 1994). These results indicate that a typical arc-like crust is preserved in the Kohistan arc.

The average velocity of the Ladakh arc is presented in Figure 3D. The mid-crustal (at 20–40 km depth) velocities are reduced by 10%–15% as compared to the Kohistan arc, possibly indicating the effect of the Baltoro granite melting (Caldwell et al., 2009). In the Ladakh arc, the HVL (Vs > 4 km/s above Moho) is relatively thin (∼10 km) beyond 50 km depth. The velocity structure in the Ladakh arc is comparable to that in southern Tibet (Yang et al., 2012). The thin HVL above the Moho beneath the Ladakh arc represents the underthrusting Indian lower crust that reaches up to varying latitudes in the Pamir and Tibet (Kumar et al., 2022). Thus, the unusually thick HVL beneath the Kohistan arc must include a density-unstable arc lower crustal root, as predicted by previous studies, overlying the underthrusting Indian lower crust.

In summary, the velocity structure beneath the Kohistan arc is comparable to the pre-collisional structure derived by laboratory measurements, which indicates that the Kohistan arc has a well-preserved arc-like crust. The anomalously dense lower crust in southern Kohistan possibly represents the density-unstable arc root, typically found in modern arcs. With time, this densified root becomes heavier than the underlying upper mantle and may undergo episodic foundering, ultimately resulting in an andesitic continental crust (Jagoutz and Behn, 2013; Jagoutz and Kelemen, 2015). In contrast, the Ladakh arc, which behaves more like southern Tibet, has undergone significant modifications after the India-Asia collision. This contrasting crustal structure of the KLA is presented in a schematic diagram in Figure 4, which highlights the variability in the nature of tectonic processes in the KLA, as envisioned by Burg (2006).

A 3-D shear wave velocity model of crust, with a lateral resolution of ∼30 km, beneath the Kohistan-Ladakh arc is derived using ambient noise and earthquake waveform analysis based on Bayesian trans-dimensional inversion. The results are summarized in Figures 3, 4 and the major findings are listed below:

1) A uniform upper and middle crust of Vs ∼3.5–3.6 km/s and an anomalous high-velocity lower crust (Vs > 4 km/s) at depths exceeding 40 km are observed beneath the Kohistan arc.