Spatiotemporal dynamics of dissolved organic matter in subtropical karst cave waters identified by optical properties

Dissolved organic matter (DOM) is composed of degradation by-products and freshly produced compounds, resulting in a heterogeneous mixture of carbohydrates, hydrocarbons, polyhydroxy phenols, fatty acids, amino acids, and quinines (Carter et al., 2012; Wei et al., 2018). Because of complicated structures and various types of attached functional groups (Li and Hur, 2017), DOM has important roles in terrestrial and aquatic environments associated with carbon and nitrogen cycles, distribution and fate of pollutants, and residence and loss of nutrients (Liu et al., 2021; Xu et al., 2018). Thus, tracing changes in DOM concentration and structure could provide a better understanding of the geochemical processes associated with the carbon cycle in ecosystems.

Karst areas play an important role in the global carbon cycle because weathering of carbonate rocks affects carbon dynamics and budgets, and underground carbon dioxide (CO₂) storage is an important part of the terrestrial carbon flux (Cao et al., 2020; Liu and Zhao, 2000; Serrano-Ortiz et al., 2010). Soil carbon, including organic and inorganic carbon, is a major driving force for carbonate dissolution and affects the hydrochemistry properties of surface water and subsurface water in karst areas (Dong et al., 2023). Soil organic carbon can be transported by percolating water from overlying soil into caves (Baker et al., 1997; van Beynen et al., 2001; Ban et al., 2008), and the fingerprint retained in different cave water sites might have spatiotemporal discrepancies due to the high heterogeneity of soil organic matter and complexity of the karst aquifer. Therefore, distinguishing these spatiotemporal discrepancies in dissolved carbon levels and compositions could provide information on organic matter transit and fate in karst cave systems.

Several characterization methods have been applied to explore DOM structure; among them, optical measurement is considered to be reliable because it provides information on chemical bonds, functional groups, and molecular spatial configurations of DOM (Ni and Li, 2020; Kellerman et al., 2014). Ultraviolet-visible (UV-vis) absorption spectroscopy is one of the most common analytical methods used to determine the spectral characteristics of DOM. Natural organic matter includes organic compounds with heavily conjugated systems that have strong UV signals in the absorption spectrum in aqueous solution (Yakimenko et al., 2016; Dodla et al., 2008; Wei et al., 2018; Weishaar et al., 2003; Yu et al., 2016). The absorption coefficient is the most common parameter for identifying DOM composition (Li and Hur, 2017). The absorbance coefficients at wavelengths 254 (a₂₅₄), 285 nm (a₂₈₅), and 350 nm (a₃₅₀) indicate aromatic, fulvic, and lignin fractions of DOM, respectively (Ni and Li, 2020; Kida et al., 2019; Piirsoo et al., 2012; Hernes and Benner, 2003). Absorption ratios, i.e., the ratios of the absorption coefficients at two wavelengths, are indices frequently used to trace the source and determine DOM quality. The E2/E3 ratio (i.e., the absorption ratio at 250 nm to 365 nm) reflects the DOM molecular size, whereas the E2/E4 ratio (i.e., the absorption ratio at 254 nm to 436 nm) indicates the origin of DOM (Li and Hur, 2017; Battin, 1998; Jaffé et al., 2004; Hur et al., 2006). The E4/E6 ratio (i.e., the absorption ratio at 465 nm to 665 nm) indicates the degree of DOM condensation, which is correlated to molecular size, carboxyl content, and O/C and C/N atomic ratios (Chen et al., 1977; Helms et al., 2008; Senesi et al., 1989). In addition to absorption ratio indices, the spectral slope coefficient is an important parameter for characterizing DOM structure features. Bricaud et al. (1981) primarily used spectral slopes to distinguish DOM sources. Helms et al. (2008) calculated the slopes of logarithmically fitted lines between 275–295 nm and 350–400 nm by converting the UV absorbance spectrum to natural logarithms. Moreover, the ratio of two spectral slopes (namely, S_R) was found to be negatively correlated with DOM molecular weight. Therefore, spectral slopes and absorption ratios can be used to investigate the geochemistry of DOM.

