Cosmogenic nuclide and solute flux data from central Cuban rivers emphasize the importance of both physical and chemical mass loss from tropical landscapes

We use 25 new measurements of in situ produced cosmogenic 26 Al and 10 Be in river sand, paired with estimates of dissolved load flux in river water, to characterize the processes and pace of landscape change in central Cuba. Long-term erosion rates inferred from 10 Be concentrations in quartz extracted from central Cuban river sand range from 3.4-189 Mg km -2 y -1 (mean = 59, median = 45). Dissolved 30 loads (10-176 Mg km -2 y -1 ; mean = 92, median = 97), calculated from stream solute concentrations and modelled runoff, exceed measured cosmogenic 10 Be-derived erosion rates in 18 of 23 basins. This disparity mandates that landscape-scale, in this environment, mass loss is not fully represented by the cosmogenic nuclide measurements. The 26 Al/ 10 Be ratios are lower than expected for steady-state exposure or erosion in 16 of 24 35 samples. Depressed 26 Al/ 10 Be ratios occur in many of the basins that have the greatest disparity between dissolved loads (high) and erosion rates inferred from cosmogenic nuclide concentrations (low). Depressed 26 Al/ 10 Be ratios are consistent with the presence of a deep, mixed, regolith layer providing extended storage times on slopes and/or burial and extended storage during fluvial transport. River water chemical analyses indicates many basins with lower 26 Al/ 10 Be ratios and high 10 Be ratios are 40 underlain at least in part by evaporitic rocks that rapidly dissolve. Our data show that when assessing mass loss in humid tropical landscapes, accounting for the contribution of rock dissolution at depth is particularly important. In such warm, wet climates, mineral dissolution can occur many meters below the surface, beyond the penetration depth of most cosmic rays, and thus the production of most cosmogenic nuclides. Our data suggest the importance of 45 estimating solute fluxes and measuring paired cosmogenic nuclides to understand better the processes and rates of mass transfer at a basin-scale.


Introduction
Cosmogenic nuclide concentrations of river sand have been used to quantify rates of landscape change (often termed erosion rates) since the 1990s (Brown et al., 1995a;Granger et al., 1996;Bierman and Steig, 1996;Portenga and Bierman, 2011;Codilean et al., 2018).Accurately establishing long-term rates of change provides an important context for understanding the effects of human activity on erosion (Reusser et al., 2015;Nearing et al., 2017), and for other common applications of cosmogenic nuclides at the basin-scale, such as quantifying the effect of tectonics (Scherler et al., 2014), climate (Marshall et al., 2017), and baselevel change (Reinhardt et al., 2007) on rates of landscape change over time.
10 Be-derived rates of landscape change at a drainage-basin scale are often implicitly assumed to reflect both physical and chemical mass loss, the sum of which is termed denudation (Regard et al., 2016).However, this assumption is only valid if all mass loss from the landscape occurs within the uppermost meter or two of Earth's surface, the penetration depth of the cosmic-ray neutrons responsible for producing most cosmogenic nuclides via spallation reactions (Bierman and Steig, 1996).Deeper mass loss by rock dissolution remains largely undetected by cosmogenic nuclide analysis.Failure to account for rock dissolution at depth and the export of mass as dissolved load below the spallationdominated nuclide production zone (~2 m) may bias cosmogenic nuclide-derived estimates of denudation (Small et al., 1999;Riebe et al., 2001a;Dixon et al., 2009a) on the low side.Incorrectly determined erosion rates can derail attempts to understand landscape evolution, soil production, and climate interaction with surface processes (Riebe et al., 2003).
Rock dissolution at depth is a major process in areas with significant groundwater-rock interactions; connecting denudation rates to landscape change requires consideration of this process.This includes any landscape where the physical removal of mass is slow, allowing for prolonged waterrock interactions, such as low-relief landscapes (Ollier, 1988).Some landscape characteristics facilitate or are the result of extensive water-rock interaction: thick saprolite (Dixon et al., 2009a), extensively jointed/fractured bedrock (Ollier, 1988), and readily soluble rocks, including carbonate (Pope, 2013) and evaporite deposits.Conditions in the humid tropics favor prolonged and extensive water rockinteraction and include the absence of recent glaciation (Modenesi-Gauttieri et al., 2011), the presence of active groundwater flow systems year-round (Ollier, 1988), and large amounts of precipitation.
Rock dissolution rates in the tropics can be among the highest globally (Pope, 2013); yet, global compilations of cosmogenic nuclide data from river sand suggest rates of landscape change in the tropics are slower than in most other climate zones (Portenga and Bierman, 2011).This dichotomy is consistent with cosmogenic rates significantly underestimating landscape denudation in areas where deep rock dissolution is ubiquitous.Only a few studies focused in the tropics compare nuclide-derived rates to measurements of dissolved load flux in streams (e.g.Salgado et al., 2006;Hinderer et al., 2013;Regard et al., 2016).As the use of cosmogenic nuclides to measure rates of landscape change in the tropics expands (e.g.Cherem et al., 2012;Barreto et al., 2013;Derrieux et al., 2014;Mandal et al., 2015;Sosa Gonzalez et al., 2016a;Jonell et al., 2017), considering the potential influence of rock dissolution at depths below the production of most cosmogenic nuclides becomes more important.
Here, we present measurements of in situ 26 Al and 10 Be in riverine quartz, along with estimates of dissolved loads, in humid, tropical central Cuba (Bierman et al., 2020).With these data, we explore the relationships between cosmogenic nuclide concentrations, dissolved load fluxes, and landscapescale parameters at a basin scale in a humid tropical location where mass is being lost from the landscape by multiple, different processes from a variety of rock types.We characterize the rates and processes by which the Cuban landscape is changing, and place these data in a global context.Our findings illustrate the importance of considering rock dissolution when using cosmogenic nuclides to assess rates of landscape change in areas with the potential for significant mass loss by solution at depth, and provide a geologic baseline for assessing the impact of human actions on the Cuban landscape.

