Abstract

Cenozoic extension of the Qinling range-Weihe Graben system has occurred in response to the uplift and growth of the Tibetan Plateau. Rapid exhumation of the northern Qinling range since the late Miocene is also regarded as resulting from the eastward expansion of the northeast part of Tibet. Tectonic evidence of this in the landscape remains unclear, but the fluvial system can provide a sensitive proxy record of tectonic forcing through space and over time scales of 105–107 a. Here, we present a study of channel profiles in the northern Qinling range, which forms a footwall highland separated from the southern Weihe Graben by active normal faults. We identify a population of knickpoints that separate river profiles with a gentle upstream gradient from steeper downstream reaches. Above the knickpoints, steepness indices increase from the central part towards the west and east, whereas channel steepness shows its highest values in the Huaxian-Huayin section. We observed no systematic changes of channel steepness pattern as a function of rock resistance, drainage area, or channel concavity. Correlation analysis between channel steepness and basin elevation and relief documents the control of tectonic forcing on regional topography. While bearing no relation to geological outcrop boundaries, the knickpoints show a strong correlation between retreat distance, catchment area, and river length. We infer that the knickpoints formed in response to an increase in mountain uplift rates and retreated as a kinematic wave. Under linear slope exponent n, we calibrated channel erodibility K~1.00±0.44×106m0.1/a and derived knickpoint ages of 5.59±1.80Ma. Combining the ages of onset of active faulting and mountain growth in the NE Tibetan Plateau (8–10 Ma, e.g., Liupan Shan, Jishi Shan, and eastern segments of the Haiyuan and Kunlun faults) and in the southwest Qinling range (9–4 Ma), we conclude that growth of the NE Tibetan Plateau began in the mid-Miocene time and expanded eastwards to the Qinling range-Weihe Graben during the late Miocene and early Pliocene.

1. Introduction

Knowledge of both spatial and temporal changes in rock uplift rates places constraints on the geometry and tectonic activity of fault systems [1, 2], provides a diagnostic test for evaluating regional landscape evolution models [3, 4], and thus constitutes a foundation for understanding large-scale lithospheric deformation dynamics [5, 6]. Traditional leveling surveys, GPS, and InSAR play important roles in determining decadal-scale surface displacements but remain challenging for mostly capturing interseismic slip and requiring continuous measurements [79]. Studies on active tectonics (103–104 a) and low-temperature thermochronology (>106 a) provide details on long-term rates and patterns, which have always been focused on deformed geologic structures, well-exposed fault outcrops, and available field samples [1012].

For an absence of such strict conditions, some geomorphic features (e.g., topographic relief and slope-length index) have been calculated to estimate relative differences in the velocity field of rock uplift [1315]. Fluvial topography evolves as a competition between the rates of rock uplift and river incision, leading to the first-order control of rock uplift patterns on both the river drainage shape and longitudinal profile morphology in active orogens [5, 16]. A number of models have been proposed to interpret various river incision processes in an effort to set theoretical and practical frameworks for understanding the response of channels to variations in rock uplift rates [1719]. Instead of focusing on specific incision processes, the well-known stream-power incision model (SPIM) treats an overall incision rate as the function of channel bed shear stress, rock strength, and flood frequency [17]. In the simple form of SPIM, flood frequency and sediment flux only scales with the drainage area [5]. Under the assumption of a balance between rock uplift and river incision, the model suggests a powerful means of channel profile steepness to map the spatial patterns of differential rock uplift [2022].

However, in most active orogens where the equilibrium assumption cannot be met, knickpoints develop along channel profiles as a response to temporary changes in tectonic uplift rates [1, 4, 23]. Work on transient landscapes is a much more complicated story. Usually, three assumptions have been mostly used in the detachment-limited system: (1) average discharge and sediment flux scales with the drainage area; (2) the incision process can be well described by the SPIM along the whole channel, i.e., one can neglect the influence of oversteepened reaches; and (3) the drainage network has remained fixed through the time scale recorded by the fluvial channels. Accordingly, the upstream migration velocity of knickpoints can be assimilated to a kinematic wave [1, 2430].

In this contribution, we present such a case study in the northern Qinling range. The mountain range is a footwall upland bounded by a normal fault along the southern Weihe Graben, Central China (Figures 1 and 2). The landscape has been undergoing intense deformation, as evidenced by active channels that cut through the steep mountain margin and by the record of historical earthquakes (Ms5) along the fault. Field investigations and work on active tectonics were carried out on the Huaxian-Huayin fault segment in relation to the largest Huaxian earthquake (Ms~8.5, 1556 CE) and to its pronounced fault scarps and fault planes [3135]. In other parts of the area, in contrast, research is limited by the scarcity of good geological exposures and obvious geomorphic markers. Despite geomorphological analysis of the tectonic activity along the southern Weihe Graben [31, 36], the fluvial response of steep channels to differential rock uplift rates remains unclear. We perform analysis on river longitudinal profiles and knickpoints. To assess the ability of the SPIM to retrieve tectonic information, we discuss the influences of potential variables (e.g., lithology, sediment flux, catchment divide migration, transient signal migration, and spatially variant rock uplift) on the evolution of the river systems. We exploit the rules of kinematic wave speed of knickpoint migration to estimate the onset timing for rapid mountain growth. Additionally, we test to what extent the active differential rock uplift is responsible for the high mountain elevation and topographic relief in this landscape.