UV-vis absorption spectral parameters are widely used to analyze DOM characteristics (Artinger et al., 2000; Carter et al., 2012; Frey et al., 2016; Harvey et al., 2015; Kellerman et al., 2015; Liu et al., 2021; Zhou et al., 2019; Zhu et al., 2017); however, the application of optical properties to trace the spatiotemporal dynamics of DOM levels, compositions, and potential linkage with inorganic carbon in karst cave systems has been rarely explored. The inhomogeneity and high permeability of karst aquifers influence the heterogeneous hydrological processes during downward infiltration, leading to large uncertainties and multiple solutions in the responses of cave waters to external environmental parameters (Zhang et al., 2023). Unraveling the dynamics of DOM in karst cave waters may provide insight into its origin and fate, allowing for a better understanding of carbon cycling in this type of terrain. The objectives of this study were to (1) trace the spatiotemporal changes of DOM levels and compositions in a karst cave system; (2) discuss the relationship between DOM dynamics and hydrological path patterns; and (3) explore the linkage between dissolved organic carbon (DOC), HCO₃⁻, and partial pressure of aqueous CO₂ (pCO₂) and further reveal the origin and fate of DOM in a karst cave system. This study hoped to obtain some information regarding DOM geochemical dynamics in highly heterogeneous hydrological processes and its potential impact on carbon cycling in karst cave systems.

Mahuang Cave is a nontourist cave located in Southwest China. The cave entrance is situated at an altitude of 720 m, and the bedrock layer (including the overlying soil) covering the cave is ~100 m thick. It is one of the primary branches of Shuanghe Cave (107°02′30″N–107°25′00″N; 28°08′00″E–28°20′00″E), which ranks as the fifth longest cave in the world (with a total documented length of 300.127 km as of 2021; Zhang et al., 2023). With the influence of the Yanshan and Himalayan tectonic movements, the Shuanghe cave system was uplifted and formed the northeast, northwest, and southwest fold-fracture zones due to regional tectonic stress (Zhang et al., 2023). The cave system developed in the Upper–Middle Cambrian Loushanguan Formation (∈2–3ls) and Lower Ordovician Tongzi Formation (O1t) with components of dolomitic gypsum, gypsum, dolomite-bearing anhydrite, dolomite-bearing secondary limestone, and dolomitic syngenetic breccia (Dong et al., 2023; Wang et al., 2020). Calcium carbonate deposits, such as stalactites, stone curtains, helictite, and stalagmites, are developed inside the cave.

The region is characterized by a typical subtropical humid monsoon climate, which is mainly controlled by the East Asian monsoon. During the monitoring period, the minimum and maximum air temperatures were 4.2 °C and 26.0 °C, occurring in December 2020 and August 2020, respectively. In addition, 80% of precipitation (~810.9 mm) occurred during the wet season (May to October). The overlying vegetation is mainly composed of evergreen and broad-leaf forests and shrubs.

The field work was conducted monthly from August 2020 to December 2021. In Mahuang Cave, five sites (MH0, MH2, MH3, MH4, MH5) were selected along cave passageways from the entrance (Fig. 1). Cave air temperatures were monitored in situ using a Kestrel 4500 portable small weather station (Kestrel, Boothwyn, Pennsylvania). Discharge rates were manually timed with measuring cylinders. Water temperature, pH, and electrical conductivity (EC) were measured in situ with a WTW Multi 3630 portable multiparameter water-quality probe (WTW, Weilheim, Germany). HCO₃⁻ levels were determined titrimetrically in situ using an alkalinity test kit (with a measuring range of 0.1–10 mmol/L; Merck, Darmstadt, Germany), with an accuracy of 0.1 mmol/L. Water samples were filtered through 0.45 μm polyethersulfone (PES) microfiber filters into 50 mL high-density polyethylene (HDPE) containers for measuring the DOC concentration and for spectrum analysis. After filtering through 0.22 μm PES microfiber filters, water samples were stored in 50 mL HDPE containers until further use. Samples for cation measurement were acidified to pH <2 by adding 1:1 nitric acid to maintain ionic activity. Headspace or air bubbles were avoided during the sealing process. Samples were stored in the dark in a refrigerator at 4 °C and transported back to the laboratory within 2 d with icepacks for protection.

