Lithospheric Flexure Controls on Geomorphology, Hydrology, and River Chemistry in the Andean Foreland Basin

Abstract

Tectonics exerts a strong control over the morphology of Earth's surface that is apparent in active mountain belts. In lowland areas, subtle processes like lithospheric flexure and isostatic rebound can impact Earth surface dynamics, hydrologic connectivity, and topography, suggesting that geomorphic and hydrologic analyses can shed light on underlying lithospheric properties. Here we examine the effect of lithospheric flexure on the geomorphology, hydrology, and river water chemistry of the Rio Bermejo fluvial system in the east Andean foreland basin of northern Argentina. Results show that proximal to the mountain front, foredeep basin subsidence causes sedimentation along a braided channel belt that is superelevated relative to the surrounding flood basin. During floods, water flows from the superelevated channel into the groundwater reservoir, causing a net loss of discharge with distance downstream. Further downstream, forebulge uplift forces channel narrowing, high lateral migration rates, and incision up to 13 m into older river deposits. This incision locally allows groundwater flow into the river, causing a ∼20% increase in river solute load. Groundwater emerges from the forebulge into the backbulge, predominantly as spring-fed channels. Here, channel migration rates decrease, suggesting a switch from net uplift to subsidence that reduces the depth to the groundwater table. This analysis shows that subtle lithospheric flexure can have significant effects on river channel morphology that determine hydrologic flow paths, and ultimately influence geochemical and ecological patterns. We suggest that these effects may elucidate lithospheric properties that are otherwise inferred from bulk geophysical observations.

Introduction

Tectonics and earth surface processes are linked. This is perhaps most apparent along actively deforming continental margins. Mountain building sets patterns of surface uplift and precipitation (e.g., Anders et al., 2008; Willett, 1999), can both drive and be modified by spatial variations in erosion and deposition (e.g., Allen et al., 1991; Koons et al., 2002; Oberlander, 1965; Wang et al., 2014; Zeitler et al., 2001), and influences weathering rates and geological carbon fluxes that modulate global climate (e.g., Raymo & Ruddiman, 1992; Chamberlin, 1899; France-Lanord & Derry, 1997; Hilton & West, 2020; Bufe et al., 2021). Similarly, the interplay of tectonics and surface processes has given rise to distinct and economically important patterns of drainage, erosion, and deposition along rifted margins (Kooi & Beaumont, 1994; Steckler & Omar, 1994; Weissel & Karner, 1989). Although most work examining feedbacks between tectonics and surface processes has focused on relatively high-relief landscapes (e.g., Burbank & Anderson, 2013; Finnegan et al., 2008; Kirby et al., 2007; Wobus et al., 2006), tectonic processes can also influence topography, hydrology, and biogeochemical cycling in low-relief landscapes. For example, lithospheric flexure and mantle flow can have a major impact on the location and direction of lowland rivers (e.g., Cox, 1989; Wickert et al., 2019). This is because low-relief rivers tend to have geometries that are finely adjusted to their discharge and sediment load (Mackin, 1948), such that small amounts (∼10¹ m) of uplift or subsidence can cause major channel instability (e.g., Burnett & Schumm, 1983; Gasparini et al., 2016).

Low-relief rivers convey much of the world's mobile sediment from upland sources to geological sinks (Milliman & Meade, 1983) and act as filters that can distort and obscure signals of tectonic or climatic perturbations during sediment transport (Jerolmack & Paola, 2010; Straub et al., 2020). Furthermore, lowland rivers construct landscapes that promote weathering and soil formation (e.g., Willenbring et al., 2013), thereby sustaining the vast majority of global human settlement and economic activity (Doyle, 2018; Hassan, 1997; Macklin & Lewin, 2015). Despite this importance, there have been relatively few studies examining linkages between tectonics and surface processes in low-relief settings (e.g., Allen & Homewood, 1986; Armitage et al., 2014; McGlue et al., 2016; Roberts & White, 2010; Ruetenik et al., 2016; Wickert et al., 2019).

Foreland basins represent one of the most commonly occurring low-relief settings where one might expect tectonics to influence surface processes. In these depositional settings, crustal thickening within the adjacent mountain belt creates a local load on the Earth's lithosphere, causing it to flex like a rigid beam (DeCelles, 2012; Karner & Watts, 1983; Turcotte & Schubert, 2002; Walcott, 1970). This creates a region of subsidence, the foredeep, immediately adjacent to the mountain range, paired with a more distal region of minor uplift (the forebulge), and a second region of subsidence (the backbulge) most distal to the mountain range (e.g., DeCelles, 2012; Figure 1a). The horizontal and vertical scales of the flexural profile are set by the local flexural rigidity, the elastic thickness of the lithosphere and the dimensions of the orogenic load (Chase et al., 2009; Turcotte & Schubert, 2002; Watts, 2001). The amplitude of subsidence and uplift along the profile decays rapidly with distance from the orogenic load, such that the amplitude of backbulge subsidence can be up to three orders of magnitude less than foredeep subsidence (e.g., DeCelles, 2012; DeCelles & Giles, 1996). Foredeep subsidence can be rapid, aided by the load of accumulating sediment from the adjacent mountain belt, to attain depths of 10³ m, and widths of 10¹–10² km. In contrast, forebulge uplift rates are slower and vertical amplitudes tend to be ∼4–7% of foredeep maxima, resulting in ∼10² m of positive flexure. Subsidence in the backbulge is typically on the order of 10⁰ m, fading into a (cratonic) far field background (e.g., DeCelles, 2012; DeCelles & Giles, 1996). Predictions of the dimensions of these foreland basin components are imprecise because lithospheric rigidity and elastic thickness cannot be measured directly (e.g., Tesauro et al., 2012). Moreover, the topographic expression of foreland flexure can be heavily modified by sediment supply (Garcia-Castellanos & Cloetingh, 2012; Sinclair et al., 1991). When sediment supply is limited, the foredeep is underfilled and the forebulge may form a notable topographic high (e.g., the European Alps; Sinclair, 1997). In contrast, under high sediment supply, the foreland becomes overfilled, forming a sediment apron which covers flexural features (e.g., the central Andean foreland; Cohen et al., 2015; McGlue et al., 2016; Garcia-Castellanos, 2002), making these features difficult to delineate, and possibly damping their influence on surface processes.

Generate image description

Figure 1. Overview map of the Rio Bermejo study area. (a) Schematic cross section of the Andean foreland basin subsurface structure and lithospheric dynamics in the region of the Rio Bermejo, modeled after Horton and DeCelles (1997). (b) South America geoid anomaly, as calculated by Chase et al. (2009), overlain on shaded relief. Inset map shows geoid anomaly map for broader region in South America. Color scale is shown in legend, with warm colors representing positive anomalies and cool colors showing negative anomalies. Also shown are major rivers, hydrologic gauging stations (PS = Pozo Sarmiento, EC = El Colorado), and water sampling sites (with downslope distance from the mountain front to the Rio Paraguay shown at select sites along the main stem Rio Bermejo). Dark shaded gray rectangles show locations of panels c, d and e, which show Google Earth satellite images of the Rio Bermejo: (c) braided reach in the foredeep at ∼80 km downstream from the Rio San Francisco (RSF) confluence, (d) meandering reach in the forebulge at ∼245 km downstream from the RSF confluence, and (e) meandering reach in the backbulge at ∼650 km downstream from the RSF confluence.