2. Regional Setting

The ~1000 km long Qinling range extends WNW–ESE and forms the middle part of the Central China orogenic belt (Figure 1). The mountains, formed by collisions between the South China and Yangtze Blocks in Paleozoic and Triassic times [37, 38], are mainly divided into three structural units [39]. The western Qinling range is constrained by the Shangdan suture to the north, Mianlue suture to the south, and Ningshan Fault to the east (Figure 1). Owing to the uplift and expansion of the Tibetan Plateau, the western Qinling has experienced strong compression deformation and evolved into parts of the northeastern plateau. Enkelmann et al. [40] suggested that the lower crustal flow has been diverted around the Sichuan Basin and beneath the SW Qinling since the late Miocene and thus caused active plateau uplift in this area, while the eastern Qinling range (east to the Ningshan Fault) is less affected by the plateau expansion and is mainly under the control of regional extension in the Cenozoic [39].

The northern Qinling range is located north to the Shangdan suture and bounded by the Northern Qinling Margin Fault (NQMF) and the Huashan Piedmont Fault (HPF) (Figure 2). Active normal faulting of the NQMF and HPF causes asymmetric topography, with moderate relief to the south of the mountain range but a remarkably steep topographic gradient to Weihe Graben. Taibai Mountain, with an elevation over 3700 m, is the highest peak of the range, while the elevation of the adjacent graben is only ~400 m (Figure 2). The northern Qinling range consists of the Archean crystalline basement, and most parts are overlain by Precambian metamorphic rocks and Mesozoic to Paleozoic sedimentary strata (Figure 3). Yanshanian granite plutons exposed in the range, e.g., the Taibai and Huashan granites, were emplaced ca. 120–160 Ma (zircon U-Pb ages). This indicates a total amount of 12–13 km of exhumation since the late Jurassic and early Cretaceous [41, 42]. Low-temperature thermochronology has revealed accelerated uplift of the mountain range during the Oligocene and since the late Miocene [4247]. Although it is difficult to link the Oligocene pulse of exhumation with the underlying geodynamic causes, e.g., far-field effect of the India-Asia collision, Pacific subduction, or opening of the Sea of Japan, the latest cooling stage is suggested to be the result of the growth and eastward expansion of the NE Tibetan Plateau [42, 43].

The uplift of the northern Qinling range and subsidence of the Weihe Graben constitute an intracontinental extension system. With ~7 km thick Cenozoic deposits, the tectonosedimentary evolution of the graben portrays polyphase exhumation of the range footwall [48, 49]. The industry borehole transects in the Lantian area and shows that the deposition began ca. 40–45 Ma [50, 51]. Deposits of late the Eocene age are about 1 km thick, resulting in an average depositional rate of <0.10 mm/a. The rapid exhumation of the northern Qinling range causes an increase in the deposition rate to be ~0.20 mm/a, which is derived from the ~2 km thick Oligocene deposits in the graben [33]. A large, late Oligocene to middle Miocene lacuna could indicate nearly no exhumation of the adjacent mountain range, which was also shown by both apatite fission track (AFT) and U-Th/He (AHe) studies in the hanging wall [42]. The final episode of rapid subsidence of the Weihe Graben began in the late Miocene time; this was inferred from accelerations in deposition rates in the basin as well as uplift rates in the northern Qinling range [43, 51, 52].

3. Methods

We extracted 55 river profiles along the northern Qinling range from the 1 arcsec SRTM DEM (digital elevation model). These rivers, which vary in length (~4–70 km) and catchment area (~2–1000 km2), cut through the steep topography and flow to the Weihe Graben.

3.1. River Profile Analysis

In active orogens, the fluvial system evolves as a competition between the rates of rock uplift (U) and river incision (E). When U=E (steady state), fluvial channels adjust local gradients (S) as a power-law scaling function of the contributing drainage area (A) [17]:
(1)dzdx=EK1/nAm/n,
where z is the elevation, x the distance, and K the channel erodibility, and m and n are constants. Equation (1) holds for A>Acr, a critical threshold for transition from debris flow to fluvial processes (typical value of 0.1–5 km2; [22, 53]). ks=U/K1/n and θ=m/n are the channel steepness and concavity indices, respectively. A linear fit to the log-transformed slope-area data can produce the two parameters. In the normalized equilibrium profile, θ is a metric to measure how concave the profile is [28] and can be an indicator of a channel bedrock substrate [20]. ks varies systematically with uplift/incision rates and thus can be used to map the spatial pattern of tectonic/climatic forcing [6, 54, 55].
Identifying variable concavities can be achieved with the help of slope-area analysis. But, a channel slope that is estimated via resampling and differentiating noisy topographic data usually contains considerable uncertainties. The integral approach (integral form to Equation (1); [56]) was proposed instead as an alternative for determining stream power parameters:
(2)z=EKA0m1/n·χ,(3)χ=0xA0Axm/ndx.
Setting the reference drainage area A0=1m2, the slope of the χz plot is ks. Nevertheless, along the whole channel, Equation (3) uses a constant θ to calculate χ. Thus, for channel profile analysis, we identified changes in channel concavities by slope-area plots, which were then used as regression limits in the χz domain (as recommended by [57]). For the dependence of ks on θ, a reference concavity (θref) is always used to calculate normalized channel steepness, ksn. Here, we chose θref=0.45 as used in most studies [3, 22, 54]. Thus, Equation (2) produces a power-law scaling function to estimate erodibility K:
(4)ksn=EK1/norK=Eksnn.