Background
Terminology referring to mass loss from watersheds has been applied ambiguously in the past and can be confusing.Here, we refer to the tempo of landscape mass loss calculated from 26 Al and 10 Be concentrations as erosion rates; these rates include all processes (physical and chemical) removing mass within ~2 m of Earth's surface.We refer to rates of landscape mass loss inferred from measurements of stream water chemistry, convolved with estimates of annual run-off volumes, as rock dissolution rates.
We use the term denudation to refer to total mass loss from sampled catchments.All of these rates are expressed in terms of mass per time per area (Mg km -2 y -1 ), which can be converted to depth over time by assuming a rock density.

Quantifying basin mass loss with cosmogenic nuclides: approaches and limitations
Landscape-scale denudation occurs through both physical removal of mass (erosion) and chemical dissolution of minerals in rocks.Sediment, produced from eroding bedrock, travels downslope towards base level, whereas rock dissolution moves mass in solution from the landscape to rivers, and then to the ocean.Measurement of cosmogenic nuclides in river sediment can be used to infer the spatially averaged erosion rate of a drainage basin (Brown et al., 1995a;Granger et al., 1996;Bierman and Steig, 1996).In a basin that is steadily eroding, the concentration of cosmogenic nuclides in a sediment sample reflects the rate at which overlying mass at and near the surface was removed as the material was exhumed, through both physical mass loss and rock dissolution (Lal, 1991).Cosmogenic erosion rates are equivalent to denudation rates if, and only if, rock dissolution only occurs within a meter or two of the surface-the depth of penetration for neutrons which produce most cosmogenic nuclides.If rock dissolution occurs below the neutron penetration depth, erosion rates calculated from measured nuclide concentrations will underestimate denudation.
Measuring multiple cosmogenic nuclides with different half lives in the same sample can provide more information on the near-surface history of surface materials, such as soil mixing depth and residence time (Lal and Chen, 2005), as well as sediment storage within the watershed (Granger and Muzikar, 2001).The production ratio of 26 Al/ 10 Be at the surface at mid-and low-latitudes is constrained by measurements and nuclear physics (Nishiizumi et al., 1989;Balco et al., 2008).If sediment that has accumulated cosmogenic nuclides is buried such that production is diminished over > 10 5 y, the production ratio decreases because 26 Al decays more rapidly than 10 Be.Vertical soil mixing intermittently buries sediment grains, suppressing the 26 Al/ 10 Be ratio in sediment shed from the landscape surface (Makhubela et al., 2019).
Paired cosmogenic isotope concentrations are visualized using a two-isotope diagram; the y-axis is the 26 Al/ 10 Be ratio and the x-axis is the concentration of 10 Be with normalization based on the production rate of nuclides at the sample site (Klein et al., 1986;Granger, 2006).Sediment samples that have experienced constant exposure with no erosion, or constant exposure under steady-state erosion, will plot within an enclosed region along the top of the diagram; samples that have experienced more complex exposure histories, including burial during or after cosmic-ray exposure, will plot below this region.Such complex histories could include development of a vertically mixed surface layer (Bierman, 1999;Lal and Chen, 2005) as well as extended burial during transport down slopes and in and along rivers.
Using measurements of cosmogenic nuclides to determine basin-averaged denudation rates requires the assumptions that mass loss from the basin is in steady-state, that the mineral used for cosmogenic nuclide measurements is uniformly distributed throughout the watershed, and that denudation occurs within the penetration depth of most cosmic rays, the upper several meters of Earth's surface (Bierman and Steig, 1996).The grain size fraction selected for cosmogenic nuclide analysis must also be representative of the grain size distribution of sediment being produced on slopes (Lukens et al., 2016) although in many landscapes cosmogenic nuclide concentrations do not vary by sediment grainsize.
Erosion rates calculated from cosmogenic nuclides can be inaccurate if these assumptions are violated.Rock dissolution can leave sediment enriched in resistant mineral phases, such as zircon, titanite, and quartz-the mineral in which 26 Al and 10 Be are most commonly measured (Riebe and Granger, 2013).Such enrichment produces underestimates of long-term denudation rates unless accounted for, because the enriched mineral will have a longer residence time relative to the surrounding regolith (Riebe et al., 2001a;Ferrier and Kirchner, 2008).Calculations of denudation rates from cosmogenic nuclide concentrations also rely on the assumption that mass loss is occurring primarily through surface lowering; however, some rock dissolution and, thus, some transfer of mass from rock to groundwater solutions occurs below the depth of most cosmogenic nuclide production (Fig. 1; Small et al., 1999;Dixon et al., 2009a;Riebe and Granger, 2013).In areas with significant rock dissolution at depth, denudation rates inferred from cosmogenic nuclides underestimate denudation because some mass loss occurs below the depth of most nuclide production.