In underfilled forelands, rivers draining the mountain belt may be forced by the forebulge topographic high into a strike-parallel course located within the foredeep (e.g., Clevis et al., 2004; de Leeuw et al., 2020). In contrast, in overfilled forelands, rivers can maintain a course perpendicular to the strike of the mountain range, crossing the major flexural elements of the foreland, while building thick deposits in subsiding areas and thinner deposits over the forebulge. In the foredeep, these rivers should elevate their bed by sediment deposition, raising the channel above the surrounding floodplain and causing frequent river avulsions (e.g., Edmonds et al., 2016). Downstream, they should incise into the uplifting forebulge, giving rise to channel belts inset below the surrounding alluvial terrain (e.g., Bufe et al., 2016; Holbrook & Schumm, 1999). We hypothesize that under these conditions, spatial patterns of river channel superelevation versus incision can reveal the locations of flexural features, permitting estimation of lithospheric mechanical properties from measurable geomorphic parameters, like channel width, sinuosity, and channel elevation relative to the floodplain.

Moreover, we expect that the vertical position of the river channel relative to the floodplain can drive changes in groundwater dynamics, with implications for river chemistry and biogeochemical cycling. Specifically, in superelevated channels the channel is a topographic high, requiring a hydraulic gradient sloping away from the channel. In this case, we expect overbank flow and outward seepage, which would cause a net water flux from the river into the groundwater reservoir (Winter et al., 1998). This discharge loss should leave the river water chemistry unaffected. In contrast, areas where the channel is incised into the floodplain the channel is a topographic low, generating a hydraulic gradient sloping toward the channel and a positive flux of groundwater into the river channel (e.g., Winter et al., 1998). In this case, we expect an increase in river dissolved solute concentrations due to the input of solute-rich groundwater (e.g., Maher, 2011). This highlights how lithospheric flexure may regulate the weathering flux exported from mountain forelands.

We explore the idea that, even in the case of a filled foreland basin, lithospheric flexure can influence groundwater-surface water connectivity by modulating the vertical position of the river channel relative to the groundwater table. Under this assumption, we can analyze changes in river discharge and water chemistry with distance downstream across the foreland basin to identify otherwise subtly expressed lithospheric dynamics. We test this hypothesis in the east Andean foreland basin of northern Argentina, where the Andean orogenic load causes lithospheric flexure radiating to the east. We first provide background on lithospheric flexure and describe our study area, the Rio Bermejo, which drains across the ∼700 km wide Andean foreland as a single channel without significant tributaries or distributaries, building a narrow fluvial megafan (McGlue et al., 2016; Repasch et al., 2020). This uniquely simple drainage pattern allows us to avoid complicated tributary mixing scenarios that arise in dendritic drainage networks, and to isolate the effects of tectonic forcing. Second, we present a topographic analysis that shows the Rio Bermejo is superelevated relative to the floodplain in areas of flexure-induced subsidence and incised into the floodplain in areas of flexural uplift, with consequences for sediment deposition and channel mobility. Third, we analyze river discharge records and dissolved load geochemistry to show that lithospheric flexure causes a net loss of river discharge moving downstream and increases the river solute load due to an influx of groundwater where the Rio Bermejo traverses the forebulge. Taken together, these metrics provide a detailed documentation of lithospheric flexure influencing fluvial geomorphology, groundwater hydrology, and river dissolved load geochemistry, as well as a means of estimating the local mechanical properties of the lithosphere.

The East Andean Chaco Foreland and the Rio Bermejo

The Andean foreland basin is a type example of a retroarc foreland basin (e.g., DeCelles, 2012; Horton & DeCelles, 1997). A combination of seismic, topographic, and sedimentological data shows that the Chaco foreland exhibits a typical four-part system consisting of a wedgetop, foredeep, forebulge, and backbulge (Figure 1a, Horton & DeCelles, 1997; McGlue et al., 2016; Cohen et al., 2015). High sediment supply from the eastern Andes has created a series of fluvial megafans in the Chaco foreland (e.g., Thalmeier et al., 2021), making this an (over-) filled foreland basin that lacks an obvious topographic expression of forebulge uplift. However, evidence for forebulge uplift (and foredeep subsidence) exists from geophysical observations and modeling. Along the central Andes of Argentina and Bolivia (between −20°S and −30°S), the foreland basin exhibits a positive geoid anomaly, corresponding to ∼550 m of positive flexure (Figure 1b), in the form of a forebulge with a width of ∼200 km, positioned ∼300–500 km from the mountain front (Chase et al., 2009).

We follow previous work (McGlue et al., 2016) in using the Rio Bermejo to examine how zones of subsidence and uplift influence Earth surface dynamics within a filled foreland basin (Figures 1a and 1b). With a drainage area of ∼120,000 km², the Rio Bermejo is an undammed river draining the eastern flank of the Andes in northern Argentina at ∼24°S. The Bermejo traverses across the East Andean Chaco foreland to join the Rio Paraguay (and, further downstream, the Rio Paraná). The Bermejo river system has been active since at least the Pleistocene, building a megafan of fluvial deposits across the lowland basin (McGlue et al., 2016) with a fan apex-to-toe relief of ∼250 m. We focus on the section of the Rio Bermejo between the mountain front and the Rio Paraguay; the high average sinuosity of the river results in a total channel length of ∼1,270 km, despite only a ∼700 km straight-line distance along this flow path. Regional precipitation occurs mostly in the austral summer (November-April), with a high rainfall rate along the mountain front (∼1,100 mm/yr), low rainfall from ∼150 to 300 km east of the mountain front (∼650 mm/yr), and high rainfall at the Rio Bermejo-Rio Paraguay confluence (∼1,400 mm/yr). The river experiences peak water discharge of ∼1,000–2,000 m³/s during February and March, while low flows during September through November can fall below 50 m³/s (Golombek et al., 2021; Sambrook Smith et al., 2016). Between 1990 and 1999, mean sediment flux recorded at the El Colorado gauging station (∼580 km downslope of the mountain front, Figure 1) averaged 110 Mt/yr, which is equivalent to 60%–90% of the total wash load of the upper Rio Paraná (Amsler & Drago, 2009).

As part of a larger sedimentary systems analysis, McGlue et al. (2016) examined the expression of lithospheric flexure within the Rio Bermejo megafan. The McGlue et al. (2016) study showed that the Rio Bermejo transitions from a wide (up to 3 km), braided channel in the foredeep to a relatively narrow (mean width ∼270 m), single-threaded meandering channel in the forebulge. Moving from the forebulge to the backbulge, the Rio Bermejo remains single-threaded, but further narrows to ∼100–200 m wide. Using satellite imagery and residual topography relative to the local riverbed elevation, previous work (Cohen et al., 2015; McGlue et al., 2016) has shown an abundance of paleochannels and topographic lows in the Rio Bermejo foredeep, indicating a high degree of channel mobility and considerable accommodation space. In the forebulge, residual topography is high and there are fewer paleochannels, suggesting decreased channel mobility and reduced accommodation space. In the backbulge, abundant topographic lows and spring-fed channels flowing parallel to the mainstem Rio Bermejo indicate an increase in accommodation space relative to the forebulge. The transition from a braided to meandering river at the foredeep-forebulge transition has also been shown to drive lateral erosion into relict sediment deposits that have been stored in the floodplain for up to ∼20 kyr (Repasch et al., 2020).