Anions and cations were measured using an ion chromatograph (VISTA MPX; Varian, Palo Alto, California) and an inductively coupled plasma–optical emission spectrometer (ICS-90; Dionex, Sunnyvale, California, USA), respectively. The DOC concentrations were measured with a C/N Analyzer (multi N/C 3100; Jena, Germany) after acidification to pH <2 with 1.5 mol/L HCl to eliminate the inorganic carbon. UV-vis absorption spectra were obtained using a 1 cm quartz cuvette at room temperature on a UV-vis spectrophotometer (Cary 300, Agilent, Santa Clara, California). A 200–800 nm wavelength range was measured at a scanning rate of 300 nm·min⁻¹ and a 1 nm interval. Ultrapure water served as a blank, which was subtracted from the optical measurements. Results are expressed as the mean of duplicate or triplicate analyses, with a coefficient of variation ≤2% and standard deviation.

The absorption coefficients a₂₅₄, a₂₈₅, and a₃₅₀ were measured at wavelengths 254 nm, 285 nm, and 350 nm, respectively. Absorption coefficients a₂₅₀, a₂₆₅, a₃₆₅, a₄₃₆, a₄₆₅, and a₆₆₅ were then used to calculate the E2/E3, E2/E4, and E4/E6 ratios. Two spectrum slopes, 275–295 nm and 350–400 nm, were applied to calculate S_R as follows:

where a_λ is the absorption coefficient at λ nm, and L is the path length of the quartz cuvette.

The pCO_2(w) of cave waters was calculated using the PHREEQC program (version 3, U.S. Geological Survey):

where K_h and K₁ are the temperature-dependent Henry’s law and first dissociation constant for CO₂ gas in water, respectively.

Spearman’s correlation coefficient was used to explore the correlation between parameters. The Mann-Whitney U test was used to verify temporal (wet and dry season) variations, while the median test was used to determine spatial (MH0, MH2, MH3, MH4, and MH5) differences. Principal component analysis (PCA) was conducted to determine variables with significant correlations and identify possible sources of variability in the karst cave system. All statistical analyses were performed using SPSS (v23.0, IBM Corp). All the figures were prepared using OriginLab Origin Pro 2022 and GraphPad Prism 9.

Throughout the observational period from August 2020 to July 2021, a notable variability in discharge rates was observed (Fig. 2). Specifically, MH0 and MH2 displayed comparatively lower average discharge rates of 0.32 and 0.10 mL/s, respectively. Conversely, MH3 and MH4 exhibited elevated average discharge rates of 2.70 and 10.05 mL/s, respectively.

MH5 is an underground stream at an elevation of 673 m. The four drip points (MH0, MH2, MH3, and MH4) were categorized based on their hydrological characteristics, using data on peak discharges and coefficients of variation (percentage standard deviation variation around the mean) collected during the monitoring period. (Fig. 3; Zhang et al., 2023; Smart and Friedrich, 1987; Tooth and Fairchild, 2003). MH0 was classified as a seasonal drip, while MH2 exhibited characteristics of both a seasonal drip and seepage flow. Both MH0 and MH2 were classified as slow-flow sites. MH3 was classified as a percolation stream, while MH4 was classified as shaft flow, both of which were high-flow sites. Discrepancies in hydrological type could lead to differences in response to the external environment.