Chemical weathering corrections to cosmogenically-determined mass loss rates
Although the importance of accounting for loss of mass by chemical weathering (rock dissolution) when calculating cosmogenic erosion rates has been recognized (Small et al., 1999;Riebe et al., 2001a;Dixon et al., 2009a;Riebe and Granger, 2013), few studies incorporate rock dissolution information or apply correction factors to cosmogenic nuclide-derived rates.In the tropics, some studies compare export rates from dissolved loads in streams to cosmogenically-derived erosion rates, but those studies have considered these two metrics of landscape change separately (Von Blackenburg et al., 2004;Salgado et al., 2006;Hinderer et al., 2013).Other studies use the measurement of insoluble elements in bedrock, saprolite, and soil to quantify quartz enrichment through the weathering process and calculate correction factors that account for the influence of rock dissolution at and near the surface (Small et al., 1999;Riebe et al., 2001a), at depth (Dixon et al., 2009b), or both (Riebe and Granger, 2013).
Of studies that do correct for the influence of chemical weathering when calculating cosmogenic nuclide-derived rates of erosion, Riebe and Granger (2013)'s chemical erosion factor (CEF) method, or earlier quartz enrichment factor method (Riebe et al., 2001a), are often used (Regard et al., 2016).
Calculating a CEF requires measurements of soil thickness and density, as well as determining the concentration of the mineral used in cosmogenic nuclide measurements (commonly quartz) and an insoluble element (commonly Zr) in numerous samples of soil, saprolite, and unweathered bedrock.The method is underpinned by the assumption that chemical mass loss is occurring exclusively in wellmixed regolith and deep saprolite (Riebe and Granger, 2013).Erosion rates calculated from cosmogenic nuclide measurements can be multiplied by the CEF to correct for the effects of chemical mass loss (Riebe and Granger, 2013).Chemical erosion factors reported in tropical environments include a CEF of 1.79 in Puerto Rico (Riebe and Granger, 2013) and 3.2 in Cameroon (Regard et al., 2016), demonstrating how significantly cosmogenic nuclide-derived estimates of erosion can underestimate total denudation rates by not accounting for the effects of deep rock dissolution.

Study area
Cuba is the largest Caribbean island (~110,000 km 2 ) and is situated along the boundary between the Caribbean and North American plates.Reflecting this complex tectonic setting, Cuban geology is varied and includes silicate, carbonate, and evaporite rocks (Pardo, 2009).Lithologies include marine deposits, accreted volcanic terrains, passive-margin sediments, and obducted ophiolite, all unconformably overlain by slightly-deformed autochthonous coarse clastic sediment and limestone (Iturralde-Vinent et al., 2016).
The Cuban landscape features a mountainous spine (600-1970 m) descending into low relief coastal plains, except along portions of the south coast where mountains meet the sea.This drainage divide parallels Cuba's east-west orientation, creating rivers that travel relatively short distances from headwaters to base level (Galford et al., 2018).Cuba's climate is tropical wet and dry, with a mean annual temperature of 24.5 °C and average annual precipitation of 1335 mm y -1 .The climate is highly seasonal; ~80% of this precipitation is delivered during the wet season from May-October (Llacer, 2012).
Centuries of agriculture have heavily altered the Cuban landscape (Whitbeck, 1922).Prior knowledge of mass loss at the basin scale is limited to measurements of suspended sediment discharge for short periods between 1964 and 1983 for 32 Cuban rivers (Pérez Zorrilla and Ya Karasik, 1989), and measurements of dissolved loads in five limestone basins with karst (Pulina and Fagundo, 1992).In central Cuba, underlying basin rock type is the primary control on surface water geochemistry (Betancourt et al., 2012), a finding supported by geochemical analyses of river waters from the same basins sampled in this study (Bierman et al., 2020).Dissolved load fluxes carried by Cuban rivers (Bierman et al., 2020), and rock dissolution rates inferred from these fluxes, are consistent with rates reported for other Caribbean islands [Dominica, Guadeloupe, and Martinique from Rad et al. (2013) and Puerto Rico from White and Blum (1995)], and high compared to global data compiled by Larsen et al. (2014a).

Field methods
We collected detrital sediment (n = 26) from the beds of active river channels in central Cuba, representing a variety of basin sizes, average slopes, and lithologies (Fig. 2; Supplement T1-2).Channel morphologies varied, but most streams were incised, and many had exposed bedrock (see Bierman et al. (2020) for photos/descriptions of select field sites).At each site we collected samples for water chemistry analysis and measured channel parameters, including width, depth, and discharge at time of sampling.

Lab methods
We prepared samples for cosmogenic analysis and extracted beryllium and aluminum following the methodology of Corbett et al. (2016).We sieved bulk sediment samples in the lab and used the 250-850 μm grain size fraction for all samples, except for CU-120, which also includes finer material (63-250 μm) due to low quartz content.Sediment samples were chemically etched to purify quartz and remove meteoric 10 Be (Kohl and Nishiizumi, 1992).Twenty-four samples yielded sufficient quartz for analysis.We measured quartz yields for all but one sample (CU-120) by recording the mass of sediment before and after dilute acid etching (Supporting Fig. 1).
We extracted 26 Al and 10 Be at the National Science Foundation/ University of Vermont Community Cosmogenic Facility, using ~5-43 g of quartz per sample (mean = 24 g).We added ~250 μg of Be to each sample using two different in-house made carriers (Supplement T5); the first batch used a low-ratio carrier made from beryl, while subsequent batches used a dilution of low-ratio commercial SPEX carrier.We added Al to samples with insufficient total Al using a commercial SPEX ICP standard in order to reach a total Al mass of ~1500 μg (Supplement T6).Samples were processed in batches of 12, each of which included at least one blank, and two batches included one quality control standard each (Corbett et al., 2019).
10 Be ratios were normalized against standard 07KNSTD3110 with an assumed ratio of 2850 x 10 -15 (Nishiizumi et al., 2007) and 26 Al/ 27 Al measurements were normalized against standard KNSTD with an assumed ratio of 1818 x 10 -15 (Nishiizumi, 2004).Full laboratory replicate sample preparations and measurements of 26 Al and 10 Be agree to within < 2% (Supplement T7; n=2).We corrected Be measurements by carrier type, since samples were prepared using different carriers; we use the average of two process blanks (1.91±1.01x 10 -15 ; 1SD) to correct 10 samples, and the average of four process blanks (4.02±1.00x 10 -15 ; 1SD) for the remaining samples (Supplement T3).We corrected Al measurements using the average of 4 process blanks (4.97±2.94x 10 -15 ; Supplement T4).We subtracted blank ratios from sample ratios and propagated uncertainties in quadrature.

