Comparable early-stage decomposition but contrasting underlying drivers between surface and cave habitats along an elevational gradient

Decomposition is a major contributor to ecosystem respiration, determining the carbon emission and nutrient cycling rates. Our current understanding of decomposition dynamics and their underlying drivers has mainly focused on surface habitats but largely ignored in subterranean environments. Here we studied abiotic and microbial drivers of early-stage litter decomposition inside and outside caves along an elevational gradient in Tenerife. We found comparable decomposition rates ( k ) and litter stabilizing factors ( S ), with contrasting drivers and elevational variation. At the surface, we observed a mid-elevational trend in k


Introduction
Nutrient accumulation and cycling are paramount for biodiversity build-up and ecosystem stability across biomes (Baldock, 2007).Particularly processes involving carbon cycling influence the productivity in an ecosystem through energy redistribution and nutrient availability (Smith andSmith, 2012, Baldock, 2007).One key component in global carbon cycling is organic matter decomposition, where decomposer organisms and climatic variables, such as temperature and water availability, transform complex organic molecules into simpler molecules (Findlay, 2013), ushering the continuation of the carbon cycle.Decomposition of plant litter contributes greatly to ecosystem respiration and is responsible for a large part of carbon emission (Djukic et al., 2018).Decomposition processes lead to break down of complex molecules (short-term carbon release), but also to the stabilization of labile molecules (transformation of a fraction of labile molecules into recalcitrant molecules).The latter is mainly driven by environmental variables and contributes to long-term carbon storage (Elumeeva et al., 2018, Prescott, 2010).
Elevational gradients provide powerful study systems with potentially steep environmental changes over relatively short distances that can be used to investigate how climatic variables and biotic drivers influence these ecological processes.Studies that have examined decomposition along elevational gradients have highlighted the importance of temperature (Wang et al., 2009, Salinas et al., 2011, Vitousek et al., 1994), humidity (Murphy et al., 1998), precipitation (Vitousek et al., 1994, Upadhyay et al., 1989), and plant community composition (Sundqvist et al., 2011, Murphy et al., 1998).However, the magnitude of these drivers and their interactions vary between gradients (Becker and Kuzyakov, 2018, Murphy et al., 1998, Salinas et al., 2011, Vitousek et al., 1994), implying that changes in climatic variables have varying consequences on decomposition and carbon storage across the globe.Moreover, most studies have focused on surface ecosystems and our understanding of decomposition in subterranean habitats such as caves along elevational gradients are lacking (Ravn et al., 2020).Cave ecosystems across the globe (Culver and Pipan, 2019) provide a window to explore the contributions to nutrient cycling in the vast dimensions of underground habitats (Mammola et al., 2019), otherwise inaccessible to humans.These ecosystems differ substantially from surface habitats as they exhibit stable diurnal and seasonal environmental conditions (Lauritzen, 2018, Castaño-Sánchez et al., 2020), limited nutrient availability (oligotrophic) (Simon et al., 2007), and invertebrate and microbial communities specialized to life in the dark (Simon et al., 2007, Ortiz et al., 2014, Gonzalez-Pimentel et al., 2018, Gonzalez-Pimentel et al., 2021, Hathaway et al., 2014b, Riquelme et al., 2015, Hathaway et al., 2014a, Simon et al., 2003).Thus, the parameters and potential drivers of decomposition in caves conceivably differ from surface environments (Ravn et al., 2020), but current insights into biological processes in caves is limited, hampering our understanding of the relative contributions and drivers of decomposition and their vulnerability to climate change.
To help fill these knowledge gaps, we examine early-stage decomposition parameters and their abiotic and microbial drivers in cave ecosystems and at the surface along an elevational gradient (from 66 up to 2300 m a.s.l.) on the Island of Tenerife, The Canary Islands (Fig. 1, Table S1).Tenerife has multiple caves from sea level to more than 3000 m a.s.l.(Bacallado et al., 1995), providing an excellent opportunity to study the dynamics of organic matter decomposition both above and below ground.We investigated the impact of a series of soil abiotic factors and bacterial community compositions on early (three months) decomposition parameters (decomposition ratek; i.e., short-term carbon cycling and litter stabilizing factor -S; i.e., long-term carbon storage) characterized using the standard Tea Bag Index -TBI (Keuskamp et al., 2013).First, we hypothesized that decomposition would be slower in caves compared to surface habitats due to the nutrient-limited stable environment within caves (Simon et al., 2007) and differences in decomposer communities (Reboleira et al., 2022, Ravn et al., 2020).Second, as decomposition tends to be positively associated with water content (Murphy et al., 1998) and temperature (Becker andKuzyakov, 2018, Vitousek et al., 1994) in surface habitats across elevational gradients, we predicted that the decomposition parameters should be positively linked to soil properties and nutrient levels in both surface and cave habitats.Finally, we hypothesized that decomposition parameters will be associated with different groups of bacterial taxa in caves and in surface soils due to marked differences in above and below ground microbiomes (Reboleira et al., 2022).