Summaries of the physicochemical characteristics of the cave waters are presented in Table 1 and Table S1.¹ Temporally, discernible variations were noted in cave air temperature, cave water temperature, pH, EC, and pCO₂ between the wet and dry seasons (Fig. 3). The cave’s internal air temperature ranged from 3.6 to 20.2 °C, displaying a narrower span in the wet season (16.1 ± 0.3 °C) compared to the dry season (11.5 ± 0.6 °C), as depicted in Figure 2. Moreover, the range of cave air temperature was inversely correlated with distance from the entrance, with MH0, situated closest to the cave entrance, exhibiting the highest variation, indicative of heightened sensitivity to external environmental factors. Similar patterns were discerned in cave water temperature, where values were 15.4 ± 0.3 °C and 11.8 ± 0.5 °C during the wet and dry seasons, respectively. The pH values demonstrated a shift from 8.32 ± 0.04 in the wet season to 8.55 ± 0.02 in the dry season. This transition coincided with rising air temperatures and precipitation, suggesting a nuanced relationship. EC displayed an opposing trend, with values of 332 ± 7 μm/cm in the wet season and 311 ± 8 μm/cm in the dry season. The pCO₂ exhibited a range of 281–1949 ppmv, characterized by higher levels in the wet season (848 ± 85 ppmv) and lower levels in the dry season (397 ± 29 ppmv). HCO₃⁻ concentrations (188.74 ± 5.30 in the wet season and 176.60 ± 6.04 in the dry season) and discharge rates (5.36 ± 2.38 mg/L in the wet season and 1.47 ± 0.46 mg/L in the dry season) followed a similar trend, albeit without statistical significance (p > 0.05), as illustrated in Figure 4.

Spatially, Ca²⁺, Mg²⁺, EC, and discharge rates exhibited significant seasonality (p < 0.05) during the monitoring period (Fig. 5). It is worth noting that variations in EC, Ca²⁺, and Mg²⁺ for MH0 and MH2 were relatively smaller in the dry season compared with the other three sites. This could be related to the different hydrological types, as mentioned in “Hydrological Type Classification” section. MH0 and MH2 were slow-flow sites with relatively longer retaining times and stronger buffering effects, and thus less sensitive to external seasonal changes.

Aquatic DOC concentrations ranged from 0.34 to 4.50 mg/L (Fig. 6). Two concentration peaks were observed during the monitoring period, including November 2020 to January 2021 and May 2021. Nevertheless, temporal analysis failed to reveal statistical significance in DOC levels between the wet and dry seasons (p > 0.05). Absorbance coefficient a₂₈₅ exhibited significant variability in the wet (10.46 ± 0.39 m⁻1) and dry seasons (9.41 ± 0.38 m⁻¹) (p < 0.05). Fulvic acid serves as a proxy for terrestrial signals derived from soil, reflecting the contribution of water-soil erosion. Consequently, the observed elevation in a₂₈₅ values during the wet season was attributed to increased seepage volume. The E2/E4 and E4/E6 ratios exhibited significant differences (p < 0.05) in the wet (4.38 ± 0.12 m⁻¹ and 1.29 ± 0.04 m⁻¹, respectively) and dry (3.78 ± 0.16 m⁻¹ and 1.06 ± 0.03 m⁻¹, respectively) seasons, indicating that the origin and DOM condensation could be associated with hydro-meteorological conditions in the karst cave system. The spectral slope ratio, S_R, was inversely correlated to DOC molecular weight, ranging between 0.10 and 10.73 with no temporal differences (p > 0.05). Spatially, both a₂₅₄ (11.71 ± 0.37 m⁻¹, 11.18 ± 0.73 m⁻¹, 9.65 ± 0.54 m⁻¹, 10.29 ± 0.73 m⁻¹, and 10.49 ± 0.86 m⁻¹ for MH0, MH2, MH3, MH4, and MH5, respectively) and S_R (0.21 ± 0.04, 1.11 ± 0.47, 0.99 ± 0.53, 1.88 ± 1.14, and 1.74 ± 0.71 for MH0, MH2, MH3, MH4, and MH5, respectively) showed significant variability (p < 0.05) during the monitoring period, suggesting that the type of hydrological pathways in the karst aquifer was a potential factor that impacted DOM fate in the cave system.