Analytical methods
We extracted drainage basins and then calculated basin slopes and effective elevations (Portenga and Bierman, 2011) using the ASTER Global Digital Elevation Model (Lp Daac), determined underlying basin rock types from the USGS Caribbean layer (French and Schenk, 2004), and utilized precipitation data from the WorldClim dataset (Hijmans et al., 2005) to estimate basin-specific mean annual precipitation (MAP).
We calculated erosion rates using version 3 of the online erosion rate calculator originally described by Balco et al. (2008) and subsequently updated (wrapper: 3.0, erates: 3.0, muons: 3.1, validate: validate_v2_input.m -3.0 consts: 2020-08-26) using the effective elevation (Portenga and Bierman, 2011) calculated for the basin upstream of the sample collection point, a sample thickness of 0 cm, a rock density of 2.6 g cm -3 , and assuming no topographic shielding across this low-relief landscape.We report erosion rates using the Stone/Lal production scaling scheme.
For samples with the highest 10 Be concentrations (n=4), we also measured the concentration of cosmogenic 21 Ne in quartz to further characterize exposure history (Supplement T10).Neon isotope measurements were made at the Berkeley Geochronology Center on aliquots of the same purified quartz samples used for 26 Al/ 10 Be analysis.They were done by vacuum degassing and noble gas mass spectrometry using the method described in Balter-Kennedy et al. ( 2020) and Balco and Shuster (2009).
We used measurements of dissolved loads in stream water (Bierman et al., 2020) and modelled annual flows from GLOH2O (Beck et al., 2015;Beck et al., 2017) to calculate rock dissolution rates for the 25 basins where we were able to collect water samples.To account for the wide range of lithologies in our upstream watersheds, including some with evaporites, we modified the approach used by Erlanger et al. (2021) (Supporting Fig. 2).We removed ions deposited as atmospheric inputs based on published data on dissolved loads in Cuban rainfall (Préndez et al., 2014).We then determined evaporite weathering rates by balancing Na with Cl and Ca with SO4.The remaining Na was used to determine the silicate contribution of Mg and Ca by using an assumed ratio of Na/Mg of 0.25 and Na/Ca of 0.35 (Erlanger et al., 2021).Silicate weathering rates were calculated as the total of SiO2 and HPO4, assumed to result from silicate weathering.Finally, we balanced the remaining Mg and Ca with bicarbonate to determine carbonate weathering rates.
Considering a variety of landscape-scale metrics, we explored the relationship between 10 Bederived erosion rates and calculated rock dissolution (total and silicate, carbonate, and evaporite) rates using linear correlations and their associated p-values.All reported means of sample populations are arithmetic.

Results
Quartz sand, isolated from central Cuban river sediment, has high concentrations of cosmogenic nuclides (0.41 to 12.6 x 10 5 atoms g -1 10 Be and 0.27 to 5.9 x 10 6 atoms g -1 26 Al). 26Al/ 10 Be ratios (Fig. 4, Table 1) vary considerably, ranging from 3.65-8.36(mean = 5.72±1.14,median = 5.83).Sixteen of 24 samples plot below the window defined by continuous exposure and steady erosion on the two-isotope diagram (Fig. 5).Because these 26 Al/ 10 Be data indicate significant burial of quartz during and/or after exposure, many central Cuban drainage basins do not meet the assumption of insignificant nuclide decay inherent in calculations of erosion rates from cosmogenic nuclide concentrations in detrital sediment (Bierman and Steig, 1996).To minimize the impact of violating this assumption, we compare erosion rates based only on the longer-lived nuclide, 10 Be, with landscape scale metrics and dissolved loads.The 10 Be rates, because they cannot properly account for loss of nuclides during burial for samples with depressed 26 Al/ 10 Be ratios, are overestimates of the true rate of erosion.Erosion rates (Supplement T8), calculated from measured concentrations of 10 Be (Supplement T7), differed considerably between sites. 10Be-derived erosion rates (Fig. 3) range from 3.4-189 Mg km - 2 y -1 (mean = 59±52, median = 45).Considered as bedrock lowering rates by assuming a bedrock density of 2.6 g cm -3 , these are 1.3-73 m/My (mean = 23±20, median = 17). 10Be-derived erosion rates in central Cuba are weakly and positively correlated with mean annual precipitation and slope (Fig. 6).
Rock dissolution rates (Fig. 3) range from 10-176 Mg km -2 y -1 (mean = 92±39, median = 97) and are higher than 10 Be-derived erosion rates in 18 of the 23 basins in which we were able to make both measurements.The median rock dissolution rate is 2.15 times higher than the median 10 Be-derived erosion rate.Rock dissolution rates and 10 Be-derived erosion rates are not correlated (Fig. 6).However, when total rock dissolution rates are partitioned into silicate, evaporite, and carbonate rates then the silicate dissolution rate is positively correlated with 10 Be-determined rates of erosion (r 2 = 0.18, p < 0.05, Table 2).Rock dissolution rates are not separable by dominant basin lithology. 10Be-inferred erosion rates are metamorphic > sedimentary, igneous > sedimentary (Fig. 4).
There is lithological dependence of 10 Be-derived erosion rates and the ratio of rock dissolution to 10 Be-derived erosion rates at the basin scale (Fig. 5). 10 Be-derived erosion rates of sedimentary rocks were lower than other rock types (p = 0.01).Samples with the lowest 26 Al/ 10 Be ratios and the highest 10 Be concentrations (Fig. 5) were collected in the north-western part of the field area in basins predominately underlain by sedimentary rocks (Fig. 2).For the most part, samples from basins dominantly underlain by igneous and metamorphic rocks plot to the left on the two-isotope diagram and have higher 26 Al/ 10 Be ratios than quartz from basins underlain by sedimentary rocks.Basins draining primarily sedimentary lithologies have the highest ratio of rock dissolution to 10 Be-derived erosion rates.Seven basins 119,120,121,122,131,and 132) stand out from the rest (Fig. 5) and are clustered in the northwestern part of our field area.These basins have much lower than average 10 Bederived erosion rates (3.4-14.5 Mg km -2 y -1 ), low 26 Al/ 10 Be ratios (3.80-5.87),and rock dissolution rates 1.7-29 times higher than the 10 Be-derived rates of erosion.All but CU-131 are underlain primarily by sedimentary rocks.
Neon isotope measurements (Supplement T10) revealed high total neon concentrations with isotope composition indistinguishable from atmosphere, so excess 21 Ne was indistinguishable from zero and could not be quantified.Expected cosmogenic 21 Ne concentrations in the samples we analyzed, calculated from observed 10 Be concentrations and the assumption of steady erosion (3-6 M atoms g -1 cosmogenic 21 Ne), would comprise less than 2% of the total amount of 21 Ne we observed and would not be detectable at typical measurement uncertainties.The neon isotope measurements are not inconsistent with the 26 Al and 10 Be data, but provide no additional information.