In addition, we also generated a χ-value map of drainage basins to determine catchment divide stability. Although cross-divide contrasts in mean gradient and local relief appear to consistently provide useful information at the present stage [58], a χ-anomaly can also indicate a possible future divide instability [59].

3.2. Catchment Concavity

In Equation (3), m/n is usually not known a priori. For a steady-state channel, concavity is what best fits the χz plot [56]. In a landscape under spatially and/or temporally variable conditions, however, channel concavity varies a lot. Especially below the knickpoints, channel segments change their gradients in response to various incision processes, causing large variations in channel concavities [24]. We thus followed Goren et al. [1] to calculate catchment concavity, which is the m/n value that minimizes scatter in all the χz curves of a drainage basin. We divided the χ domain into 100 bins and calculated the elevation scatter as
(5)Scatter=1NMi=1Nziz2,
where N is number of pixels for all the channels, M=100 is number of bins, zi is the elevation at each node, and z^ is average elevation of each bin.

Without assuming a linear relationship between χ and elevation, the colinearity between the trunk and all the tributaries of a drainage basin (Equation (4)) holds under three scenarios: (1) steady-state catchments with spatially uniform uplift patterns; (2) block uplift conditions, i.e., uplift rates variable through time but invariant in space; and (3) tectonic uplift rates following the dominant spatial distribution of χ space [1, 56, 60].

3.3. Slope-Break Knickpoint

As uplift differs from incision, we can observe knickpoints along river longitudinal profiles. According to the different forms observed on slope-area plots, knickpoints can be classified into two end-member morphologies: vertical-step and slope-break [61]. The former has abrupt increases only in local channel gradients, which is often anchored in space and relevant to small-scale heterogeneities, e.g., changes in bedrock resistance [62] and steep cascades caused by debris-flow-, earthquake-, or landslide-related deposits [5, 63]. By contrast, slope-break knickpoint, which represents a break in the power-law scaling function of the channel slope and drainage area, usually sets a mobile boundary for dividing channel segments under equilibrium with the previous (above) and present (below) tectonic uplift rates. Thus, it indicates a transient landscape, which is accompanied by significant variations in steepness indices up- and downstream [27, 28, 64]. Slope-break knickpoints are usually diagnostic of a passing incisional wave. Assuming the simple case of an abrupt increase in uplift rates, from Ui to Uf, the wave speed (horizontal velocity vH) can be written as [29]
(6)vH=Am/nK1/nUfUiUf1/nUi1/n.
Assuming that the channel segments below and above the knickpoint are under steady state with respect to present and previous conditions, respectively, we can derive
(7)Uf=Kksn_fn,(8)Ui=Kksn_in,
where ksn_f and ksn_i are downstream and upstream channel steepness, respectively. Thus, we can derive another form:
(9)vH=KAm/nksn_fnksn_inksn_fksn_i.
Accordingly, the response time (or knickpoint age), τ, for knickpoints propagating from the river outlet to a point along the channel profile is
(10)τ=0x1vHdx=ksn_fksn_iksn_fnksn_in0x1KAm/ndx=ksn_fksn_iksn_fnksn_in1KA0m/nχ.
Assuming a linear relationship between river incision and channel gradient (i.e., n=1), Equations (9) and (10) regress to a simple form:
(11)vH=KAm.(12)τ=1KA0mχ

Equations (11) and (12), a widely used linear form, show an area-dependent kinematic wave speed of the knickpoint retreat [65]. Equation (10) additionally indicates a tectonically independent recession rate and a constant χ value for knickpoints that formed at the same time [1, 24].

Under the assumptions (1) that the upstream reach of the knickpoint has been at equilibrium under a set of previous conditions and (2) that the drainage network is stationary without a long-wavelength tilt, one can reconstruct the “relict” channel profile using the χz scaling function (Equation (2)). Projecting the paleochannel to the river outlet, the difference (Δz) between the present-day elevation and the reconstructed elevation constrains the lower bound of the incision magnitude since the onset of the transient signal [28]. Dividing Δz by τ, we can derive the minimum rates of channel incision.

4. Results