Fig. 1.
Measuring decomposition parameters across caves and their surface habitats along an elevational gradient on Tenerife.A. The location of the island of Tenerife among the seven main islands of The Canary Islands, Spain.B. Tenerife is a volcanic island, reaching 4000 m a.s.l. and is consequently characterized by extraordinary variation in climatic conditions, as illustrated in the schematic that also indicates the mid-elevation cloud layer that experiences most precipitation and most diverse plant communities.The illustration further gives the locations and elevations of the eight caves we sampled to obtain a gradient of environmental conditions across both elevations and intercardinal NE and SW directions.C. To estimate the decomposition rate (k) and stabilizing factor (S), we buried three Green and three Rooibos teabags at 8 cm depth in the surface soil, following the standard Tea Bag Index (TBI) protocol.Given the thin layer of sediment within caves, the teabags were only buried under the available sediment layer (approximately 2 cm deep).Tea bags were left for three months, after which they were recovered for analyses.(For interpretation of the references to colour in this figure legend, the reader is referred to the web version of this article.)

Sampling sites
Fieldwork was conducted in eight caves on the island of Tenerife (Fig. 1, Table S1), Canaries, Spain, between November 2019 and February 2020.The caves were located along an elevational gradient from 66 m a.s.l up to 2,300 m a.s.l (Fig. 1, Table S1).