Discussion
In central Cuba, erosion rates inferred from the concentration of 10 Be in river sand vary by more than an order of magnitude.The lowest 10 Be-inferred erosion rate (3.4 Mg km -2 y -1 ; 1.3 m/My) is less than those measured in tectonically stable arid landscapes including Namibia and Australia (Bierman and Caffee, 2001;Codilean et al., 2021).The highest 10 Be-inferred rate (189 Mg km -2 y -1 ; 73 m My -1 ) exceeds those measured in temperate, humid, tectonically stable areas such as the southern Appalachian Mountains (Portenga et al., 2019;Duxbury et al., 2015;Linari et al.) and is similar to or less than rates measured on other Caribbean Islands including Puerto Rico and Dominica (Quock et al., 2021;Brocard et al., 2015;Brown et al., 1995b).
Variability in 10 Be concentration, and thus inferred rates of erosion, between central Cuban drainage basins, many within just a few tens of kilometers of each other with similar basin slope, suggests significant landscape-scale controls on 10 Be concentration and thus mass loss.Indeed, we find that 10 Be-determined erosion rates are positively correlated with slope (R 2 =0.20, p=0.03), mean annual precipitation (R 2 =0.22,p=0.02), and rates of silicate dissolution (R 2 =0.17, p=0.05) (Fig. 6, Table 2).
Erosion rates are lowest for basins underlain dominantly by sedimentary rocks and highest for basins underlain by metamorphic rocks (Fig. 3).However, accurately quantifying rates of denudation (total mass loss) from central Cuban landscapes is complicated by significant export of mass in solution and near-surface quartz enrichment.In the sections that follow, we discuss the 10 Be data in the context of dissolved load export in river water and the landscape-scale insight on active processes provided by dual-isotope measurements ( 26 Al and 10 Be) made in riverine quartz.

Cosmogenic erosion rates underestimate landscape scale mass loss in Cuba
Our data show clearly that significant, landscape-scale mass loss is occurring by solution in central Cuba.Rock dissolution rates exceed, some by more than an order of magnitude, corresponding 10 Be-derived mass loss rates in central Cuba, demonstrating that the cosmogenic nuclide measurements are an incomplete assessment of total mass loss from the landscape.Rock dissolution rates are greater than cosmogenic erosion rates for 18 of the 23 basins we analyzed and the median rock dissolution rate in Cuba is 2.15 times higher than median cosmogenic nuclide-derived rate (Table 1, Fig. 7).Rock dissolution rates that significantly exceed corresponding 10 Be-inferred rates have also been reported from Uganda (Hinderer et al., 2013) and Cameroon (Regard et al., 2016), where they were attributed to the influence of easily weathered volcanic tephra and deep weathering associated with thick regolith, respectively.Most other studies that compare rock dissolution rates and 10 Be-derived erosion rates in the tropics documented rock dissolution rates within the range of cosmogenic nuclide-derived rates (Von Blackenburg et al., 2004;Salgado et al., 2006;Cherem et al., 2012;Sosa Gonzalez et al., 2016b;Quock et al., 2021).Cuba is different.
The discordance between high rock dissolution rates and low 10 Be-derived erosion rates in central Cuba suggests that significant rock weathering is occurring below the depth of most cosmogenic nuclide production (Bierman and Steig, 1996; Fig. 1).The discordance, along with high rates of carbonate and evaporite dissolution in some basins, suggests that many lithologies in our field area are highly susceptible to dissolution.Bierman et al. (2020) attribute high rock dissolution rates and the relationship between stream water chemistry and bedrock type in central Cuba to extensive rockgroundwater interaction along subsurface flow paths, controlled by ongoing bedrock uplift and associated rock fracturing.The prevalence of rock dissolution at depth in Cuba is consistent with findings from other humid, tropical landscapes, including Puerto Rico (White et al., 1998;Kurtz et al., 2011;Chapela Lara et al., 2017;Moore et al., 2019), Guadeloupe, Martinique, Dominica (Rad et al., 2007), and Hawaii (Schopka and Derry, 2012).
We observed no correlation between 10 Be-derived erosion and rock dissolution rates in central Cuba (Fig. 7), in contrast to other studies in the tropics that have observed generally positive correlations (Salgado et al., 2006;Cherem et al., 2012;Sosa Gonzalez et al., 2016b) The lack of correlation suggests that mass loss below several meters, the depth at which most cosmogenic nuclides are produced, is an important component of denudation in Cuba.Discordance between high rock dissolution rates and low 10 Be-derived erosion rates observed in Cuba occurs in basins with different underlying lithologies (Fig. 6).Such widespread discordance suggests deep chemical weathering is occurring throughout central Cuba.
Carbonate weathering dominates river water geochemistry in central Cuba.Our analysis of Cuban water composition suggests that the rate of carbonate dissolution varies widely and in most basins we sampled, exceeds by several-fold the rate of silicate dissolution (Fig. 4).Silicate dissolution rates are low (<25 Mg km -2 y -1 ) and similar between all lithologies.Export rates of elements calculated to reflect the presence of evaporite minerals are also generally low (<35 Mg km -2 y -1 ), except in four basins dominated by sedimentary rocks 121,122,132).Water geochemistry data from four of these basins suggest the presence of significant evaporite deposits due to high concentrations of Cl, SO4, Br, and Na (Bierman et al., 2020).Together these data imply that lithologies underlying the basins we sampled are not uniform and that silicate rocks do not account for most of the dissolved mass loss in at least some, and likely many, of the basins we sampled.
Together, underlying lithology and topography are important controlling factors in how and how rapidly the Cuban landscapes we studied are losing mass by both physical and chemical weathering.
Lowland basins, primarily underlain by sedimentary rocks, on average have low rates of 10 Be-inferred mass loss and high rates of dissolution.Six basins underlain by sedimentary lithologies  have the highest 10 Be concentrations and lowest erosion rates indicating longer total exposure durations for the quartz we analyzed.All are low slope (0.5 to 1.6 degrees).These six basins also demonstrate the greatest disparity between high rock dissolution rates and low 10 Bederived erosion rates (5.7-29X).One basin underlain by igneous rocks (CU-131) has similarly low slope (0.6 degrees) and high 10 Be concentration but a much lower ratio of dissolution to erosion rates (1.7) likely reflecting the paucity of readily soluble minerals.As a result, 10 Be-derived erosion rates are weakly and positively correlated with average basin slope (R 2 = 0.20, p=0.03), but rock dissolution rates are not correlated (R 2 = 0.04, p=0.36) with slope.