We identified 24 rivers with smoothed concave-up profiles and 31 rivers with slope-break knickpoints that separate channels with different levels of steepness. Although the river shape may not be diagnostic of equilibrium conditions, we hypothesize nonetheless that these profiles document steady and transient states (stream-power parameters for these rivers are given in Table S1). The steady-state rivers (θ~0.42±0.17) and upstream channel segments (θ~0.41±0.20) have similar concavities (Figure 4(a)). They satisfy the normal distribution ~ N(0.42, 0.192), i.e., meanvalue=0.42 and variance=0.192. Below the knickpoints, θ varies widely from −5 to 4 (mean value 0.27±1.20) but mostly falls in a narrow range of 0.2–0.6. Furthermore, the concavities of these catchments vary from 0.2 to 0.7 and also meet N(0.42, 0.102) (Figures 4(b) and S1). Most channels thus share concavity values similar to steady-state rivers developed under uniform settings of lithology and climate [20, 22].

However, channel steepness indices vary systematically below and above the knickpoints (Figures 4(c) and 4(d)). The ksn indices of steady-state channels and downstream reaches vary like a sine wave, according to which the southern Weihe Graben could be divided into six segments. In Segment 1 (mean ksn119±20m0.9), from the west to Baoji, the steepness decreases from ~150 to 105 m0.9. In Segment 2 (150±33m0.9), i.e., Baoji to Zhouzhi, ksn increases to >180 m0.9 and then decreases to ~120 m0.9. In Segment 3 (150±25m0.9), i.e., Zhouzhi to Lantian, ksn increases to ~200 m0.9 and then decreases to ~120 m0.9. Along the NQMF, average ksn indices of each segment are nearly the same, whereas along the HPF, ksn varies a lot. Segment 4 (145±44m0.9), from Lantian to Huaxian, displays a sharp increase from ~120 to 250 m0.9 and then a sudden drop to 150 m0.9 followed by a gradual decrease to ~120 m0.9. In Segment 5 (189±47m0.9), i.e., Huaxian to Huayin, ksn reaches the highest value of 272 m0.9 (river No. 33) but decreases toward both sides. In Segment 6, from Huayin to Lingbao, ksn increases to 150 m0.9 and then decreases to ~50 m0.9. The average ksn value of steady-state streams (81±20m0.9) is lower than that of the downstream reaches (121±25m0.9), whereas above the knickpoints, ksn shows a different pattern. The ksn indices increase from the links (~35 m0.9) between NQMF and HPF towards both sides and reach the high values in Segments 2 (103±20m0.9) and 5 (100±24m0.9). Segments 3 (84±10m0.9) and 6 (85±18m0.9) have similar average ksn indices.

The knickpoints, which separate channel segments with different steepness values, are typically of the slope-break type. We estimated knickpoint elevations relative to the river outlet, which varies widely from ~200 m to 1200 m (Figure 5(a)). Minimum incision magnitude also ranges widely (~70–670 m). Such variations may indicate spatial variations in an uplift pattern. Nevertheless, the χ values, most of which lie within a relatively narrow range of 5.59±1.80m, meet the criteria of a normal distribution (Figure 5(b)). We can thus consider that the χ values are almost constant, implying that these knickpoints might have begun to form and migrate at the same time.

5. Discussion

5.1. Tectonic Controls on Channel Steepness Patterns and Knickpoint Migration

The steepness indices are higher in the Huaxian-Huayin section than in other parts, which appears to indicate focused tectonic uplift around that area. The slope-break knickpoints, which separate the steeper downstream from the more gentle upstream reaches, might imply a temporary increase in uplift rates. However, to relate the geomorphic features to tectonic control, detailed discussions on the influences of other potential variables, such as rock resistance, sediment flux, and drainage network reorganization, are needed.

5.1.1. Bedrock Resistance

Rock resistance to stream erosion can exert important effects on channel erodibility. Spatial variations in substrate properties can cause different river incision processes and thus considerable scatter in the slope-area plot [66]. Along channels with spatially variable incision processes, vertical-step knickpoints may develop below (and very close to) slope-break knickpoints, resulting in abnormal concavities, e.g., negative values or much higher-than-average concavity. The knickpoint on river No. 33, for example, is situated at the transition between the Archean gneiss and granite (Figure 3), with an abrupt increase in the downstream channel gradients and negative concavity. This example shows the influence of lithology on changes in channel attributes. However, we did not find any obvious knickpoints at other lithologic boundaries along the upstream and downstream reaches of rivers, e.g., Nos. 9 and 25. For most rivers, knickpoints furthermore do not correlate with mapped lithologic contacts and exist in a variety of rock types (e.g., Nos. 4 and 45, Figure 6), suggesting that rock type is not the dominant factor controlling knickpoint formation.

To assess the potential role of lithologic resistance on channel gradients, we compared the channel steepness pattern with mapped variations in regional lithology (Figure 3). Along the southern Weihe Graben, many rivers flow through areas of granite (especially in the following segments: Baoji-Zhouzhi, Lantian-Huaxian, and Huaxian-Huayin). Channel steepness nonetheless varies a lot (~30–270 m0.9) despite the uniform lithology. The granite is characterized by sparse jointing or foliation, which should promote low channel erodibility (K) and thus high steepness indices [5]. However, some rivers on granite also exhibit low steepness values. For example, steepness indices of river Nos. 2–4 (~110 m0.9, granite) are much lower than that of Nos. 14–16 (~130–180 m0.9, Cambrian limestone, sandstone). Thus, in this landscape, lithologic resistance to river incision appears to be limited to local channel reaches.