Assessment of decomposition parameters
We examined decomposition of organic matter in caves using the standard Tea Bag Index (TBI), which has been validated from studies across biomes in surface ecosystems (Keuskamp et al., 2013, Djukic et al., 2018).Accordingly, we used tea bags of the slowly decomposable Rooibos tea and the faster decomposable Green tea, as pre-made "litterbags" for the determination of decomposition parameters (Fig. 1).The use of two tea types with contrasting decomposability allows estimating a decomposition trajectory using a single measurement in time (Keuskamp et al., 2013).Furthermore, the use of the commercially available pre-made "litterbags" (tea bags) reduces variation caused by differences in preparation methods and material used.Based on weight measurement of tea bags before and after the three-month incubation period, we calculated parameters describing decomposition rate per day (k) and litter stabilization factor (S) using the well-established and publicly available protocol from the global initiative TeaTime4Science (http://www.teatime4science.org).S indicates the fraction of labile material remaining after the incubation period, while k represents the rate of breakdown of the labile fraction.S can be greatly influenced by environmental variables and also indicates the fraction of recalcitrant compounds that are resistant to decomposition (Duddigan et al., 2020, Keuskamp et al., 2013).Thus, higher S implies higher level of carbon sequestration, with limited breakdown of labile material (Keuskamp et al., 2013, Prescott, 2010, Becker and Kuzyakov, 2018).
We buried three pairs of Rooibos and Green tea bags in the sediment at each surface and cave sampling locality (Fig. 1).At the surface, we ensured that we buried teabags in natural habitats, with the minimum impact of anthropogenic activities (except for in San Marcos: Table S1).Tea bags were installed at the same time (November 2019) at all locations and collected after three months (February 2020).Subsequently, the tea bags were dried for 48 h at 70 ℃ before weighing.The mass loss of tea bags during the three-month period was used to calculate k and S (Keuskamp et al., 2013).

Environmental variables and soil properties
During the experiment, we recorded temperatures every two hours within caves and at the surface using dataloggers TidbiT v2 Temp UTBI-001.Temperature data were downloaded through an Optic USB Base Station (BASE-U-4) and HOBOware Software.These temperatures were used to calculate the average temperature during the three months of the experiment.Surface dataloggers from San Marcos and Perdiz caves were lost during the experiment.
During the installation of the tea bags, soil samples were collected in triplicates (near each pair of teabags) from each surface and cave habitat, at each of the eight elevations, and frozen until further analyses.In the laboratory, soil was sieved (2 mm), carefully mixed and any roots were removed.Each fresh soil sample was divided into five subsamples.One subsample of 10 g was suspended in demineralized water (ratio 1:5) and used for determination of pH (pHM240 MeterLab) and conductivity (SevenCompact Conductivity).
To obtain soil and microbially-bound nutrient and carbon concentrations, we used water as the extractant.A range of different extractants (e.g., H 2 O, K 2 SO 4 , CaCl 2 , KCl, NaHCO 3 ), and of different molar strengths, are used for extraction of soluble and microbial fractions of C, N, P and other elements in the soils (Rennert et al. 2007;Clemmensen et al. 2008;Ravn et al. 2017;Mclaren and Buckeridge, 2019;Schwalb et al. 2023).While it should be noted that salt extractants may often lead to higher element concentrations, we chose water as the extractant in our study for two reasons.This was firstly because water reflects actual uptake conditions relevant for both plants and microbes in many soils, and secondly because it enabled us to lyophilize the extract for subsequent combustion in the elemental analyzer, as was done by e.g., Ravn et al. (2017) and Reboleira et al. (2022).While the choice of extractant may affect element amounts, this was unproblematic for our experiment as the extractant was the same for all sites and soil types (cave vs. surface).The total microbial biomass C and P was estimated using the chloroform fumigation method (Brookes et al., 1982, Vance et al., 1987).To achieve this, a subsample of 20 g soil was suspended right away for one hour in demineralized water (ratio 1 g soil:5 mL H 2 O) and filtered (Whatman GF/D).Another 20 g subsample was incubated for 24 h in a vacuum desiccator with chloroform before extraction and filtration.All filtrated extractions were kept frozen until analysis.Samples were thawed, centrifuged for 10 min at 4300 rpm prior to further analyses.Before freeze-drying we added 100 and 50 μL 2 M HCl to 100 mL extracts of non-fumigated soil and 50 mL extracts of fumigated soil, respectively.All material from each freeze-dried sample was packed individually in tin capsules and analyzed on an isotope ratio mass spectrometer (IRMS; Isoprime) connected to an Eurovector CN elemental analyzer to determine total dissolved C (TDC) and total dissolved N (TDN).Phosphate (PO 4 3− ) in fumigated and non-fumigated extracts were measured using flow injection analysis (FIAstar 5000 Analyzer).Extractions of non-fumigated soil were analyzed for nitrate (NO 3 -) and ammonium (NH 4

+
) content using flow injection analysis (FIAstar 5000 Analyzer).Estimation of microbial biomass C and P was based on the difference between fumigated and non-fumigated samples using an extractability factor of 0.45 for C and 0.40 for P (Brookes et al., 1982, Clemmensen et al., 2008, Ravn et al., 2017).To determine soil water content, another soil subsample was dried at 60 • C for three days.Subsequently the dry soil was ground in a ball mixer and approximately 10 mg of soil was packed in tin capsules and analyzed on IRMS coupled to an elemental analyzer to determine C and N concentration.Soil organic matter (SOM) was determined through loss on ignition by burning of another sample of soil at 550 • C overnight.

Characterization of soil microbiomes
DNA was extracted from 0.25 g of soil from each location on the fifth subsample of soil using the Qiagen Power soil kit (Qiagen, Germany) following the manufacturer's guidelines.Initial PCRs to identify samples with bacterial DNA were conducted using two primers targeting the v4 region of the 16S rRNA gene ('SB711 and 'SA504) and following a well establish protocol for the primers (Bodawatta et al., 2020).All the samples amplified positively and were sent to the Microbiome core at the University of Michigan for amplicon MiSeq sequencing on an Illumina platform.
Amplicon sequences were analyzed in QIIME2 (Bolyen et al., 2019) using the DADA2 pipeline (Callahan et al., 2016) and assigned to ASVs at 100% similarity.Consequently, ASVs were assigned to taxonomy using the SILVA 132 bacterial reference library (Quast et al., 2013).Subsequently mitochondrial, chloroplast and Archaeal sequences were removed from the data set.We generated a bacterial phylogeny using a align-to-tree-mafft-fasttree function in QIIME2.For further analyses we removed communities with <1500 sequences (three samples).Prior to analyzing the communities, we rarefied the dataset using the sample with the lowest number of sequences (12,181 sequences per sample) using reafy_even_depth function in phyloseq package (McMurdie and Holmes, 2013).

Statistical analysis
All subsequent statistical analyses were conducted in R 4.1.0(R Core K.H. Bodawatta et al.Team, 2020).Prior to analyzing soil properties, we investigated the Pearson's correlations between different soil parameters using the ggpairs function in GGally package (Schloerke et al., 2021) in the two habitat types.Multiple parameters revealed significant positive correlations with one another (Fig. S1).For example, at the surface, total N, C, C:N ratio, PO 4 3-and microbially-bound elements were significantly correlated with SOM, while NO 3 -, inorganic nitrogen, and microbiallybound elements were significantly correlated with soil pH.We also observed similar correlations in caves, where SOM was significantly correlated with total N, C, NO 3 -, NH 4 + , and inorganic nitrogen, while soil water content was significantly correlated with C:N ratio, PO 4 3-and microbial bound elements.Thus, we selected the following variables that did not correlate with one another for the subsequent analyses: soil pH, water content, SOM, and DOC (Fig. S1).We investigated the influence of habitats (surface and cave), different sites at different elevations (cave identity) and the interaction between these two variables on soil properties using generalized linear models (glms).The distributions of dependent variables were visualized with histograms and for each model, the link function was adjusted accordingly to identity or log.We explored elevational trends of soil properties, decomposition parameters, and temperature variation using linear and polynomial regressions in the dplyr package (Wickham et al., 2019).The best fit model was identified using AIC values and comparing two models statistically.We conducted these analyses separately for surfaces and caves.To investigate how k and S associated with soil properties, we conducted similar linear and polynomial regression analyses.Alpha diversity indexes of bacterial communities (chao1 richness estimate and Shannon's diversity index) were calculated using the diversity function in the microbiome package (Lahti and Shetty, 2017).We also calculated Faith's phylogenetic diversity index for bacterial communities using the pd function in the picante package (Kembel et al., 2010).We examined the elevational relationships between bacterial alpha diversity indexes using linear and polynomial regressions.Similar analyses were conducted to investigate the associations between elevation and chao1 richness estimates of the 10 dominant bacterial phyla, that collectively accounted for 94.4% bacterial sequences.
Bacterial community level differences were visualized using nonmetric multidimensional scaling (NMDS) analysis in vegan package.Differences in surface and cave microbial communities (measured with Bray-Curtis and weighted UniFrac distances) were examined conducting Fig. 2. Contrasting elevational trends in abiotic variables and decomposition parameters.The results of our evaluation of the effect of elevation (x-axis) for both abiotic factors (A-D) and decomposition parameters (E-F).Grey shaded area indicates standard errors of data points, where each data point represents a single replicate for a given site (n = 3 per site).Soil pH (A) and organic matter content (B) exhibited significant mid-elevational peaks in surface soils but with negative trends (significant only for pH) within caves.Soil water content (C) and three-month average temperature (D) were not associated with elevation on the surface, contrasting caves that exhibited a significant mid-elevational drop in water content (C) and a significant negative association with temperate (D).Elevation associated with decomposition rates in surface soil, where there was a significant mid-elevational peak, contrasting the pattern in caves (E).Stabilizing factor was not associated with elevation in either habitat (F).For the full statistical analyses, see Table S4.
K.H. Bodawatta et al. permutational multivariate analysis of variance (PERMANOVA) with 10,000 permutations using the adonis2 function in the vegan package (Oksanen et al., 2019).Subsequently we investigated the influence of elevation, pH, water content, SOM, and DOC on community level differences using the envfit function in the vegan package.To test if specific bacterial genera were associated with k and S, we performed Pearson's correlations between relative abundances of bacterial genera with the decomposition parameters using the microeco package (Liu et al., 2021).Significant values were adjusted using false discovery rates (FDR).Data was visualized using ggplot2 (Wickham, 2016) and viridis packages (Garnier, 2018).