Low 26
Al/ 10 Be ratio evidence for a deep mixed surface layer and possible quartz enrichment The 26 Al/ 10 Be ratios suggest that most sediment we collected from central Cuban rivers does not have a simple exposure history. 26Al/ 10 Be data in 16 of 24 sampled basins are inconsistent with steady surface erosion (Fig. 5).Many of the basins with the lowest 26 Al/ 10 Be ratios drain predominantly marine sedimentary lithologies, and have low average basin slopes (0.5-0.7°); the remaining basin drains primarily igneous rocks and has an average basin slope of 0.6°.These are the same seven basins discussed in the section above, all but one of which have high ratios of dissolved load to 10 Be-inferred erosion rates.There is a significant (R 2 = 0.34, p = 0.003) relationship between average basin slope and 26 Al/ 10 Be.
Observed 26 Al/ 10 Be ratios in most of the low-ratio samples are consistent with bioturbation and prolonged near-surface exposure (Struck et al., 2018).We suspect that at least some of the inconsistency between measured 26 Al/ 10 Be ratios and those predicted by a simple, steady-state surface erosion model is due to (deep) soil mixing.Typically, the lower boundary of the simple exposure region of a two-isotope diagram (Fig. 5) is constructed based on the assumption that all grains move monotonically towards the surface at the rate that the surface is eroding (Granger, 2006).Vertical mixing, due to bioturbation or other soil processes taking place in the upper layers of soil, violates this assumption.Within a mixed soil layer, grains circulate at a higher velocity than the erosion rate, and therefore experience an average production rate lower than the surface rate and spend time below the surface where the rate of nuclide decay may exceed the rate of nuclide production.During burial, (Fig. 5).
Rapid chemical mass loss due to the presence of readily soluble evaporite and marine or igneous lithologies in some basins likely enriches the remaining sediment in quartz.The combination of mass loss by rapid rock dissolution and the retention of weathering residuum favored by subdued topography in low-relief basins allows less-soluble material (e.g., quartz) to accumulate at and near the surface, creating thick regolith.Extensive vertical mixing of near-surface soil, as is expected for flat, forested landscapes where the rate of bioturbation is likely high in relation to slow erosion rates, leads to longer residence times for these residual mineral grains, and therefore a lower 26 Al/ 10 Be ratio, in a mixed surface layer compared to a surface eroding at the same rate without vertical mixing.This assertion is supported by the consistency between measured 26 Al/ 10 Be ratios and expected nuclide concentrations and ratios calculated assuming the presence of a mixed surface layer (per Lal and Chen (2005), equation 12).Expected 26 Al/ 10 Be ratios calculated assuming a mixed layer depth of 40-160 cm agree well with measured low 26 Al/ 10 Be ratios from basins CU-119, CU-122 and 132 (Fig. 5).This mixed layer depth range is consistent with the soil depths of 90-150 cm reported for the location of these basins (Bennett and Allison, 1928).In deeply weathered tropical soils, bioturbation can extend to depths of several meters (Von Blackenburg et al., 2004) so it is plausible that mixing depths are even greater than the model suggests, providing an explanation for the 26 Al/ 10 Be ratios measured in CU-120, 121, 131, and 106.We were not able to measure regolith depths in the drainage basins we sampled.
The 26 Al/ 10 Be ratio in other samples (e.g., 118,and 110) is too low to be attributed solely to the effects of a deep mixed surface layer, and requires that some fraction of the sample has experienced both surface exposure and a significant period of burial well below the surface where cosmogenic nuclide production is negligible.Factors that could lead to this low ratio include the incorporation of previously deeply buried sediment through channel avulsion (Wittmann et al., 2011) or incision into terraces (Hu et al., 2011).We conclude terrace storage, along with a combination of quartz enrichment due to high chemical weathering rates of soluble marine rocks in combination with very low slope basins and a deep mixing layer, combine to generate detrital quartz with high concentrations of 10 Be and lower than expected 26 Al/ 10 Be ratios.