5.1.2. Drainage Network Evolution

In active orogens, river networks evolve as a response to river capture and drainage divide migration [59, 67]. In the upper reach of a channel gaining drainage area by river capture, marked features like a steep gorge will develop above the capture point in response to the rapid incision [59, 68]. Usually, such a feature will be of the vertical-step knickpoint variety and is easily identified from the river longitudinal profile and χz plot [59, 68]. No such geomorphic features, however, were encountered along these rivers in our study area.

Willett et al. [59] proposed that the χ-value difference across the catchment divide could indicate disequilibrium river networks even without precise constraints on the rock uplift pattern. Taking a threshold for drainage area of 1 km2, we calculated the χ-value map of the catchments to gauge the horizontal motion of drainage divides in the region (Figure 7). We found a disequilibrium drainage network in this landscape. For example, along the Baoji-Zhouzhi segment, catchment Nos. 6 to 9 are likely growing at the expense of their neighbor No. 10. At the eastern termination of the Huayin-Lingbao segment, the main channel of the victim catchment No. 55 is losing area to drainage basin Nos. 44 to 54.

During drainage divide migration, the conquering (aggressor) catchments will have higher incision rates, while the incision rates for those losing areas (victim) are expected to become lower [59]. Both could exert significant influence on bedrock channel morphology and even upset interpretations of the rock uplift process from river profiles [3]. Thus, if under the first-order control of the disequilibrium drainage network, then we would observe a significant correlation between channel steepness and catchment area. However, weak correlation between channel steepness and catchment area indicates limited influence of changes in the drainage area (Figure 8; upstream, R2=0.11, P value = 0.06; downstream, R2=0.00, P value = 0.84; steady-state rivers, R2=0.01, P value = 0.71). Thus, we exclude the dominant control of the evolved drainage network on the ksn pattern.

5.1.3. Sediment Flux

Sediment flux (Qs) exerts an influence on channel morphology either by protecting the channel bed from being more eroded (Qs>Qc, Qc is the sediment carrying capacity) or by increasing available tools (Qs<Qc). The former will lower channel gradients and river incision rates, imparting a high concavity. Conversely, the latter can induce higher incision rates and gradients, leading to lower or negative concavity values [5]. Furthermore, big earthquakes, massive landslides, glacier melting, and/or high magnitude, low-frequency flood events could deliver large amounts of sediment to river system, resulting in a series of cascades at local channel segments. River profiles would depart from a simple detachment-limited system, characterized by large scatters in channel slope and a wide range of channel concavities. In this area, some channel reaches with a wide range of concavities or large scatter in channel gradients, e.g., Nos. 9 and 25 (Figure 6), might serve as an example.

We found at least three reasons for ruling out any major control on channel profiles by sediment flux. Firstly, all of the reaches with large concavities are only situated below knickpoints. Above knickpoints, in contrast, channel gradients exhibit a good power-law scaling, with nearly uniform concavity indices. There is no evidence to demonstrate that the influence of sediment flux should be focused exclusively on the downstream channels. Secondly, the slope-break knickpoints retreat as incision waves and are thus not static (see details in Section 5.1.4). In some cases, channels with distinct variations in local gradients develop immediately below slope-break knickpoints as a response to variation in incision processes ([5]; Whipple et al., 2011). Although all the downstream reaches have much higher steepness indices than the upstream sections, the concavity indices do not correlate with systematic variation in channel steepness. That means variations in river incision rates cause changes in steepness indices rather than in channel concavities. In addition, channels with larger drainage areas are usually likely to receive more sediment [6, 69, 70]. We thus turned to the drainage area to make a simple estimate of the effects of sediment flux. In this landscape, however, we found no correlation between the drainage area and channel steepness (Figure 8). The influence, if any, of sediment flux on knickpoints and channel steepness pattern is thus overall minor.

5.1.4. Tectonic Control of Rock Uplift

Because well-exposed fault outcrops and fault planes are scarce in the study area, limited data concerning fault slip rates have been reported. This limits our knowledge of the quantitative relationship between channel steepness and rock uplift rate. Several lines of evidence nonetheless support the conclusion that the channel steepness and knickpoint patterns reflect tectonic control.