Aqueous DOC concentrations in Mahuang Cave were relatively lower in August 2020. This was attributed to regional abnormal precipitation (monthly average: 28.8 mm), causing less DOC flushing from overlying soil to the cavity and mitigation, as driven by the piston effect. According to a previous study (Dong et al., 2023), DOC stored in soil pores and fissures could be flushed out with rising rainfall in the months following the late wet season. Such a trend was observed during the late dry and early wet periods. DOC concentrations experienced a reduction during the early wet season, primarily attributed to the dilution effect resulting from elevated discharge rates induced by high-intensity rainfall events. The E2/E3 and E4/E6 ratios were lower in the initial dry season, indicating that DOM in cave water had a larger molecular size and higher condensation degree in this period. The E2/E4 ratio was used as a proxy for estimating the relative composition of autochthonous versus terrestrial DOM (Hur et al., 2006), which also exhibited lower values, suggesting more terrestrial DOM input (such as lignin) and weakened microbial activities during this period.

Multiple linear regressions were used to identify the relationships between DOC concentrations and DOM absorption coefficients. As shown in Equations 4 and 5, opposite correlations were observed between the wet and dry seasons:

DOC had a strong positive association with a₂₅₄ in the wet season and with a₃₅₀ in the dry season. This indicates that different compositions dominated DOC concentrations under different hydro-meteorological conditions in the cave system. A possible reason could be that with increasing rainfall, DOM with relatively higher aromaticity, which was previously stored in soil pores and bedrock fissures during the dry period, was flushed due to the piston effect, resulting in a strong positive association between DOC and a₂₅₄ values in the early wet season. A positive relationship between the aromatic fraction and DOC was also observed in the Florida Coastal Everglades (Jaffé et al., 2004) from March to May. In September, a₂₅₄ showed a higher average value, which was in accordance with the results observed for the Longchuan River (Ni and Li, 2020), possibly due to the rainfall after the relatively dry month of August. In addition, the microorganisms alternatively consumed labile DOM fractions and caused the accumulation of refractory DOM (Ni and Li, 2022; Hu et al., 2018), such as the high aromatic fraction a₂₅₄. This result is in accordance with the positive correlation between DOC and the E2/E3 ratio (Fig. 7A).

After entering late autumn, a large amount of lignin from decomposed plant debris was flushed into the cave system with rainwater and slowly filtrated downward. The a₃₅₀ serves as an indicator of terrestrial inputs and vascular plants (Ni et al., 2019; Kolic et al., 2014), which exhibited a strong relationship with DOC. Freshly flushed labile DOM would be preferentially utilized by microorganisms, resulting in the accumulation of recalcitrant DOM during the wet season. However, this process exhibited variations in response to meteorological characteristics, and the contribution of recalcitrant DOM became somewhat more prominent during the dry season. The transition process from high-molecular-weight DOM flushed into the karst aquifer to low-molecular-weight DOM occurs gradually, “step by step” (Ni and Li, 2023; Ni et al., 2023). It implies that refractory DOM may undergo incomplete degradation, resulting in the formation of labile DOM. In the late dry season (March and April), the E2/E3 and E4/E6 ratios increased, indicating that DOM with a smaller molecular size and lower condensation degree was generated. The E2/E4 values also increased because the higher temperature caused increased microbial activities (Fig. 6). Overall, the meteorological characteristics served as important constraints of the predominant DOM composition.