Constraining total rates of landscape denudation
The disagreement between high rock dissolution rates and low 10 Be-derived erosion rates raises questions about how to best characterize total landscape denudation rates.It is clear from our data set that neither cosmogenic nuclide measurements nor stream solute flux are capturing all or even, in some cases, the majority of landscape denudation in central Cuba.Evidence for deep rock dissolution presented in section 6.1 suggests that sediments and solutes are being sourced at least partially from different depths in the landscape.Because most mass loss in much of central Cuba occurs in solution (rock dissolution rates are higher than 10 Be-derived erosion rates in 18 of 23 basins), rock dissolution rates typically represent a minimum bound on total landscape denudation.
Treating the removal of mass in solution and through physical erosion as entirely discrete processes happening at different depths in the landscape sets an upper limit on total landscape denudation: the sum of inferred rock dissolution rates and 10 Be-derived erosion rates.Summing 10 Bederived erosion rates and chemical denudation rates increases estimates of total landscape denudation across study basins by a factor of 1.4-30 (mean = 6.3, median = 2.7) above 10 Be-derived erosion rates.
Disregarding the six basins with evidence of evaporite deposits  leads to an average increase of a factor of 2.6 (median = 2.5) above 10 Be-derived erosion rates.
These mean and median values for the basins without evaporites are between the reported CEF of 1.79 for the Luquillo Critical Zone Observatory in humid, tropical Puerto Rico (Riebe and Granger, 2013) and the CEF of 3.2 for the thickly saprolite-mantled, tropical environment of south Cameroon (Regard et al., 2016).These comparisons suggest that for landscapes with a significant proportion of total denudation occurring through deep rock dissolution, summing rock dissolution rates and cosmogenic nuclide-derived rates provides a reasonable estimate of total landscape denudation.
In landscapes like central Cuba, total denudation rates may be difficult to predict based on landscape metrics.Summed chemical denudation rates and cosmogenic nuclide-derived erosion rates are not correlated with rock type, as rock type appears to have opposing influences on these rates (i.e., basins underlain by sedimentary rocks had the highest rock dissolution rates but lowest cosmogenic nuclide-derived rates, Fig. 3).Summed rock dissolution rates and cosmogenic nuclide-derived rates do increase with mean basin slope (R 2 =0.19, p=0.04) and mean basin elevation (R 2 =0.17, p=0.05) (Fig. 5), but those relationships are confounded because 10 Be-derived rates are highest in high elevation, steep basins and rock dissolution rates are highest in low slope, low elevation basins-relationships that are primarily controlled by the influence of rock type on these two different mass loss processes.
In central Cuba, the lack of correlation between rock dissolution rates and 10 Be-derived erosion rates (Fig. 7) suggests a possible mechanism for limiting total reduction in landscape relief.While global data demonstrates significant, positive correlations between 10 Be-derived erosion rates and basin slope and relief (Portenga and Bierman, 2011), accounting for the influence of rock dissolution may alter this dynamic.The possibility of combined physical and chemical processes limiting reductions in relief has significant implications for the study of deeply weathered, high-relief tropical landscapes.The dual importance of rock dissolution in low-lying areas and physical erosion in steeper terrain could explain the relationship behind sustained high relief topography and low 10 Be-derived erosion rates common across some tropical landscapes, such as Brazil (Vasconcelos et al., 2019) or Sri Lanka (Von Blackenburg et al., 2004).Landscapes with high rock dissolution rates and low physical erosion rates appear to be relatively common (Larsen et al., 2014b).As lowlands weather primarily through rock dissolution and high relief areas by physical erosion, total relief would change more slowly than 10 Beestimated rate differentials would suggest.
Regardless of rock type, however, cosmogenic nuclide-derived erosion rates are positively correlated with MAP (R 2 = 0.22, p = 0.02).While MAP does not vary widely across our study basins in central Cuba (956 to 1555 mm y -1 ), this correlation suggests a climatic control on denudation rates across this landscape.This finding is contrary to other studies in the humid tropics (Von Blackenburg et al., 2004), and beyond (Riebe et al., 2001b;Portenga and Bierman, 2011), that have found no correlation between climate variables and cosmogenic nuclide-derived long-term erosion rates but similar to recent findings in humid, temperate Tasmania (Vanlandingham et al., 2022).Since in Cuba 10 Be-derived erosion rates are positively correlated with MAP but chemical denudation rates are uncorrelated with MAP, this trend likely highlights the importance of rainfall in allowing for the physical export of sediment from a drainage basin that is transport-limited rather than weatheringlimited.
Our data clearly demonstrate that cosmogenic nuclide measurements can underestimate total denudation in landscapes with significant rock dissolution at depth, particularly in the tropics.This suggests that similar underestimates of total denudation rates produced by relying on measurements of cosmogenic nuclides may be a factor in other tropical landscapes.While rock dissolution rates in the tropics have been documented as among the highest globally (White and Blum, 1995;Rad et al., 2013;Larsen et al., 2014a), a global compilation of 10 Be-derived erosion rates demonstrated that such isotopically determined rates of erosion are lower in the tropics than in all other climate zones, apart from arid regions (Portenga and Bierman, 2011).The contrast between these two depictions of tropical denudation suggests that 10 Be-derived erosion rates for tropical areas may be incomplete representations of total mass loss from these landscapes because dissolved loads are incompletely accounted for by measurements of 10 Be in river sand.This discrepancy highlights the need for more studies that compare rock dissolution rates and cosmogenic nuclide-derived rates to provide more accurate estimations of total landscape denudation (Vanlandingham et al., 2022).