Decomposition parameters and soil properties differ along the elevational gradient
Overall, we observed similar average k (caves: 0.0125 day − 1 ± 0.007; surface: 0.0105 day − 1 ± 0.003) and S (caves: 0.5156 ± 0.2431; surface: 0.5462 ± 0.0654) between surface and caves (Table S2), indicating that decomposition parameters in general are similar in above and below ground habitats.However, the decomposition parameters differed notably between elevations (Fig. S1, Tables S2  and S3), where differences in levels of k and S were higher between surfaces and caves in lower than higher elevations.There was a midelevational peak in k at the surface, but not in caves (Fig. 2E).In contrast, S did not show an elevational trend in either habitat (Fig. 2F).
Soil properties, such as pH, water content, and SOM differed significantly between caves and surfaces (Fig. S2, Tables S2 and S3), and their levels and magnitude differences between cave and surface habitats varied along the elevational gradient (Fig. S2, Table S3).At the surface, we observed significant mid-elevational peaks in pH and SOM (Fig. 2A, B, Table S4).Within caves, pH decreased significantly with the elevation, while SOM did not (Fig. 2A, B).There was a non-significant negative association with elevation in soil water content at the surface, while cave water content dropped significantly at mid-elevation (Fig. 2C, Table S4).Average temperature of the three months decomposition period was not associated with elevation at the surface but decreased significantly with increasing elevation in caves (Fig. 2D, Table S4).