Aquatic CO₂ within karst cave water, partially originating from the decay of organic matter, descends into the deeper unsaturated zone. This process may share common cycling pathways with DOM within the aquifer or be influenced by the dynamics of DOM. As highly heterogeneous compounds, variations in DOM levels and composition in a diverse karst aquifer space could be controlled by hydro-meteorological conditions. Under increased temperature and precipitation conditions, HCO₃⁻ and associated hydrochemical parameters, including EC, Ca²⁺, and Mg²⁺, exhibited a positive correlation with a₂₈₅ (Fig. 7A). The smaller molecular size of a₂₈₅ suggests a potential association with water-rock processes and the introduction of terrigenous organic matter. Fulvic acid, serving as a terrestrial signal originating from soil, displayed consistent variations alongside HCO₃⁻, EC, Ca²⁺, and Mg²⁺, potentially influenced by water-soil erosion resulting from rock weathering. These results are consistent with a previous study at Mahuang Cave, which demonstrated that carbonate dissolution was the predominant process in the wet season (Zhang et al., 2023). However, DOM fractions showed an insignificant association with pCO₂. During this period, discharge rates were much higher, and karst fissures and pipes were filled with freshwater, which might have inhibited CO₂ degassing.

As discharge rates decreased with lower rainfall, the pipes and fissures could have been relatively open and well ventilated (Fairchild et al., 2000; Banner et al., 2007), facilitating CO₂ degassing. During this period, heightened a₃₅₀ values primarily governed DOM dynamics and exhibited a notably positive correlation with pCO₂ (Fig. 7B). This suggests a coordinated trajectory between CO₂ and terrestrial DOM inputs. Our findings were similar to the study conducted in the Longchuan River, Southwest China (Ni et al., 2019), which demonstrated that the generated riverine CO₂ had a significantly positive correlation with lignin compounds and aquatic respiration of DOM, explaining the pCO₂ levels in the dry period.

The highly diverse nature of the karst zone determines the heterogeneous hydrological processes, with most underground flow exhibiting mixed behaviors (Zhang et al., 2023; Markowska et al., 2015; Worthington et al., 2016), thus affecting the organic matter fingerprint from the overlying soil. Consequently, exploration of the drip-water hydrology types of Mahuang Cave could provide insight into the regulators controlling the dynamics of DOM. The four drip points were categorized based on their hydrological characteristics (see “Hydrological Type Classification” section) into seasonal drip (MH0), seepage flow/seasonal drip (MH2), percolation stream (MH3), and shaft flow (MH4). The variations of discharge rates revealed that MH0 and MH2 were located in zones with a low-flow subsystem, while MH3 and MH4 were located in zones with a fast-flow subsystem. To further analyze their hydrologic paths in heterogeneous aquifers, the hydrological paths of the four drip-water points in Mahuang Cave were simulated based on the responses of discharge rates to rainfall. As shown in Figure 8, the aboveground subsurface space in the karst region was divided into the atmosphere, soil, epikarst, and vadose zone (including cave cavity) components, with the hydrological paths of the cave water sites indicated. During the wet season, higher recharge flow rates filled the fissures and pipes with freshwater, resulting in faster discharge rates. In contrast, reduced rainfall caused lower recharge flow rates and discharge rates, and the pipes and fissures remained relatively open, with intermittent flow occurring in several fissures possibly.

The a₂₅₄, S_R, EC, Ca²⁺, and Mg²⁺ parameters exhibited spatial differences (p < 0.05) during the monitoring period, suggesting a considerable impact from heterogeneous hydrological processes on DOM variation and water-rock interaction in the karst cave system. MH0 was identified as a seasonal drip, characterized by a higher coefficient of variation, which was more sensitive to seasonal rainfall changes with associated time lags (Tooth and Fairchild, 2003). MH2 exhibited a similar pattern in hydrological behaviors. The difference was that MH0 may have fewer supply channels than MH2, resulting in higher seasonal variations of discharge rates. These two slow-flow drip sites had relatively higher mean a₂₅₄ values than the other three sites, possibly due to different residence times and recharge sources. Extended residence time depleted the labile fraction, resulting in increased accumulation of refractory DOM due to incomplete microbial degradation in fissures. Furthermore, the EC, Ca²⁺, and Mg²⁺ variations were higher in these slow-flow sites, especially EC in the dry season, which could be explained by a longer residence time and, consequently, more dissolution.