The observations of Cohen et al. (2015), McGlue et al. (2016), and Repasch et al. (2020) are consistent with lithospheric flexure influencing the geomorphic expression of river networks. Below we expand on these analyses to develop a deeper, process-based explanation of how this lithospheric flexure directly influences geomorphic processes and in turn modulates the hydrology and dissolved load geochemistry of the Rio Bermejo.

Topographic Expression of Lithospheric Flexure

We begin by examining the topographic and river morphometric trends along the Rio Bermejo and evaluating links between these trends and foredeep subsidence and forebulge uplift.

Methods

We extracted the longitudinal profile of the Rio Bermejo by hand from the TanDEM-X DEM (12 m resolution) in ArcMap. Widths of the channel and active channel belt were calculated by first manually tracing the edges of these features as identified in 2015 Landsat 7 imagery (30 m resolution). We then used the ChanGeom algorithm (Fisher et al., 2013) to determine the channel centerline and calculate channel and channel belt widths at 10 m intervals.

We quantified the superelevation and incision of the Rio Bermejo relative to the surrounding topography by making a residual elevation map relative to the active channel (sensu Cohen et al., 2015; McGlue et al., 2016). We created a series of 1 km wide swaths set perpendicular to the Rio Bermejo flow direction, mapped the banks of the Rio Bermejo by hand from aerial imagery overlaying the TanDEM-X DEM, and extracted the minimum elevation of the active Rio Bermejo channel within each swath from the TanDEM-X DEM. Variability in DEM values for water surfaces created ∼2–4 m noise in elevation values of the active channel, which were smoothed by fitting a fourth order polynomial (R² = 0.999) to the extracted active channel elevations to estimate the downstream change in elevation. We subtracted the polynomial fit from the DEM to yield a DEM that was de-trended from the Rio Bermejo channel slope, with values representing the elevation of the surrounding terrain relative to the river channel. We note that the TanDEM-X DEM is sensitive to vegetation, resulting in levee elevations that appear amplified due to tree height.

To measure lateral migration rates, we mapped the active channel from Landsat imagery for four different years, yielding three time-intervals (1988–1993, 1993–2002, and 2002–2015). For each interval, we compared the locations of the two channel centerlines. Between successive crossing points of the two centerlines, we measured the area between the two lines, which represents the area traversed by the channel centerline over the given timespan. This area, divided by the along-channel distance and the timespan between images, yields a lateral migration rate for that channel segment. We averaged these rates over 15 km long channel segments. This length scale was chosen so that, along the meandering reach of the river, at least three individual meander bends compose a single channel migration value.

Results and Interpretation

Our analysis of channel morphometrics along the Rio Bermejo is shown in Figure 3 (and reported in Tables S1 and S2 in Supporting Information S2), where we plot data against straight line distance downslope from the mountain front. These morphometric data reveal distinct patterns of geomorphic change across the Andean foreland basin, which we summarize in Table 1.

Table 1 Summary of Geomorphic and Topographic Expressions of Lithospheric Flexure

Distance downslope from mountain front	0–100 km	100–200 km	200–500 km	500–700 km
River planform	Braided	Transitional	Meandering	Meandering
River vertical position	Superelevated ∼0.8 m above floodplain	Transitional	Incised ∼4.9 m below floodplain	Incised ∼2.6 m below floodplain
River slope	Decreasing from 5 × 10−4 to 3 × 10−4	Decreasing from 3 × 10−4 to 2 × 10−4	∼1 to 2 × 10−4	∼1 to 2 × 10−4
River width	1–3 km	Narrowing	270 m	170 m
Chanel belt width	5 km	5–10 km	2 km	1.2 km
Paleochannel evidence for past avulsions	Abundant paleochannels which branch out from the modern Bermejo channel		Paleochannels both branching out from and running parallel to the modern channel	Paleochannels running parallel to, but not stemming from the modern channel
Lithospheric flexure	Negative	Transitional	Positive	Near zero
Interpretation	Foredeep subsidence creates accommodation space, triggering sediment deposition, channel superelevation, and avulsions	Transitional regime between foredeep subsidence and forebulge uplift	Forebulge uplift forces fluvial incision, channel narrowing, and a transition from braided to meandering. Paleochannels indicate past phases of an over-filled basin with river avulsions	Decreased amplitude of flexure reduces the amount of channel incision. Emergence of springs and standing water suggests a reduced depth to groundwater table

Near the megafan apex (∼0–150 km downslope from the mountain front) the river channel occupies an up to 3 km-wide braid plain positioned within a wider channel belt (Figure 2b). The channel belt is flanked on both sides by a 3–5 km wide zone of nearly contiguous positive relief, forming the levees that contain the active channel (Figure 3). Here, the river is elevated ∼0.8 m on average above the surrounding flood basin beyond the levees, suggesting that the channel is aggrading and building toward a future avulsion (Figures 2c and 3). Channel migration rates within the channel belt are high, sustained at 40–80 m/yr for decades. Channel slope is 2–5 × 10⁻⁴ through the braided reach, as the river descends into the Chaco Plain (Figure 2a). Examination of the residual relief map (Figure 3) and satellite imagery (Figure 1c) show abundant paleochannels and paleo-levees, which fan out from the current main stem Rio Bermejo, providing additional evidence for frequent avulsions between 0 and 200 km downslope of the mountain front.

Generate image description

Figure 2. Rio Bermejo channel morphometrics. (a) Mainstem channel slope (20 km running average). Gray circles show locations of water samples collected from the mainstem Rio Bermejo (see Figure 1b for locations of tributary, groundwater, and precipitation samples). (b) Active channel belt width (black) and active channel width (gray). (c) Elevation of floodplain surface relative to active channel (positive numbers represent the depth of channel incision into the floodplain surface, negative numbers represent channel super-elevation above the floodplain). (d) Rio Bermejo mainstem channel lateral migration rate. (e) Magnitude of lithospheric flexure in the Rio Bermejo study area, as modeled by Chase et al. (2009).

Generate image description

Figure 3. Residual topography of the Rio Bermejo megafan derived from TanDEM-X data. (a) Elevation of the floodplain calculated relative to the left (northeast) bank of the Rio Bermejo over swaths extending 25–40 km from the channel centerline. Positive numbers represent mainstem channel incision into the floodplain surface and negative numbers represent channel super-elevation above the floodplain. Red boxes show the aerial extents of the detailed views shown in panels (b) (foredeep), (c) (forebulge), and (d) (backbulge). Pink lines in panels (b)–(d) show the locations of cross-section profiles used to generate the relative elevation profiles in panels (e) foredeep, (f) forebulge, and (g) backbulge. Gray lines in panels (e)–(g) show the relative elevation along individual cross-sections and the black line represents the mean of those cross-sections.