Abiotic drivers of k and S differ between surface and caves
Soil pH was neither associated with k or S in both habitats.However, at the surface, k was significantly positively associated with water content and SOM (Fig. 3A and 3C, Table S5), and negatively associated with the average three-month temperature (Fig. 3B, Table S5).Moreover, dissolved organic carbon (DOC) was not linked to k at the surface.Fig. 3. Abiotic factors affect decomposition rate only in surface soil and stabilizing factor mainly in caves.At the surface, decomposition rates were significantly positively associated with soil water content (A) and organic matter (C) and negatively associated with temperature (B).None of the abiotic parameters were associated with decomposition rate in caves.Within caves, stabilizing factor was significantly negatively associated with soil water (E) and temperature (F), as well as the content of organic matter (G) and dissolved organic carbon (H).This contrasted surfaces, where dissolved organic carbon level was the only parameter significantly associated with stabilizing factor (H). Colors and labels as in Fig. 2.
In caves, we did not find an association between k and abiotic properties (Fig. 3, Table S5).Only DOC demonstrated a significant polynomial association with S at the surface, where lower values were associated with moderate levels of DOC (Fig. 3H).We found the opposite trend in the caves, where higher S was associated with moderate levels of DOC (Fig. 3H) and significantly negatively associated with soil water content (Fig. 3E, Table S5).SOM levels were positively associated with S in caves but this trend was mainly driven by the cave Honda de Güímar (Fig. 3G, Table S5).Moreover, we observed a negative association between S and three-month average temperature in caves, except in Honda de Güímar (Fig. 3F, Table S5).Overall, this points to levels of organic matter and temperature as important drivers of decomposition rate at the surface, Fig. 4. Surface and cave bacterial community compositions are shaped primarily by soil pH, elevation, and soil organic matter content.A. Average relative abundance of major bacterial phyla identified in cave and surface soil zones (n = 3) depicted similar relative abundances of phyla between surface and caves, with Proteobacteria (surface: 23.1% ± 7.6%; caves: 27.8% ± 6.9%), Acidobacteria (surface: 16.9% ± 5.2%; caves: 19.2% ± 7.4%), Actinobacteria (surface: 14.9% ± 4.7%; caves: 10.7% ± 7.3%), Chloroflexi (surface: 14.9% ± 16.5%; caves: 8.3% ± 4.4%), and Planctomycetes (surface: 8.5% ± 2.9%; caves: 7.5% ± 2.8%) being the most dominant.Bars only represent the relative abundance of the 16 most abundant bacterial phyla.B. However, the microbial community composition differs significantly between caves and surfaces (Permutational multivariate analysis of variance: PERMANOVA (Bray-Curtis) : F = 3.308, R 2 = 0.0538, p < 0.0001; PERMA-NOVA (UniFrac) : F = 3.101, R 2 = 0.0481, p = 0.0021).The bacterial community composition also differed between caves along the elevational gradient implying that community composition of these microbial communities is influenced by elevational changes (PERMANOVA (Bray-Curtis) : F = 2.453, R 2 = 0.2392, p < 0.0001; PERMANOVA (UniFrac) : F = 2.321, R 2 = 0.2521, p < 0.0001).Furthermore, the variation in community composition was larger in caves than surface communities, as evident from the larger 95% confidence intervals (dashed ellipses) of the NMDS ordination plot (stress = 0.1465) based on Bray-Curtis distances.Envfit analyses (Table S7) indicated that this pattern was affected significantly by differences in soil pH, elevation, and organic matter (dashed grey arrows) (each data point represents a single sample from a given site; n = 3 per site).C. NMDS plots for each of the two habitats (caves: top panel; stress = 0.1212; surfaces: bottom panel; stress = 0.1307) indicated that variation in community composition within caves was significantly affected by soil pH, water content, soil organic matter, and elevation, while surface microbiomes were only influenced by soil pH and elevation (dashed grey arrows) (for the full envfit analyses results, see Table S8).
while temperature and water availability influence stabilization in caves.

Soil bacterial communities differ markedly between surfaces and caves
Overall, we acquired 1,306,074 bacterial sequences (average ± SD: 29,023 ± 7,805) from 22 cave samples and 23 surface soil samples.These sequences were assigned to 26,136 amplicon sequence variants (ASVs) (Table S6) from 35 phyla (Fig. 4A).At the bacterial community level, we observed significantly different compositions between surfaces and caves (Fig. 4B), aligning with previous findings (Reboleira et al., 2022).Community composition also differed between caves along the elevational gradient (Fig. 4B), implying that environmental changes associated with elevation strongly influence bacterial communities.This inference was confirmed by envfit analyses of the influence of abiotic parameters on bacterial community composition, which revealed significant effects of elevation, pH, and SOM levels on microbiome compositions (Fig. 4B, Table S7).
Individual envfit analyses at surface and cave communities revealed that elevation and pH were significantly linked to bacterial community composition in both caves and surfaces (Fig. 4C, Table S8).However, in caves, soil water content and organic matter availability (SOMwith Bray Curtis distances and DOCwith UniFrac distances) also significantly impacted bacterial community composition (Fig. 4C, Fig. S3, Table S8).These factors were not associated in surface bacterial communities.This suggests that organic matter is more important for structuring cave than surface microbiomes, where elevation-associated soil pH is more important.