Firstly, there is a clear segmentation along the strike, with downstream steepness indices alternating between high and low values every 100 or 200 km (Figure 4(d)). At a larger scale, this segmentation correlates with the relief of the mountain range and somehow corresponds to varying fault strike directions (Figure 3). At Baoji, a boundary city between Sections 1 and 2, the fault changes its orientation from W-E to NW-SE (Figure 3). At Zhouzhi, the fault nearly tends to be W-E despite some curves near Huxian (a boundary identified by [31]). Lantian city could be the division between NQMF and HPF. From Lantian to Huaxian, the piedmont fault strikes SW-NE. Then, it tends to be ENE in the Huaxian-Huayin segment. In the eastern segment, Huayin-Lingbao, the fault nearly stretches W-E. Observing variations of ksn along the southern Weihe Graben (Figure 3), we notice that fault linkages where the strike direction changes occur always have lower channel steepness. Chen et al. [31] also identified tectonic segment boundaries based on other geomorphic indices, e.g., hypsometric integral, catchment asymmetry factor, stream-length gradient, and valley width/height ratio, which are generally consistent with our study. Although the entire Qinling range has been recording major strike-slip-related deformation (Figure 1), the rapid uplift of the northern Qinling range is mainly controlled by a pure normal-fault system [39]. There are many high fault scarps along the mountain front, but with barely any geomorphic evidence for lateral-strike-slip motion such as offset stream channels, terrace risers, or mountain ridges [33, 35]. The long faults (NQMF and HPF) were characterized by a succession of shorter normal segments during the Eocene-Oligocene [33] and then linked together in the late Miocene [39]. Thus, it seems likely that changes in fault strike may not induce a variable amount of strike-slip motion or mountain uplift. We propose that the periodic variation in channel steepness might indicate that, for each segment, fault throw rate decreases from the fault center to each of its tips. Despite the lack of a more detailed fault trace mapping, we suspect that the fault throw rates (or channel steepness indices) could be more strongly correlated with the fault dip angle. For example, along the Huaxian-Huayin segment (with the highest ksn values), the angle attains the highest value of 65–80°, whereas in the east or the west portions, the dip angle is no more than 68° (the data of fault angle are from [71]).

Secondly, channel steepness indices are much higher in the Huaxian-Huayin segment than in others. Among the earthquakes recorded by historical and paleoseismological studies, the 1556 Ms~8.5 Huaxian earthquake, caused by the slip of Huaxian-Huayin fault segment, was the most devastating. The intense activity of this fault segment can also be indicated by the pronounced fault scarps and bedrock fault planes, which have not been encountered in other segments along the southern Weihe Graben [31, 35]. Additionally, the faulted terraces with exposure ages of 2–3 ka [33] were conspicuously active during the late Holocene. The steepness pattern is thus consistent with existing studies of local paleoseismology and active tectonics.

Thirdly, we found strong correlation between channel steepness and mountain topography. We plotted ksn indices against average catchment slope, topographic relief (a series of local scales, 1.5, 2.5, 3.5, 4.5, and 5.0 km radius, and catchment-wide), and mean elevation (Figure 9). For upstream channel steepness, we obtained moderate correlation (R2 from 0.28 to 0.47) when the scale radius was <4.5 km. The largest value is obtained (R2=0.47) for a radius of ~1.5 km. Catchment-wide relief (R2=0.14) and average elevation (R2=0.07), meanwhile, were uncorrelated with upstream ksn. However, Figure 9 shows a strong correlation between downstream steepness and topographic relief and basin elevation. Except for the average slope (R2=0.35) which is a proxy for local relief within a rectangle of 3×3 pixels, nearly all the coefficients (R20.47 and mostly ≥0.49) indicate high correlation. Especially for average elevation, the highest coefficient (R2=0.7) indicates a strong geographic association between downstream steepness and the tectonically active range front.

In addition, the identified knickpoints develop as a response to an increase in tectonic uplift rates and retreat as an incision wave. To assess knickpoint recession behavior, we followed three different deterministic models [72], i.e., knickpoint retreat distance vs. (1) drainage area of the catchment, (2) total length of the main stream, and (3) length-area exponent (also known as Hack’s exponent; [13]). Figure 10(a) shows a power-law relationship between knickpoint retreat length and drainage catchment area (R2=0.79, P value ≤ 0.01) rather than with the present positions of knickpoints (R2=0.17, P value = 0.08), indicating a commonly reported case where catchment area ranks as a first-order control on how far the knickpoint has migrated [21, 24, 73, 74]. The fittedexponent~0.49 is similar to case studies (0.4–0.5) in active orogens where the knickpoints had all started at the same time [72, 74, 75]. Nevertheless, the limitation of the model is that it cannot account for the losses of the drainage area. He et al. [36] pointed out that the drainage network was not stationary and recognized disequilibrium water divides, which might reduce the ability of using the drainage area as a key parameter to determine knickpoint migration.

As proposed by Castillo et al. [72], stream length, a proxy for water discharge, could be appropriate for predicting knickpoint retreat, because the model is insensitive to losses of the drainage area. In this study, both area- and length-dependent models exhibit the same high correlation coefficients. Furthermore, the knickpoint retreat distance depends on stream length in a quasilinear form (R2=0.79, P value ≤ 0.01; Figure 10(b)). We were thus able to infer (1) that the area gain or loss caused by drainage evolution and knickpoints migrating through tributaries exerts little influence on knickpoint migration and (2) that the drainage area can also be a good proxy for water discharge (an important assumption of the stream-power incision model). The model based on Hack’s exponent allows the dependence of knickpoint migration on the degree of elongation of a river basin to be explored. The poor correlation coefficient (R2=0.11, P value = 0.10; Figure 10(c)), however, indicates that the upstream migration of knickpoints does not relate to whether the basins are broad or narrow. Accordingly, although variations in channel concavities might indicate nonuniform local steepness indices along the downstream reaches, these slope-break knickpoints are likely records of the transition from less intense to more intense tectonic activity.