As previously observed (McGlue et al., 2016; Repasch et al., 2020), between ∼100 and 200 km downslope of the mountain front, the Rio Bermejo exhibits a change in planform channel morphology from braided to meandering, a roughly three-fold decrease in channel width, and a sharp break in channel slope (Figure 2). Over this same area, the Rio Bermejo transitions from a superelevated position above the surrounding flood basin to being incised up to 13 m below the alluvial surface (Figures 2c and 3).

Between ∼200–500 km downslope of the mountain front, the channel slope (∼1–2 × 10⁻⁴), channel width (271 ± 111 m), channel belt width (1.98 ± 0.98 km), relative incision (4.85 ± 0.79 m), and planform morphology of the Rio Bermejo stay approximately constant. In this foreland segment, the Rio Bermejo is consistently incised (an average of 4.85 ± 0.79 m) into older fluvial deposits. The Rio Bermejito (Figure 1b), which occupies a paleo-channel belt from which the Rio Bermejo avulsed in 1870 (Page, 1889), emerges in this segment of the foreland basin. The residual elevation map shows sinuous paleochannels up to ∼50 km away from the modern Bermejo channel. These paleochannels are slightly inset within the elevated alluvial surface surrounding the modern channel belt (Figures 3a and 3c), documenting temporally sustained river incision across this foreland segment and periods when the river occupied different positions within the megafan. These paleochannels have different upstream entry points, suggesting that upstream avulsions can determine the location of drainage across this higher segment of the foreland (Figure 3a). At multiple locations within the raised alluvial surface, paleochannel features with lengths exceeding 50 km are superimposed on the megafan, splaying out in a fan shape. This indicates that there are phases when the river is not stable within an incised channel belt, but is instead free to migrate across this segment of the foreland. We suggest that this most likely occurs shortly after the main river has switched to a new location within the uppermost foreland basin. Under these conditions, an extensive wetland may form and persist until the river can incise through the elevated portion of the foreland basin. The Rio Pilcomayo, traversing the east Andean foreland immediately to the north of the Rio Bermejo may currently be in this state (see Section 5.3).

Starting at ∼500 km downslope of the mountain front, we observe decreases in channel width (173 ± 42.7 m), channel belt width (1.2 ± 0.6 km) and meander migration rate (3 ± 3 m/yr) (Figure 2b). The river channel remains incised but is on average less incised than in the forebulge, and the incision depth becomes more variable (2.61 ± 1.41 m) (Figure 2c). Low and irregular levees flank the channel. They become more pronounced in the lowermost part of the megafan, where they contain the active channel in a super-elevated position. At the transition into this foreland segment, paleochannels coming off the elevated alluvial surface upstream have fanned out, leaving sinuous sedimentary ridges elevated up to a few meters above the surrounding terrain. Further downslope these ridges are more distinct, and higher, running sub-parallel to the current channel belt. Another salient feature of this foreland segment is the emergence of spring-fed channels near −59.4° longitude (Figure 3; Hartley et al., 2013; McGlue et al., 2016).

These findings (summarized in Table 1) are consistent with our hypothesis of channel response to lithospheric flexure in a filled foreland basin. Within the foredeep, subsidence creates accommodation space and encourages sediment deposition within the Rio Bermejo channel belt. As observed in other depositional landscapes (e.g., Ganti et al., 2014; Hajek & Edmonds, 2014; Jerolmack & Mohrig, 2007), we suggest that subsidence-induced sediment deposition can cause the channel to become locally superelevated above the surrounding topography, creating an unstable situation where large floods can cause a river to avulse to a lower position. Repeated avulsion from a node near the mountain front has caused fanning of the Rio Bermejo within the foredeep. At times, the presence of a topographic barrier downstream may have caused ponding of river water and sediment, and wetland formation within this foreland segment.

Changes in channel planform, relative incision/superelevation, and avulsion frequency of the Rio Bermejo occur in concert with the inferred switch from negative to positive lithospheric flexure at ∼140 km downslope from the mountain front (Figure 2; Chase et al., 2009). This suggests a link between the foredeep-forebulge transition and river morphodynamics. Rock uplift can trigger channel narrowing and incision (e.g., Amos & Burbank, 2007; Duvall et al., 2004; Schumm, 1985), and the observation that the onset of incision is coincident with the inferred change from negative to positive lithospheric flexure suggests a process link. Uplift-induced incision can provide a mechanistic basis for the transition from a braided river with a wide channel belt to a single-thread meandering river confined to a narrow channel belt. This interpretation of uplift-induced incision driving the braided-to-meander transition is compatible with classic stability analyses where alluvial rivers transition from braided to meandering following increases in the depth to width ratio and/or decreases in the ratio of slope to Froude number (Parker, 1976). In the case of the Bermejo, we posit that forebulge uplift drives incision and associated decreases in channel width and increases in flow depth. This causes an increase in the depth to width ratio, which promotes the planform change from braided to meandering (Text S1, Figure S1 in Supporting Information S1). Furthermore, while naturally formed levees constrain the active channel in the braided section of the Rio Bermejo (e.g., Figure 3b), incision in the forebulge cuts into older floodplain sediment up to ∼20 ky in age (Scheingross et al., 2021). This aged, cohesive material likely has significantly higher bank shear strength than the relatively young levee deposits upstream, which can further facilitate the onset of meandering (e.g., Braudrick et al., 2009; Ielpi et al., 2022). In both of these cases, it is uplift-induced incision that is key to initiating the suite of geomorphic processes that cause a transition from braided to meandering. A high sediment load, combined with the significant reduction of accommodation space may drive the high meander migration rates through the forebulge (e.g., Constantine et al., 2014). The subtle but sustained incision of the Rio Bermejo across the forebulge, which otherwise has an alluvial cover, indicates that this foreland basin segment is teetering between being under-filled and over-filled. In this state, the river has managed to find multiple pathways across the forebulge, leaving evidence of shallow incision of older meandering channel belts, as well as phases of channel instability and splaying.

The forebulge to backbulge transition is marked by an abrupt change of channel characteristics, commensurate with the decreasing amplitude of flexure with increasing distance from the mountain load. At this transition, springs emerge at the surface (Hartley et al., 2013) and standing water becomes more prevalent, suggesting a decreased depth to the water table that is driven by minor subsidence. The abrupt decrease in channel migration rates at ∼500 km downstream may be driven by the change from uplift and incision to subsidence and increased accommodation space. Decreased lateral channel mobility may also be driven by reduced water and sediment discharge, as groundwater emerges at the surface and is distributed in small channels and lakes across the backbulge. However, this mechanism is more likely to cause a gradual reduction in migration rates moving downstream, rather than the abrupt change we observe. Nonetheless, this downstream pattern of channel migration rates is consistent with the downstream patterns of channel width, relative elevation, and planform change. Overall, these patterns align with predicted changes in lithospheric flexure, as modeled by Chase et al. (2009) (Figure 2), suggesting that subtle lithospheric dynamics can have profound impacts on fluvial geomorphology.