Elevational trends and nutrient associations of bacterial diversity differed between surface and caves
Bacterial richness (chao1) estimates and Faith's phylogenetic diversity (PD) of surface communities showed mid-elevational peaks (Fig. 5, Table S9), and these alpha diversities were significantly positively associated with pH (Fig. 5, Table S9).In contrast, cave microbial diversity was not associated with elevation, but all cave microbial diversity indexes were significantly negatively associated with soil DOC levels (Fig. 5, Table S9), implying that caves with reduced nutrient availability harbor more diverse bacterial communities.This aligns with previous findings from bacterial communities of sub-arctic caves, where nutrient poor caves harbored more diverse and complex bacterial communities (Reboleira et al., 2022).We only found significant polynomial associations between chao1 richness estimates and decomposition parameters in surface communities, where both k and S were lower when bacterial richness was moderate (Table S10, Fig. S4).This indicates that overall community level bacterial diversity is not strongly associated with decomposition parameters.
To determine whether individual bacterial phyla direct decomposition, we investigated the associations between decomposition parameters and richness of the 10 dominant bacterial phyla (accounting for 94.4% of the sequences) individually.At the surface, we observed significant mid-elevational peaks in the richness of the most abundant phyla, such as Proteobacteria, Acidobacteria, Actinobacteria, and Bacteroidetes (Fig. S5, Table S11), which contribute strongly to the Fig. 5. Abiotic factors associate differently to bacterial diversities between surfaces and caves.Elevation and soil pH were strongly associated with surface alpha diversity indexes, while cave alpha diversity was associated with levels of dissolved organic carbon.This implies that elevation and soil acidity impact the diversity in surface communities, while increased soil nutrients negatively impact bacterial diversities in caves.Labels as in Fig. 2. observed overall mid-elevational peak (Fig. 5).In contrast, within caves, only the richness of Bacteroidetes and Verrucomicrobia showed significant mid-elevational trends (Fig. S5, Table S11).This suggests that the non-significant associations between overall alpha diversity indexes and elevation in cave communities is driven by high variation in phylumspecific trends.
Associations between richness of the dominant phyla and decomposition parameters revealed trends that were specific to surface and cave environments.At the surface, the richness of Proteobacteria, Planctomycetes, and Verrucomicrobia were significantly positively associated with k, indicating that higher richness of these phyla was linked to higher decomposition rates (Fig. 6A, Table S12).Only Planctomycetes associated positively with S in surface communities (Table S12).Caves lacked any associations between phylum-level richness and k.In contrast, richness of Actinobacteria and Firmicutes significantly increased with S, while Rokubacteria richness decreased with S (Fig. 6B and Table S12) in caves.This indicates that only a subset of bacterial taxa is associated with different decomposition parameters in above and below ground habitats.These results also align with the observations of association between abiotic parameters and decomposition rates, where abiotic drivers mainly associated with k in surface, while these drivers were associated with S within caves (Fig. 3). .k and S were significantly associated with the richness of different habitatspecific bacterial phyla in caves and surfaces, driven by their correlations with specific bacterial genera.The richness of Proteobacteria, Planctomycetes and Verrucomicrobia was positively associated with decomposition rates in surface soils (A).In contrast, stabilizing factor significantly positively associated with the richness of Actinobacteria and Firmicutes, but negatively with Rokubacteria (B).Aligning with this, Pearson's correlation analyses (false discovery rate (FDR) adjusted p values; 0.01 < * < 0.05; 0.001 < ** < 0.01) revealed significant positive correlations between the relative abundance of 11 surface soil bacterial genera with decomposition rate (C).This sharply contrasted in caves, where none of the bacterial genera were significantly associated with decomposition rates, but the relative abundance of 25 bacterial genera correlated positively and three negatively with stabilizing factor (D).
Overall, this implies that specific bacterial taxa are associated with decomposition rate in the surface, while in caves bacterial taxa are associated with stabilizing factor.

Associations between specific bacterial genera and decomposition parameters differ between surfaces and caves
Phylum level associations with decomposition parameters alone do not reveal which bacterial taxa/genera are linked to decomposition in these habitats.Thus, to explore this, we determined correlations of k and S with relative abundances of bacterial genera separately in caves and surface.These analyses revealed that the relative abundance of 11 bacterial genera at the surface, of which six belonged to the Proteobacteria, were significantly positively associated with k (Fig. 6C).This attests to the significant association we observed between richness of the phylum Proteobacteria and k in the surface (Fig. 6A), suggesting that members of this phylum play a marked role in decomposition at the surface.Within caves, there were significant correlations between 28 bacterial genera and S (Fig. 6D), of which 25 were positive and three were negative.Eight of the positively associated genera belonged to the Firmicutes and five to Actinobacteria, for both of which richness was also significantly positively associated with S (Fig. 6B), indicating a positive relationship between these phyla and recalcitrant compounds in cave environments.In caves, none of the bacterial genera were associated with k.Overall, bacterial genera-level analyses further support that bacterial taxa within caves are associated with longer-term storage of C, while taxa at the surface are associated with short-term carbon cycling.