MH3 is a percolation stream feature, exhibiting relatively stable water recharge throughout the year. We can infer a stable aquifer recharge for this site because of the low coefficient of variation, and no additional intermittent flow was observed. For MH3, although channels were well developed in the karst aquifer, the hydrochemical variation of underground water was noticeable due to the low storage and recharge capacity (Zhang et al., 2023). The relatively more continuous and higher-speed flow possibly reduced the accumulation of refractory DOM in fissures, leading to smaller fluctuation of a₂₅₄ in MH3. The shaft flow of MH4 had the highest drip-water volume, coefficient of variation, and fluctuation among the monitoring sites and could rapidly respond to external changes. A storage reservoir at MH4 could be inferred according to the previous study (Zhang et al., 2023), suggesting that a particular threshold of water input must occur before flow becomes preferentially routed to this site (Tooth and Fairchild, 2003). In addition, the S_R trend was consistent with these results. MH5, which had the lowest fluctuation among the sites, might have been affected by multiple recharge sources, leading to a relatively stronger buffering capacity.

To further reveal the effects of different hydrological path patterns on the dynamics of DOM, PCA was used to obtain relationships among DOM concentration and composition variables to determine their possible roles in Mahuang Cave; the results are summarized in Tables 2 and 3. In the wet season, MH0 was correlated with a₂₈₅, HCO₃⁻, and pCO₂ on PC1 (57.32%), indicating a geochemical driver linking HCO₃⁻ and pCO₂ with DOM, while MH2 was related to HCO₃⁻ and pCO₂ (PC1, 33.26%). DOC composition at MH0 was mainly controlled by fulvic and lignin fractions (PC2, 41.25%), while DOC at MH2 was dominated by lignin fractions (PC3, 33.07%). Although both MH0 and MH2 were slow-flow sites, different hydrological patterns determined their distinct dissolved carbon dynamics. For example, MH2 flowed much slower, possibly resulting in higher consumption of labile DOM and enhanced water-rock interaction. MH3 was correlated with a₂₅₄, a₂₈₅, HCO₃⁻, and pCO₂ on PC1 (50.84%), while PC2 explained 29.79% in DOC and a₃₅₀ variation, indicating the fractions that were preferentially utilized. Both the aromatic fraction and fulvic fraction exhibited strong association with dissolved inorganic carbon, while the lignin fraction dominated DOC composition. As a percolation stream, MH3 was characterized by a rapid flow rate with stable aquifer recharge from several channels and exhibited mixed behaviors. MH4 was related to a₂₅₄, a₂₈₅, a₃₅₀, and pCO₂ on PC1 (61.59%), indicating common cycling pathways for DOM fractions and CO₂, while PC2 (28.53%) was correlated with DOC and HCO₃⁻, implying the coupling of inorganic and organic carbon. This discrepancy could be explained by faster carbon exchange at MH4 during the wet season. Large recharge could mitigate the piston effect, resulting in a rapid runoff when heavy rainfall events occurred (Fig. 8). MH5 was correlated to a₂₅₄, a₂₈₅, and a₃₅₀ on PC1 (51.34%), while PC2 (41.60%) was associated with DOC, HCO₃⁻, and pCO₂. This suggests that organic carbon was coupled with inorganic carbon; however, an insignificant association between inorganic carbon and specific organic carbon fractions was observed due to the complex sources. In the dry season, MH0 was correlated with DOC and a₂₅₄ on PC1 (49.66%), while PC2 (49.01%) was associated with a₂₈₅, a₃₅₀, HCO₃⁻, and pCO₂. MH2 was related to a₂₅₄, a₂₈₅, a₃₅₀, HCO₃⁻, and pCO₂ on PC1 (46.21%), while PC2 (39.38%) was correlated to DOC. These results indicate the impact of hydrological patterns on subsurface reworking. With less rainfall input and longer residence time in the aquifer, metabolism enhanced the refractory DOC accumulation (Ni and Li, 2022; Catalá et al., 2015). As a result, the lignin and fulvic fractions were consumed at MH0, while a cotrajectory could be inferred at MH2 between DOM fractions and inorganic carbon. MH3 was associated with a₂₅₄, a₃₅₀, and pCO₂ on PC1 (41.49%), while PC2 (32.93%) correlated with a₂₈₅ and HCO₃⁻. This indicates the common pathways from overlying soil to cavity among the aromatic fraction, lignin fraction, and CO₂ in the dry season. PC2 showed a strong correlation with a₂₈₅ and HCO₃⁻, as discussed in the “Potential Linkage Among Dissolved Carbon Species” section. At MH4, DOC consisted of a₂₅₄, a₂₈₅, and a₃₅₀ (PC1 53.57%); however, HCO₃⁻ and pCO₂ exhibited insignificant associations with the organic fraction. Nevertheless, the fingerprint of the overlying soil organic matter could have been retained in the cave water, since MH4 is a rapid material-exchange site. At MH5, inorganic carbon was related to a₂₅₄ and a₂₈₅ (PC1, 44.59%), while a₃₅₀ correlated with the DOC level (PC2, 36.17%).