Hydrologic and Geochemical Signatures of Lithospheric Flexure

Next, we test our hypothesis that the geomorphic changes associated with lithospheric flexure give rise to changes in the hydrology and river water geochemistry of the Rio Bermejo. We posit that the superelevation of the Rio Bermejo in the foredeep causes seepage of river water into groundwater, thereby reducing surface water discharge with distance downstream. Conversely, incision of the river can allow subsurface flow from surrounding topography into the channel, possibly increasing the river dissolved load concentrations. Testing our hypothesis requires discharge records both upstream and downstream of the foredeep and forebulge, estimation of water evaporation along the Rio Bermejo, and measurements of the dissolved load geochemistry. This is perhaps the first test of how foreland flexure can influence hydrology and geochemistry in the Rio Bermejo or elsewhere.

Methods Water Discharge Estimates

We use daily water discharge data collected from 1970 to 2016 by the Argentina National System for Hydrologic Information (https://snih.hidricosargentina.gob.ar) to assess downstream trends in water discharge along the Rio Bermejo. Daily mean discharge data from a gauging station in El Colorado, Argentina (Figure 1b), ∼580 km downstream of the mountain front and ∼100 km up-valley of the confluence with the Rio Paraguay, is our measure of water discharge in the downstream portion of the Rio Bermejo, Q_down. Estimating the water discharge that enters the foredeep at the mountain front is more difficult, because there is no gaging station on the mainstem Rio Bermejo immediately downstream of the Rio Bermejo-Rio San Francisco confluence. We estimate this upstream water discharge using two different approaches. First, we estimate a minimum upstream discharge by a simple combination of the two upstream gauging stations (Rio San Francisco (RSF) at Caimancito and Rio Bermejo at Pozo Sarmiento (Figure 1b)), [Image Omitted. See PDF]where Q_up is the estimated discharge from the Rio Bermejo headwaters into the foredeep, Q_PS is the Rio Bermejo discharge at Pozo Sarmiento, and Q_SF is the RSF discharge at Caimancito. Second, we apply a drainage area correction to the combination of these measurements resulting in a larger, but likely more realistic, discharge estimate, [Image Omitted. See PDF]where A_PS = 24,128 km² is the contributing drainage area at Pozo Sarmiento, A_SF = 21,917 km² is the RSF drainage area at Caimancito, and A_BER@CNF = 25,366 km² and A_SF@CNF = 26,508 km² are the drainage areas for the Rio Bermejo and RSF immediately upstream of their confluence, respectively. This drainage area correction assumes no evaporative or groundwater losses between the gauging stations and the Rio Bermejo-Rio San Francisco confluence and yields an upstream discharge that is ∼8% greater than the combined Q_PS and Q_SF without drainage area correction (Equation 1). We calculated long term daily mean discharges for each day of the calendar year using daily discharge data recorded at the respective gauging stations over the duration of the record. We removed years in which more than 20% of the daily discharge records were not included in our measurement to avoid underestimation of the total annual discharge. As a result, 38% of the 47-year discharge record is excluded from our discharge analysis.

Water Sampling and Analytical Measurements

We collected water samples from the Rio Bermejo and its headwater tributaries, adjacent groundwater wells, spring-fed channels, and precipitation events that occurred during our field campaigns. We analyzed samples for major elements to detect any spatial changes in solute concentrations driven by groundwater influx into the channel and stable isotope composition to evaluate the potential influence of evaporation on downstream water loss (Text S2 in Supporting Information S2). We collected water samples during the rainy season at multiple stations along the length of the Rio Bermejo in three years (March 2016, March 2017, March 2020) and during the transitional (May 2015) and dry (November 2019) season in one year each (Figure 1). Mainstem river water samples were complemented by collection of groundwater at five accessible groundwater wells at locations upstream, within, and downstream of the forebulge, in addition to samples from select headwater tributaries and spring-fed channels (Table S3 in Supporting Information S2, Figure 1). Rio Bermejo water samples were collected either with a high-density polyethylene (HDPE) bottle, an aluminum bucket, or an 8-L van Dorn-type sampling bottle, and groundwater samples were taken directly from water taps on wells. Sampled groundwater wells are used for local drinking water, and it was not possible to measure the depth to the water table at these locations. After collection, we filtered water samples through a 0.22 μm polyethersulphone (PES) filter and pre-rinsed sample vessels with filtered sample water prior to filling sample bottles. Samples for water stable isotope and major element analyses were collected and stored in glass vials and HDPE bottles, respectively. Mainstem water samples in 2017 were collected from multiple depths below the river water surface, while water samples in all other years were collected from the water surface. Samples for cation measurements were immediately acidified to pH < 2 with nitric acid. All samples were refrigerated upon return from the field.

Cation concentrations were measured with a Varian 720 inductively coupled plasma optical emission spectrometer at the GFZ Helmholtz Laboratory for the Geochemistry of the Earth Surface (HELGES), using SLRS-5 (Saint-Laurent River Surface, National Research Council-Conseil National de Recherches Canada) and USGS M212 and USGS T187 as external standards. We corrected for instrument drift by measuring an internal standard (GFZ-RW1) every 10 samples and we determined measurement uncertainty using calibration curve uncertainty. Anion concentrations were measured with a Dionex ICS1100 Ion Chromatograph, using USGS standards M206 and M212 as external standards for quality control, with uncertainty determined from triplicate analysis.

We corrected cation concentrations for cyclic salt inputs following Bickle et al. (2005). We set the concentration of Cl in rainwater to Cl_rain = 19.7 μmol/l using the lowest Cl concentration measured in the headwater tributaries of the Bermejo (Table S3 in Supporting Information S2). We estimate Mg_rain, Na_rain, K_rain, and Ca_rain as X_rain = (X/Cl)_rain * Cl_rain, where X is the solute of interest and the subscript ‘rain’ denotes rainwater. Data on solute concentrations in precipitation are rare in our study region and our precipitation samples were not large enough to measure solutes directly, therefore we estimate cation to Cl ratios by compiling precipitation data from the Amazon Basin (Forti et al., 2000; Moquet et al., 2011) and Uruguay (Zunckel et al., 2003) to give average molar ratios of Mg/Cl, Na/Cl, K/Cl, and Ca/Cl for rainwater of 0.2, 1.05, 0.34, and 0.97, respectively. The above procedure results in cyclic corrections that are typically less than 10% of measured concentrations. After correcting for cyclic inputs, we estimate $${\text{HCO}}_{3}^{-}$$ using a charge balance, and sum all the measured anions and cations (excluding Si) as a measure of total dissolved solids (TDS) concentration.

Results and Interpretations Downstream Decrease in Water Discharge

Water discharge in the Rio Bermejo is highly seasonal; however, apart from April through July when river discharge is rapidly decreasing, the upstream water discharge entering the foredeep (calculated by Equation 1 or Equation 2) is always greater than the downstream water discharge measured at El Colorado when averaging over the 46-year record (Figure 4). The long-term mean and standard deviation of upstream river discharge over the water year is 523 ± 525 m³/s using Equation 1 (or 567 ± 568 m³/s with the drainage area correction in Equation 2, Figure S2 in Supporting Information S1) and 432 ± 382 m³/s at the downstream gauging station, representing a 17%–24% average discharge loss along the 1,086 km river flow path between gauging stations.

Generate image description

Figure 4. Comparison of Rio Bermejo long-term daily mean water discharge from gauging stations at the mountain front (Qup, Equation 1) and at 580 km downstream (Qdown, measured at El Colorado). Long-term daily mean discharge was calculated as the mean discharge recorded for each calendar day over the 46-year gauging record from 1970 to 2016, excluding years with missing data (see text for details). Ratios of Qup/Qdown > 1 represent a mean water loss between the mountain front and 580 km downstream, whereas Qup/Qdown [less than] 1 represents a mean gain of water discharge between these points. This threshold is marked by the horizontal gray line at Qup/Qdown = 1.

The reach of the Rio Bermejo traversing the foreland basin has no notable tributaries or distributaries, with two exceptions. The Rio Bermejito paleochannel adds <2% of the main stem Rio Bermejo water discharge (Orfeo, 2006), and a diversion canal at ∼265 km downslope, established in 2013, removes <6% of river water based on its intake rate of 30 m³/s at high flows (Rohrmann, 2021). Therefore, the only plausible mechanisms for the 17%–24% downstream decrease in water discharge are transfer of river water to groundwater, either directly from the channel or via inundation of the flood basin, or evaporation. Analysis of stable water isotopes (δ¹⁸O and δ²H) measured in groundwater and river water throughout the catchment suggests that evaporation is not a major term in the water balance of the Rio Bermejo, and especially not during the rainy season (Text S2, Figure S3 in Supporting Information S1). Both the hydrologic and water isotope data are consistent with a net loss of discharge due to seepage out of the channel. We propose that the greatest potential for this exists where the Rio Bermejo traverses the foredeep, due to its superelevated position above the flood basin. We expect that groundwater flows into the river channel where the river cuts into the elevated surface in the forebulge area. Based on the observed upstream and downstream water fluxes (Figure 4), conservation of mass implies that this groundwater influx in the forebulge must be (on average) less than the loss of river water to groundwater in the foredeep.

The most thorough test of this hypothesis would be to measure the hydraulic gradient across the flood basin via groundwater monitoring wells. Unfortunately, the groundwater wells in the Bermejo floodplain were designed for local drinking water and are not monitored for hydraulic head. However, our observations are consistent with geomorphology and corresponding hydrogeologic theory. Because the channel belt in the foredeep is a local topographic high, the groundwater table must either rise to meet the river bed, requiring the hydraulic gradient to be flowing away from the channel, or the groundwater table must be lower than the river bed, such that water flows from the river channel to the groundwater reservoir, but not vice versa (sensu Winter et al., 1999; Woessner, 2000). In both cases surface water infiltrates into the groundwater reservoir, resulting in a losing stream. In contrast, in the forebulge wetlands are common (Benzaquen et al., 2017), indicating a relatively high groundwater table, and the channel is a topographic low. Together these observations in the forebulge suggest that the hydraulic gradient slopes toward the channel and groundwater flows into the Bermejo. Furthermore, moving from the foredeep to the forebulge the channel gradient and width-to-depth ratio decrease, while sinuosity and incision increase, all of which are consistent with a general transition from predominantly channel-parallel groundwater flow in the foredeep to channel-perpendicular flow in the forebulge (Larkin & Sharp, 1992). In this case, surface water lost to groundwater in the foredeep is unlikely to return to the channel until it reaches the forebulge. We test this idea by evaluating changes in river water solute concentrations across this transition.

Changes in River Solute Concentrations Across the Forebulge

If there exists an influx of groundwater into the Rio Bermejo in the forebulge, we expect the river solute load to increase as the river crosses the forebulge. Total dissolved solids (TDS) concentrations measured in Rio Bermejo water in the wet (March 2016, March 2017, and March 2020) and transitional (May 2015) seasons range from 180 to 290 ppm, whereas the dry season (November 2019) has systematically higher TDS concentrations ranging from ∼360 to 670 ppm (Figure 5a, Table S3 in Supporting Information S2). We normalized TDS concentrations to TDS measured 116 km downslope of the mountain front to make it easier to evaluate the pattern of downstream changes in solute loads between years. In each year, TDS concentrations during the wet and transitional seasons showed similar patterns of low concentrations in the foredeep, a progressive increase to maximum concentrations near the forebulge/backbulge transition (∼500 km downslope of the mountain front), and approximately steady concentrations at distances further downstream (Figure 5b). Samples collected in the rainy and transitional season showed a 20%–30% relative increase in TDS across the forebulge, whereas there was a ∼60% increase in TDS concentration over the same distance for the dry season samples (Figure 5b). In all cases, this 20%–60% increase significantly exceeds the analytical uncertainty of individual data points.

Generate image description

Figure 5. (a) Total dissolved solids (TDS) concentrations for Rio Bermejo river water and groundwater wells. Similar TDS values for samples collected at the same location over multiple field campaigns results in symbols plotting on top of one another, thus hiding some symbols. Individual sample values are reported in Table S3 in Supporting Information S2. Error bars are smaller than symbols where not shown. Groundwater samples with upward pointing arrows are minimum estimates for samples where Na concentration exceeded the calibrated range for a given analytical run. (b) Rio Bermejo river water TDS concentration normalized by the concentration measured at 116 km downslope from the mountain front. Error bars represent propagated error of analytical measurements. Gray shading in both panels denotes the approximate forebulge extent.

The pattern of TDS concentration as a function of distance along the Rio Bermejo also exists in the individual concentrations of most major elements (Figures S4 and S5 in Supporting Information S1, Table S3 in Supporting Information S2). Ca, Mg, Na, F, and SO₄ concentrations all increased by a factor ∼1.2 to ∼2.0 across the forebulge, before approximately stabilizing in the backbulge (Figures S4 and S5 in Supporting Information S1). Si concentrations are an exception to this general pattern. Si concentrations decreased progressively along the length of the Rio Bermejo from ∼7 ppm near the mountain front to ∼4.5 ppm near the outlet (Figure S4f in Supporting Information S1), perhaps due to biological uptake by benthic diatoms or secondary mineral formation, the latter of which has been observed in the Rio Bermejo (McGlue et al., 2016).

These changes in downstream solute concentrations cannot be fully explained by evaporative water loss. We observe a consistent trend of increasing solute concentrations across the forebulge in all 5 years of our sampling campaign (Figure 5), including three years where water stable isotope analyses indicate no significant evaporative losses (March 2016, 2017, and 2020, Text S2 in Supporting Information S2). Evaporation from the forebulge and backbulge (Figure S3 in Supporting Information S1) may partially explain the larger increase in normalized TDS concentration in November 2019 (Figure 5), but this increase could also be due to the proportionally low water discharge in the dry season, allowing groundwater infiltration to represent a larger fraction of total discharge. Evaporative losses in March 2015 occur upstream of the forebulge (Figure S3 in Supporting Information S1), suggesting that increased TDS concentrations within the forebulge in that year are not a product of evaporation. A concentration discharge analysis (Text S2, Figure S6 in Supporting Information S1) provides further evidence against evaporation as a driver of our observed downstream changes in solute concentrations.

Lacking evidence for a consistent major role of evaporation in the Rio Bermejo, we suggest the most likely cause of downstream increases in solute concentrations is groundwater flow into the river channel within the forebulge region. In this region, the Rio Bermejo is systematically incised into elevated, older river sediments (Scheingross et al., 2021). Incision to a depth below the water table allows for increased groundwater flow into the river. While groundwater table data are not accessible for our study area, the existence of wetlands throughout the forebulge (e.g., Benzaquen et al., 2017) indicates a high groundwater table in this region, such that any amount of incision is likely to increase groundwater flow into the river channel. Groundwater samples in this region, and elsewhere in the megafan, consistently have TDS concentrations ∼2–10 times higher than those in Rio Bermejo mainstem water (Figure 5a) and high concentrations of individual major elements (Table S3 in Supporting Information S2), thereby providing a source for input of concentrated solute loads. This suggests that lithospheric flexure may play a critical role in regulating the dissolved load geochemistry of the Rio Bermejo.

Synthesis and Implications

Our analysis of the Rio Bermejo basin highlights how feedbacks between tectonics and surface processes can have profound impacts on low-relief landscapes. Our data show how tectonic loading and associated lithospheric flexure can drive changes in river planform morphology, channel morphodynamics, groundwater-surface water interactions, and the geochemical composition of river water. Even in a low-relief landscape where the surficial expressions of lithospheric dynamics are imperceptible in the field, the patterns of geomorphic and geochemical change align with model predictions of changes between subsidence and uplift across the foreland basin. We show how these changes can help delineate the extent of the foredeep, forebulge, and backbulge (allowing independent estimates the mechanical properties of the lithosphere) and discuss broader implications for landscape evolution in the Andean foreland basin and interpretation of weathering signals in foreland basins.

Geomorphic and Geochemical Delineation of the Foredeep, Forebulge, and Backbulge

Based on our analysis, the foredeep is a zone of subsidence extending from the mountain front to ∼150 km downslope. Here, subsidence creates accommodation space, forcing the river to deposit sediment, aggrade its bed (to the point of super-elevation relative to the surrounding surface) and periodically avulses when sediment deposition in the channel belt exceeds a critical height (e.g., Jerolmack & Mohrig, 2007). We suggest that the foredeep is the most likely location for the observed water exfiltration (Figure 4) from the channel due to the superelevated position of the channel belt (Figures 3 and 6). Water lost from the river to groundwater is likely to flow downslope away from the mountain front and subsequently increase in solute concentration due to weathering in the floodplain. Groundwater enriched in dissolved solutes may re-enter the incised mainstem channel or re-emerge in springs beyond the forebulge (Figure 6).

Generate image description

Figure 6. (a) Schematic conceptual model showing how lithospheric flexure influences the geomorphology and hydrogeology of the Rio Bermejo fluvial system in the Andean foreland basin. Schematic channel cross sections show the groundwater response to (b) foredeep subsidence causing super-elevation of the channel above the floodplain, (c) forebulge uplift causing channel incision below the surrounding floodplain, and (d) backbulge subsidence resulting in the development of a spring-fed distributary system. Dark blue arrows depict the relative direction of groundwater flow.

The forebulge is a zone of upward lithospheric flexure spanning ∼150–500 km downstream from the mountain front, where the river incises up to 13 m below the surrounding surface to maintain grade. Uplift-induced incision into the elevated paleo-floodplain surface and the lack of accommodation space in the forebulge drive the development of a narrow, meandering channel. Channel migration rates are high through the forebulge, averaging ∼14 m/yr, with some meander bends migrating at rates up to 50 m/yr. This rapid channel migration may be due to the combination of a large river sediment load and the absence of accommodation space in a channel belt that is incising to keep up with the rate of uplift. Incision of the channel below the groundwater table (Figure 6b) promotes an influx of groundwater that increases river solute concentrations by 20%–30% in the rainy season; however, the lack of discharge data in this reach of the river makes it difficult to estimate the volumetric influx of groundwater. The increased river solute concentrations mark the forebulge as a geochemically important segment of the foreland, where weathering can impact downstream water chemistry due to the fluvial export of groundwater.

The backbulge extends from ∼500 km downstream of the mountain front to a point beyond the confluence of the Rio Bermejo with the Rio Paraguay. On the downslope side of the forebulge, spring-fed channels signal groundwater emergence, as previously observed (Hartley et al., 2013). We posit that the emergence of springs is linked to the transition from relative uplift to subsidence, reducing the depth to the groundwater table (Figure 6c). Subsidence in the backbulge also drives ponding across the area (McGlue et al., 2016); however, the mainstem river channel and adjacent spring-fed channels remain incised through the backbulge (Figure 6c), likely a result of the lower base level of the Rio Paraguay. This is currently reflected by a ∼10 m elevation drop over the lower 15 km of the Rio Bermejo before it joins the Rio Paraguay. River water chemistry remains steady through the backbulge, suggesting that either groundwater and surface water have similar solute concentrations, or that any exchanges between the groundwater and surface water are small enough to not significantly affect the river water chemistry.

Estimates of Lithospheric Mechanical Properties

In regions where direct geophysical measurements are not available, model predictions of lithospheric dynamics are uncertain, primarily due to difficulty in measuring the elastic thickness and rigidity of the underlying plate (e.g., Tesauro et al., 2012). Our observations from the Rio Bermejo allow us to delineate the extent of the foredeep, forebulge and backbulge. These delineations can be used to make independent estimates of the mechanical properties of the underlying lithosphere (Text S3 in Supporting Information S2) and thus provide a comparison against model predictions.

Following established theory (Turcotte & Schubert, 2002; Walcott, 1970), we find that the flexural rigidity of the underlying lithosphere is ∼3.77 × 10²³ N and we calculate the elastic plate thickness to be 39.3 km for the east Andean foreland basin at 16°S–18°S (Text S3 in Supporting Information S2), which is within the range of estimated elastic thicknesses of 20–90 km for foreland basins globally (DeCelles, 2012; Jordan, 1981; Lyon-Caen & Molnar, 1985; Roddaz et al., 2005; Watts, 2001) and similar to predictions for the East Andean Foreland of an ∼50 km elastic thickness (Chase et al., 2009). While there are uncertainties in this approach, our analysis shows that geomorphic measurements can complement geophysical constraints to improve our estimates of lithosphere mechanical properties.

Regional Implications

The effects of lithospheric flexure are not limited to the Rio Bermejo, as they exist regionally along the Andean mountain front. For example, the Rio Pilcomayo lies immediately north of the Rio Bermejo and flows eastward from the Andean mountain front across the foreland basin. A noteworthy feature of the Pilcomayo is its discontinuous mainstem channel. Like the Rio Bermejo, the upstream segment of the Pilcomayo is braided, superelevated, and avulses periodically (Ritter, 1977) before transitioning into a single-thread, narrow meandering channel ∼200 km downslope from the mountain front. In contrast to the Bermejo, the Pilcomayo then spills out into a system of wetlands, the spatial location of which is coincident with the region of predicted forebulge uplift. Moving further east, several spring-fed channels emerge at ∼500 km downslope from the mountain front, draining the Pilcomayo system to the Rio Paraguay (Hartley et al., 2013). We propose that these geomorphic and hydrologic patterns in the Rio Pilcomayo are regulated by lithospheric flexure.

While both the upstream segments of the Bermejo and Pilcomayo systems display similar channel geometry and autocyclic avulsion in the foredeep, the Pilcomayo does not currently incise across the region of forebulge uplift. One hypothesis for this difference is that the Bermejo and Pilcomayo are teetering between being overfilled and underfilled basins, and thus their ability to incise through the forebulge may depend on time since the last major river avulsion. For example, after a foredeep avulsion, the new river path will likely transverse the forebulge at a different location relative to the previous channel path. Even a small topographic high in the uplifting forebulge (which should exist if the basin is even slightly underfilled) can cause water and sediment to pool against the upstream side of the forebulge, progressively filling the foredeep basin until water spills over the topographic barrier. As water spills across the forebulge, shallow but extensive wetlands may develop on this topographic high because the channel is not yet incised there. With sustained river discharge, the river begins incising into the forebulge in the now overfilled system. Assuming a steady rate of forebulge uplift across the foreland basin, this is an autocyclic phenomenon where channels avulse in the foredeep, aggrade against the upstream side of the forebulge until the point of spilling over, incise through the forebulge and eventually avulse again in the foredeep. The rate and duration of each cycle should be determined by the avulsion frequency, which is set by discharge and sediment supply of the river system (Jerolmack & Mohrig, 2007). We posit that the Pilcomayo may be responding to a recent avulsion and is thus in the aggradation stage prior to spillover, whereas the Bermejo is in a phase of forebulge channel entrenchment and building toward a new foredeep avulsion. If this hypothesis is true, wetlands ponded against the forebulge in the Pilcomayo system are ephemeral.

Implications for Chemical Weathering in Lowland Floodplain Deposits

Our measurements of increased river solute concentrations crossing the forebulge not only serve as a proxy for groundwater addition to the river in areas of upward lithospheric flexure, but may also reveal the weathering of foreland river deposits. Chemical weathering plays an important role in regulating Earth's climate, for example, via CO₂ drawdown during silicate mineral dissolution and subsequent precipitation of carbonate minerals in the oceans (Berner et al., 1983; Gaillardet et al., 1999; Kump et al., 2000; Walker et al., 1981). Erosion of the Andes generates large amounts of sediment, composed mostly of granular silicates. More labile minority constituents such as carbonates and sulfides may weather quickly in the mountain headwaters and do not reach the mountain front in significant quantities. The remaining sediment is delivered to the foreland basin where transient storage, orders of magnitude longer than in mountain hillslopes (e.g., Nakamura et al., 1995; Repasch et al., 2020), affords time for the slow reactions involving carbonic acid and silicate minerals (cf. Bufe et al., 2021). As a result, a significant amount of silicate weathering can occur in lowland basins (e.g., Bouchez et al., 2012; Lupker et al., 2012), especially where the groundwater level is far beneath the surface. Moreover, slow and steady percolation of water in the floodplain weathering zone (Maher & Chamberlain, 2014) and seasonal water table fluctuations may promote dissolution and precipitation reactions (Ramos et al., 2022).

In the Rio Bermejo, we observe a significant downstream increase in river dissolved solute concentrations, collocated with the flexural forebulge. This could indicate high weathering intensities in this part of the foreland basin, but it could also reflect the inflow of groundwater-borne weathering products from locations farther upslope in the foreland. Freshly eroded sediment aggrades in the foredeep basin, an area that is prone to seasonal inundation, making this a likely hotspot of silicate weathering. Our data indicate that the net transfer of surface water to groundwater in the foredeep prevents any foredeep weathering from being represented in the Rio Bermejo dissolved load until at least the distal end of the foredeep. At depth, foredeep groundwater is likely to flow parallel to the river, across the flexural strike and away from the mountain front. This may add to groundwater generated locally in the forebulge by rain infiltration and drive the increased river solute concentrations across the forebulge, where groundwater flow is diverted into the channel because of incision into the uplifted floodplain. We observed no further increase in dissolved solutes through the backbulge. This may be due to a lack of systematic convergence of groundwater into the river channel, but it could also suggest that sediment delivered to the distal part of the fan is chemically depleted and unlikely to continue weathering. These observations may indicate that the size of lowland basins does not determine the degree of chemical weathering and potential for CO₂ drawdown, but rather that lowland basin size sets the available accommodation space for sediment trapping and storage, and weathering may be modulated by flexural subsidence and uplift, and by regional climate that allows for seasonal water table fluctuations.

Conclusions

The Rio Bermejo fluvial system currently sits within a filled foreland basin with no obvious topographic expression of lithospheric uplift, yet here we show that aspects of river hydrology, topography, and geochemistry are sensitive recorders of subtle tectonic forcing of the landscape. The morphodynamics of the Rio Bermejo and changes in its water chemistry across the Andean foreland basin reveal that the residual topography and landscape segmentation created by lithospheric flexure have profound effects on the geomorphology and hydrology of the river system. Results of this study suggest that measurements of channel morphometrics and dissolved load chemistry may help elucidate the properties of lithospheric flexure in other low-relief landscapes. Lithospheric dynamics can also determine where sediment can be trapped within the landscape and the relative position of the groundwater table, two key factors that drive weathering and organic carbon cycling in lowland landscapes, with potential ecological consequences. Ultimately, lithospheric flexure shapes landscapes and water flow paths, which determine land use and river management options and govern human settlement of foreland basins.

Acknowledgments

This study is a salute to Philip A. Allen, who recognized the links between surface processes and tectonics in foreland basin records. We acknowledge funding from an Alexander von Humboldt Foundation Postdoctoral Fellowship (J.S.S.), the Helmholtz Association (N.H.), and the StRATEGy international research training group funding by the German Research Foundation (DFG) and the State of Brandenburg (D.S., M.R., S.D.). We thank R. Lopez, H. Hassenruck-Gudipati, F. Latosinski, and N. Golombek, and H. Wittman for field assistance, D. Postler for preliminary hydrologic analyses, and B. Purinton for initial meander migration rate analyses. D. Frick, O. Rach, J. Schuessler, B. Zimmermann, and C. Zorn assisted with geochemical sample preparation and measurements. We benefited from constructive reviews from Editor M. Bayani Cardenas, Dr. Ying Fan, and an anonymous reviewer.

Conflict of Interest

The authors declare no conflicts of interest relevant to this study.

Data Availability Statement

All data generated for this study have been archived in the PANGAEA data repository (Repasch et al., 2023) and are available for download at https://doi.org/10.1594/PANGAEA.959938.