Discussion
Our findings document comparable decomposition (k: short-term carbon cycling and S: long-term carbon storage) (Keuskamp et al., 2013) in surface and cave habitats, underlining the importance of cave environments to carbon cycling in these regions.The magnitude of differences in k and S between habitats changed along the elevational gradient, with some caves having even higher decomposition rates than their respective surfaces.This contradicts our first hypothesis that presumed oligotrophic caves (Simon et al., 2007) would exhibit reduced decomposition compared to surfaces.Consistent with this, but in contrast to general belief, we did not find nutrient-poor conditions within caves but rather levels of organic matter comparable to -and at times even higher than -surfaces.Nutrient levels varied across caves, likely leading to the observed similarities and disparities in decomposition between surfaces and caves.The influx of organic matter into caves depends on multiple variables, including the volume of water seepage from the surface (percolating organic matter from surface into the cave), width and dimension of the cave entrance, and the movement of surface taxa that can introduce organic matter, such as bats and pigeons (Ravn et al., 2020).Taken together, this implies that caves, depending on their nutrient influx and environmental conditions, may significantly contribute to carbon cycling and influence regional carbon budgets.
Despite similar average levels of decomposition parameters, shortterm carbon cycling was strongly linked with both abiotic and biotic factors at the surface, while they were associated with long-term carbon storage in caves.The mid-elevation peak in decomposition rates at surfaces aligns with the trend observed in another tropical elevational gradient (Becker and Kuzyakov, 2018), but contradicts surface studies where decomposition rates either decreased (Salinas et al., 2011, Vitousek et al., 1994) or increased (Murphy et al., 1998) with elevation.However, some of these studies investigated small elevational spans (Murphy et al., 1998) that may not have fully captured elevational trends in decomposition (c.f.Nogues-Bravo et al., 2008).The effect of elevation on surface decomposition rates in our study appears to be governed by the combination of environmental variables, where higher decomposition is associated with lower temperatures, higher water availability, and higher levels of SOM.This contrasts the positive effect of temperature on k that was observed along a Peruvian elevational gradient (Salinas et al., 2011), but agrees with the strong influence of high water availability at higher elevations in Arizona, USA (Murphy et al., 1998).Mid-elevation sites on Tenerife experience higher levels of moisture and harbor diverse laurel forests (del-Arco andDelgado, 2018, Bryant et al., 2008), which in combination may account for higher levels of water availability and SOM.Moreover, SOM levels were also positively correlated with multiple other soil nutrients, including nitrogen.Higher nitrogen levels in less acidic soils have been found to positively affect decomposition rates in multiple biomes (Liu et al., 2010, Cusack, 2013, Averill and Waring, 2018).Thus, the presence of more nitrogen in the more alkaline mid-elevational soils of Tenerife likely increased rates of short-term carbon cycling in surfaces.Overall, temperature, precipitation, and nutrient availability has different impact on decomposition across different elevational gradients (Vitousek et al., 1994, Salinas et al., 2011, Murphy et al., 1998, Becker and Kuzyakov, 2018), implying that changes to these conditions due to climate change will have unpredictable consequences for decomposition in surface habitats across regions.
Similar to abiotic variables, bacterial community characteristics were mainly associated with decomposition rates at the surface.We observed positive correlations between richness within specific bacterial phyla (e.g., Proteobacteria, Planctomycetes, and Verrucomicrobia) and decomposition rates.Members of these phyla encode a wide range of carbohydrate-active enzymes (Dedysh and Ivanova, 2019, Bao et al., 2019, Tao et al., 2020), attesting to their ability to degrade plant-derived compounds.The likely importance of particularly Proteobacteria on decomposition was supported by correlations to relative abundances of multiple bacterial genera belonging to the phylum.Overall bacterial community richness is thus not associated with short-term carbon cycling at the surface, but only a subset of taxa is linked to decomposition, as seen in other habitats (Daebeler et al., 2022).However, to understand how these bacterial phyla and genera influence decomposition rates, investigation of microbial functions is needed.Furthermore, it is important to investigate associations between decomposition and other soil microbes, such as fungi that play critical roles in degradation of complex carbon molecules (van der Wal et al., 2013).
The mid-elevational trend in surface bacterial community richness and phylogenetic diversity suggests that soil properties associated with different elevations tend to have a stronger effect on bacterial communities than elevation alone.This mid-elevational trend aligns with patterns observed along a Japanese elevational gradient (Singh et al., 2012) while other studies have found a mid-elevational drop (Singh et al., 2014, Shen et al., 2020), a negative (Singh et al., 2014, Shen et al., 2019), or a positive (D'Alo et al., 2022) association in bacterial richness with elevation.The strong association between soil pH and alpha and beta diversity of bacterial communities indicates that pH is the major driver of bacterial community composition, as has been found in other elevational (Shen et al., 2019, Shen et al., 2020, Singh et al., 2014) and soil depth (Kim et al., 2014) gradients.This importance of pH on governing soil bacterial communities, further iterate the vulnerability of these complex communities to changes in soil acidity due to anthropogenic activities.
Our insights into the long-term carbon storage in the environment, revealed that S was not strongly associated with abiotic parameters at the surface.This contradicts patterns observed in a few other studies, including in tropical (Becker and Kuzyakov, 2018) and alpine (Elumeeva et al., 2018) elevational gradients, where S was found to positively associate with cold and dry conditions.In our elevational gradient, colder elevations experienced higher water availability (Fig. 2), potentially indicating sub-optimal conditions for litter stabilisation.Moreover, S is also influenced by other factors such as forest stratification and soil clay content (de Godoy Fernandes et al., 2021), which we did not test.These disparities between studies again underline how regional properties in elevation-associated abiotic drivers impact decomposition parameters.Use of the stabilizing factor to investigate long-term carbon storage is still at its early stages (Kwon et al., 2021, Djukic et al., 2018).Thus, to better understand the drivers of long-term carbon storage along elevational gradients, and the interaction between different decomposition parameters, we need more studies across climatic and environmental zones (e.g., latitudinal and elevational gradients).And we need to interpret the data on stabilizing factor based on teabag index with caution as the assumptions related to calculation of S may not always be met in nature (Mori et al., 2022).
In contrast to surfaces, S in caves was linked with multiple abiotic variables, such as reduced water availability and cooler temperatures.This aligns with previous findings from surface habitats across elevational gradients (Becker and Kuzyakov, 2018;Elumeeva et al., 2018).Conditions within caves are more stable across seasons and throughout the day than surfaces (Lauritzen, 2018), implying that caves with reduced water availability and temperatures are characterized by reduced decomposition year around, turning cave ecosystems with the right environmental conditions into carbon sinks.The significant associations between SOM and S in caves was mainly driven by the Honda de Güímar cave, which may experience a particularly high influx of organic matter due to the wider cave entrance and the presence of a large pigeon roost.The consistently high S and low k within this cave supports that decomposition is inhibited, which can be caused by abiotic factors such as the soil clay content (de Godoy Fernandes et al., 2021) and iron content (Jin et al. 2022).
We observed positive associations of reduced decomposition (high S) between the richness of certain bacterial phyla and relative abundance of bacterial genera within caves.However, it is difficult to deduce the potential mechanisms underlying these trends.For example, untested abiotic parameters that inhibit decomposition (Xiaogai et al. 2013) may also influence certain bacterial taxa, without a direct interaction between the two variables.This can lead to correlations between certain microbial taxa and decomposition parameters, without direct associations.However, thorough examination of microbial functions and interactions between bacterial taxa are needed for proper assessment of the association between microbes and S within cave environments.

Conclusions
Our comparison of decomposition parameters between surface and cave ecosystems provides the first glimpse of the differences in organic matter cycling in these habitats along an elevational gradient.Despite differences in associations of decomposition rates and stabilizing factors with tested variables between surface and cave habitats, decomposition parameters in both habitats tend to be sensitive to similar environmental variables, such as temperature and water availability.This implies that changes in temperature and precipitation due to climate change will strongly influence different aspects of decomposition in above and below ground habitats.These changes may have cascading effects across trophic levels, ultimately disrupting natural homeostasis of ecosystems.The influence of stable conditions within caves on long-term carbon storage underlines their potential important contributions to regional carbon budgets.However, we remain on the verge of decoding the vulnerability of biological process within cave ecosystems to aboveground anthropogenic activities, which can threaten the stable conditions within these underground habitats.Therefore, it is timely to delve further into cave ecosystems to understand the potential effects of global change and anthropogenic pressures on disturbing natural nutrient cycles in these less explored yet globally distributed habitats.
Measurements of soil parameters can be found in the Table S2.Microbiome sequences are deposited at the Sequence Read Archive (SRA) database in GenBank (Accession PRJNA800241) and microbiome data can be found in the Table S6.

Declaration of Competing Interest
The authors declare that they have no known competing financial interests or personal relationships that could have appeared to influence the work reported in this paper.
Fig.6.k and S were significantly associated with the richness of different habitatspecific bacterial phyla in caves and surfaces, driven by their correlations with specific bacterial genera.The richness of Proteobacteria, Planctomycetes and Verrucomicrobia was positively associated with decomposition rates in surface soils (A).In contrast, stabilizing factor significantly positively associated with the richness of Actinobacteria and Firmicutes, but negatively with Rokubacteria (B).Aligning with this, Pearson's correlation analyses (false discovery rate (FDR) adjusted p values; 0.01 < * < 0.05; 0.001 < ** < 0.01) revealed significant positive correlations between the relative abundance of 11 surface soil bacterial genera with decomposition rate (C).This sharply contrasted in caves, where none of the bacterial genera were significantly associated with decomposition rates, but the relative abundance of 25 bacterial genera correlated positively and three negatively with stabilizing factor (D). Overall, this implies that specific bacterial taxa are associated with decomposition rate in the surface, while in caves bacterial taxa are associated with stabilizing factor.