Based on the discussions above, we conclude that variations in steepness index for both downstream and steady-state channel reaches are suitable for mapping the differential tectonic uplift pattern. The fault segment of Huaxian-Huayin is characterized by the highest vertical slip rate. These knickpoints develop as a response to rapid acceleration in the mountain uplift rates along the southern Weihe Graben. Knickpoints migrate as a transient kinematic wave, allowing the onset timing of rapid river incision (or active mountain uplift) to be constrained. Considering the correlation analysis shown in Figure 9, we conclude that the modern topographic relief and average elevation is likely under the control of currently active tectonic forcing, while the previous imprints on the landscape might be mostly erased by the recent tectonic activity.

5.2. Knickpoint Ages and Implications for Mountain Growth

In many cases, river channels increase their gradients to accommodate the increase in the uplift and/or denudation rates. Usually, as drainage basins with an average slope of 5-25°, one can observe a power-law scaling between the rates of the base-level fall with topographic metrics, e.g., average slope and channel steepness indices [76]. This produces constraints on both slope exponent (n) and channel erodibility (K), which is key to determining knickpoint ages. But, when the basin slope is close to or above the threshold angle (25–30°), erosion is stochastic in nature, and any correlation between denudation rates and catchment metrics falls apart [77]. In the northern Qinling range, owing to the steep topography (average basin slope of ~25–32°; Figure 9(a)), He et al. [36] found little correlation between denudation rates and either basin slope or channel steepness.

We used different n values (1, 1.5, and 2) to calculate erodibility (K) and knickpoint ages. When n<1, an increase in tectonic uplift rates would cause a stretch zone instead of slope-break knickpoints [56]. In the χ domain, the stretch zone exhibits a curve lacking tectonic information, which is thus not compatible with the geology of this landscape. Thus, n<1 was ruled out. Based on the catchment-wide denudation rates (from [36]) and average steepness indices (Figures S2 and S3) and the linear n, we estimated erodibility values for each basin and derived a plot of K values with distance along the boundary of the northern Qinling range (Figure 11(a)). The low coefficient (R2=0.10) between K and distance indicates that channel erodibility exhibits nearly no explicit spatial variations and thus allows an estimate on the average value of erodibility to be made. Here, K~1.00±0.44×106m0.1/a. Combing the χ values, we derived knickpoint ages and obtained a mean value of 5.59±1.80Ma (Figure 11(b)). This timing defines the onset for rapid river incision or mountain uplift of the northern Qinling range. Taking the minimum incision magnitude of 311±162m (Figure 5(a)), we derived an incision rate of ~0.06±0.04mm/a which is lower than the catchment-wide average denudation rates (mostly 0.1–0.2 mm/a; [36]).

We propose two reasons why the estimated rate is low. Firstly, denudation rates documented by 10Be dating could be overestimates because of landslides in the catchments. For example, many landslides have occurred along the northern flank of the Huashan, which would cause higher erosion rates [33]. Secondly, the local channel may not have sufficient time to evolve as responses to the rapid throw of the active fault, especially in the Huaxian-Huayin fault segment. Through field investigation along nearly the entire Huaxian-Huayin segment and 14C and OSL dating of faulted terrace, Du et al. [78] and Li et al. [33] constrained the throw rate of the fault segment to be 0.61±0.15mm/a since the middle Pleistocene. During the Holocene, fault slip rates increased to 2–3 mm/a [33, 34]. Abrupt increases in such a short time scale (from the middle Pleistocene to the Holocene) generated pronounced fault scarps and bedrock fault planes, with fluvial incision failing to balance the rapid tectonic uplift. Therefore, although we derived the onset of acceleration in river incision from knickpoint positions, it remains difficult to retrieve details of uplift rate variations through time just from the channel profiles.

We also used the ksn data of He et al. [36] and nearly produced the same K and knickpoint ages (Figures 11(a) and 11(b)). According to Figures 11(c)–11(f), a higher n causes a lower magnitude of K (if neglecting its unit) and younger knickpoint age, e.g., n=1.5 and n=2 produce late-Pliocene to early Quaternary ages (~2–4 Ma). This is much younger than age constraints provided by low-temperature thermochronology. Thus, based on the metrics of knickpoints, we suggest that a linear exponent of n=1 would be better. Equation (8) shows that when the river incision rate is a linear function of the local channel gradient, knickpoints that started to migrate at the same time should have a constant χ value (also demonstrated by [56]). Additionally, we found little correlation between knickpoint retreat distances with channel steepness (R2=0.00, Figure 10(d)). In the xz coordinate, the horizontal velocity is dependent on the contributing area rather than the channel gradients (when n=1, Equation (10)). Given that fluvial systems adjust channel gradients in response to changes in uplift rates, the distance over which a knickpoint has migrated upstream does not relate to downstream channel steepness (a proxy for recent tectonic uplift). Accordingly, we conclude that the slope exponent could be linear in this landscape.