Conclusions
The first cosmogenic nuclide measurements from the island of Cuba provide insight into how mass is lost from landscapes in humid, tropical settings.Solution plays a large role in total mass flux, and significant mineral dissolution occurs along weathering fronts meters below the landscape surface.
Rock type exerts the primary control on the pace of denudation, and precipitation influences rates of landscape change.We find evidence for thick mixed surface layers in lowland basins and river water chemistry data suggest that deep rock dissolution dominates denudational processes in low slope basins where weathering products remain near the surface for long periods of time.
These findings highlight the necessity of accounting for mass loss by solution at depths below the penetration of most cosmic rays when interpreting cosmogenic nuclide-derived rates in landscapes with the potential for significant rock dissolution.The discrepancy between high rock dissolution rates and low 10 Be-derived erosion rates observed in central Cuba emphasizes how relying on cosmogenic nuclide measurements alone to determine total rates of mass loss from landscapes can lead to considerable underestimates of denudation.Summing mass loss rates in solution with mass loss rates inferred from cosmogenic nuclides provides an upper limit for total mass loss from landscapes when significant rock dissolution occurs below the penetration depth of cosmic ray neutrons.These findings suggest that estimating rock dissolution rates is important when applying cosmogenic nuclides in landscapes, especially those which are humid, tropical, have soluble rocks, and/or have deep weathered regolith.5: 26 Al/ 10 Be two-isotope plots.To permit comparison of data from different locations and elevations on the same plot, nuclide concentrations have been normalized by dividing measured concentrations by calculated mean production rates in the respective drainage basins, using production rate calculations from version 3 of the online exposure age calculator described by Balco et al. (2008) and subsequently updated.In both plots, uncertainty ellipses denote 68% confidence regions for the normalized nuclide concentrations, and the black lines are the boundaries of the simple exposure region (Lal, 1991) calculated using the conventional assumption of steady block erosion without vertical mixing.Panel A shows that 26 Al/ 10 Be ratios from basins with high nuclide concentrations, implying low erosion rates, are systematically lower than predicted by steady erosion without vertical mixing.Panel B shows that this inconsistency can, at least in part, be explained by the presence of a mixed layer of at least 160 cm.Circles show expected steady-state nuclide concentrations in a fully mixed surface layer calculated according to Lal and Chen (2005) for a range of erosion rates and mixed-layer thicknesses, which highlights that sediment derived from a deep mixed layer has lower nuclide concentrations and lower 26 Al/ 10 Be ratios than would be expected if the mixed layer were absent.Sample ID as in panel A.

Figure 1 :
Figure 1: Conceptual diagram showing cosmogenic nuclide production and concentration in a column 550 of soil, saprolite, and rock.The dashed line shows decreasing production of cosmogenic nuclides with depth; solid line shows nuclide concentration with depth, and the white arrow represents mass loss by solution below the depth of significant cosmogenic nuclide production.

Figure 2 :
Figure 2: Maps showing underlying basin geology (panel A; French and Schenk, 2004), elevations 555 (panel B; Lp Daac), and location of study area within the island of Cuba (panel C).Legend for panel A includes the category of the rock units in terms of sedimentary (S), igneous (I), and metamorphic (M).Note that the two marine units in A are separated because they have different chemical load signatures.

Figure 3 :
Figure 3: Maps showing rates of landscape change and isotopic data for each study watershed.A,B.560 Rock dissolution rates.C,D. 10 Be-derived erosion rates.In both maps, darker colors in the basins indicate faster rates of landscape change.E,F.Ratio between rock dissolution rates and 10 Be-derived erosion rates; darker colors indicate larger ratios.G,H. 26 Al/ 10 Be ratio.Darker colors are higher ratios.Box plots show the maximum and minimum values in the lines extending from the box; the upper side of the box represents the upper quartile, the line inside of the box represents the median value, and the 565 bottom of the box represents the lower quartile.Y axis units are the same as shown in the corresponding map legend.

Figure 4 :
Figure 4: Dissolution rates partitioned by lithology.A,B.Silicate.C,D.Carbonate.E,F.Evaporite.570 Darker colors represent higher rates.Box plots show the maximum and minimum values in the lines extending from the box; the upper side of the box represents the upper quartile, the line inside of the box represents the median value, and the bottom of the box represents the lower quartile.
Figure5: 26 Al/ 10 Be two-isotope plots.To permit comparison of data from different locations and elevations on the same plot, nuclide concentrations have been normalized by dividing measured concentrations by calculated mean production rates in the respective drainage basins, using production rate calculations from version 3 of the online exposure age calculator described byBalco et al. (2008) and subsequently updated.In both plots, uncertainty ellipses denote 68% confidence regions for the normalized nuclide concentrations, and the black lines are the boundaries of the simple exposure region(Lal, 1991) calculated using the conventional assumption of steady block erosion without vertical mixing.Panel A shows that 26 Al/ 10 Be ratios from basins with high nuclide concentrations, implying low erosion rates, are systematically lower than predicted by steady erosion without vertical mixing.Panel B shows that this inconsistency can, at least in part, be explained by the presence of a mixed layer of at least 160 cm.Circles show expected steady-state nuclide concentrations in a fully mixed surface layer calculated according toLal and Chen (2005) for a range of erosion rates and mixed-layer thicknesses, which highlights that sediment derived from a deep mixed layer has lower nuclide concentrations and lower 26 Al/ 10 Be ratios than would be expected if the mixed layer were absent.Sample ID as in panel A.

Figure 6 :
Figure 6: Relationship between measured 10 Be-derived erosion rates, chemical denudation rates, and the sum of chemical denudation rates and 10 Be-derived erosion rates to basin characteristics.Different shaped/colored points represent the dominant underlying rock type in that basin.Plots with p > 0.05 are shown with a gray background.Small numbers in upper right are R 2 values; * indicates p ≤ 0.05, ** p < 595 0.01.

Figure 7 :
Figure 7: Scatterplot of rock dissolution rates vs. 10 Be-derived erosion rates.A. Data shown for 600 calculations of carbonate, silicate, and evaporite dissolution rates.B. Data shown for dominant lithologies underlying each sampled basin.Horizontal lines extending from the points demonstrate the uncertainty associated with the calculation of 10 Be-derived erosion rates.Histograms on axes show distribution of data.Dashed line is 1:1.