In a karst area, the internal linkage between DOC and dissolved inorganic carbon is essential but difficult to quantify, especially in these systems. In contrast to other natural water bodies, karst water geochemistry is well known for inorganic-organic carbon transformation mechanisms and allochthonous inputs (Ni and Li, 2022; Sun et al., 2021). Freshly produced low-molecular-weight DOM preferentially degrades to CO₂ (Mayorga et al., 2005). The produced CO₂ can interact with carbonatites, influencing karstification and the carbon sink. This transformation not only establishes a linkage between organic and inorganic carbon but also identifies an important and elusive part of the carbon cycle in karst cave systems. According to our findings, DOC contributed to the hydrochemical features in the cave waters, which has been previously undervalued in karst cave research. In addition, due to the highly diverse nature of the karst zone and varied hydrological processes, the fate of DOM could be impacted by different hydrological paths, as discussed in section 4.3. Divergent subsurface reworking processes result in distinct DOM fingerprints originating from overlying soil at various drip-water sites. These variations can contribute to differences in stalagmite organic laminae formation, emphasizing the need for caution in interpreting these formations due to the uncertain microbial contributions. Further work is needed for high-frequency and long-term monitoring to capture the DOC and dissolved inorganic carbon coupling dynamics more accurately.

Karst areas have a unique carbon cycle with special DOC and dissolved inorganic carbon coupling dynamics. Therefore, the temporal and spatial dynamics of DOC and the linkage with karstification deserve attention. The present study demonstrated spatiotemporal variations in DOC levels and composition, explored the potential relationship with dissolved inorganic carbon, and determined the effect of different hydrological path patterns on DOM dynamics. The results demonstrated that DOC in cave waters was constrained by both terrigenous and biogenic drivers. Meteorological conditions significantly impacted the variation in DOC levels and compositions in karst cave waters. On the one hand, meteorological conditions controlled the DOC input from the overlying soil and flushing by the piston effect. On the other hand, at higher air temperatures and precipitation, DOC had a lower molecular size, lower degree of condensation, and higher microbial activity. However, as a result of the different recharge patterns and residence times determined by different hydrological processes, discrepancies were observed for the consumption of DOC fractions and dissolved inorganic carbon coupling, impacting hydrochemical features. The special linkage between DOC and dissolved inorganic carbon under meteorological and hydrological controls should be considered in future studies of carbon cycling in karst cave systems.

This research was cofunded by the National Natural Science Foundation of China (42161048) and the Science and Technology Program of Guizhou Province (ZK [2021] General 191). We are grateful to all the staff of Shuanghe Cave National Geopark for their assistance with the field work.