We now discuss whether the knickpoint age derived from the assumption of linear n is reasonable. Studies on the sedimentology and paleogeography reveal that, during the Eocene-Oligocene, the Weihe Graben was separated by a series of small rift basins bounded by short normal faults [33]. These faults started to link together in the late Miocene, resulting in a large-scale normal fault system that controls the prominent topographic relief. Fault linkage and active normal faulting caused rapid mountain growth and basin subsidence. Subsequently, in the early Pliocene, those small basins started to connect to one another to form Lake Sanmen, which was much more extensive than the present Weihe Graben. A drill core near Huxian city has revealed an ~3 km thick depositional stratigraphy of late Miocene-early Pliocene age, which amounts to about one-half of the entire Cenozoic basin fill [52, 79]. During the middle Miocene (16–11 Ma), deposition did not exceed 1 km [39, 50]. Late Miocene mountain growth of the northern Qinling range can thus be inferred from this evidence. The knickpoints are thus likely a response of the channel profiles to active fault growth and basin subsidence.

Forward modelling on vertical transects of low-temperature thermochronology (AFT and AHe) revealed that the onset of rapid mountain growth began around 7 Ma [42], i.e., a little earlier than our estimates. Inversion on channel profiles can retrieve tectonic and/or climatic details relevant to 105–107 yr time scale, but these details would be difficult to match up with results from low-temperature thermochronology. Uncertainties may arise from the diversity of knickpoint χ values, from errors in the estimated stream-power parameters, as well as from the thermal modeling strategies. However, both our work and previous studies [42, 43] indicate rapid uplift of the northern Qinling range since the late Miocene and early Pliocene.

Owing to intracontinental extension in the north Qinling range-Weihe Graben area, the Qinling Mountains have experienced two episodes of rapid exhumation: one during the Oligocene and the other from late Miocene to the Present. The latest stage of rapid exhumation is attributed to the growth and eastward expansion of the NE Tibetan Plateau. Along the NE margin of the plateau, thermochronology results have revealed middle to late Miocene strike-slip fault displacement (5–8 Ma along the eastern segment of the Kunlun fault, 8–10 Ma at the eastern tip of the Haiyuan fault; [80]) and rapid mountain uplift (~10 Ma of the northern and eastern Qilian Shan, ~8 Ma at the Liupan Shan; [12, 81, 82]).

Thus, compared with the initial ages of rapid topographic growth along the NE plateau margin, the timing for the Qinling Mountains and other areas resulting from the eastward expansion of the NE Tibetan Plateau should be a bit later. Based on analysis of fission track ages of detrital apatite samples in the Linxia Basin, where sediment provenance is from the western Qinling range in the south and the Jishi Shan in the west, Zheng et al. [83] documented two stages of uplift in the Jishi Shan: one at 14 Ma (attributed to convective thinning of the lithosphere in the NE Tibetan Plateau) and another at 5–8 Ma (the eastward expansion of the plateau). A much younger age of 3–4 Ma was also identified by Zheng et al. [83], which could indicate uplift of the Jishi Shan as late as the Pliocene time. Again based on AFT modelling, Enkelmann et al. [40] also described a somewhat younger onset of rapid cooling (9–4 Ma) for the SW Qinling range. Thus, the derived knickpoint ages of 5.59±1.80Ma appear to be a reasonable estimate of the initial age of the Qinling Mountain growth. In reference to Section 5.1.4, we conclude that the metrics of mountain topography (average elevation, local, and catchment-wide relief) have been mainly controlled by tectonic forcing since the late Miocene and early Pliocene.

6. Conclusions

In this contribution, we analyzed longitudinal profiles of streams draining the northern Qinling range. We identified slope-break knickpoints and calculated concavity and steepness indices. Concavity values for upstream channel reaches and whole catchments mostly lie within a narrow range of 0.2–0.6, which matches a normal distribution and thus indicates a similar shape under the normalized equilibrium profiles. However, channel steepness values differ substantially, with higher values below the knickpoints as opposed to more gentle upstream segments. For the downstream reaches, ksn attains its highest values regionally between Huaxian and Huayin. Such a pattern indicates that active tectonic uplift is mostly focused along this fault segment. These knickpoints develop as a response to a sustained increase in mountain uplift rates and migrate upstream as an incision wave, the retreat velocity of which is dependent on the catchment area and stream length. Based on power-law scaling properties between channel steepness and catchment-wide denudation rates, we determined stream-power parameters and estimated knickpoint ages of 5.59±1.80Ma. From this, we conclude that the different ksn patterns above and below existing knickpoints indicate rapid mountain growth since the late Miocene and early Pliocene.

Data Availability

The DEM data used in this study is publicly free (http://earthexplorer.usgs.gov/).

Additional Points

Key Points. There is tectonic control on river profiles and the present topography in the northern Qinling range. Knickpoints develop as a response to temporary increases in the uplift rates and retreat as a kinematic wave. There has been an acceleration of the extension in the northern Qinling range-Weihe Graben since the late Miocene and early Pliocene.

Conflicts of Interest

The authors declare that they have no conflict of interest.

Acknowledgments

We thank He Chuanqi (Zhejiang University) for providing the 10Be data used in their paper [36]. We are grateful for the grants from the National Science Foundation of China (41802227, 41772198).

This is an open access article distributed under the Creative Commons Attribution License, which permits unrestricted use, distribution, and reproduction in any medium, provided the original work is properly cited.