Secondary tropical forests recover dung beetle functional diversity and trait composition

Secondary forests dominate some human‐modified tropical biomes, and this is expected to increase via both abandonment of marginal agricultural land as well as forest and landscape restoration programmes. A key question is whether promoting the recovery and protection of secondary tropical forests will return invertebrate functional diversity and associated functional traits. Dung beetles are ideal for assessing functional diversity as they play vital roles in several ecosystem functions, including seed dispersal, nutrient cycling and bioturbation. We examined how taxonomic and functional diversity, and the functional trait composition of native dung beetle species recovers in naturally regenerating secondary forests in comparison to both cattle pastures and primary forest in the Colombian Choco‐Andes, a global hotspot of threatened biodiversity. Using a space‐for‐time approach, we found that taxonomic and functional diversity recovered to levels comparable to primary forest within approximately 30 years of secondary forest regrowth. Functional richness and FD, measures of the diversity of traits present in a community, were similar in secondary and primary forest, but significantly lower in pasture. Rolling dung beetle species were positively associated with forest habitats, particularly primary, while dwelling species were more common in pasture. Thus, the functional trait composition of secondary forests was more similar to primary forest than to pasture. The ability of secondary forests to rapidly accumulate primary‐forest dung beetle functional diversity, and a representative suite of functional traits, provides an opportunity to protect biodiversity and ecosystem functioning, especially in regions where marginal agricultural land allows cost‐effective conservation actions.


Introduction
Agricultural expansion is driving tropical land-use change, resulting in the conversion of over 150 million hectares of tropical forest between 1980 and 2012 (Gibbs et al., 2010;Hansen et al., 2013). This habitat loss and subsequent fragmentation effects are the leading cause of tropical biodiversity decline . Extensive reductions in species richness (Gibson et al., 2011) and the replacement of forest specialist species with habitat generalists (Clavel et al., 2010) are driving large-scale biotic homogenization (Socolar et al., 2017).
While primary forests remain under significant threat in some regions, secondary forests have become dominant features of human-modified tropical landscapes in others. For example 36.2 million hectares of secondary forest regrew between 2000 and 2010 in Latin America and the Caribbean (Aide et al., 2013), especially in the tropical Andes, Brazilian Caatinga and Costa Rica (Nanni et al., 2019). Farmland is abandoned due to complex socioeconomic and biophysical drivers, especially steep topography and related agricultural marginality, climate, declining rural populations and urbanization (Lambin & Meyfroidt, 2010;Nanni et al., 2019), allowing secondary forests to naturally regenerate.
This trend of land abandonment may be expected to continue. Forest and Landscape Restoration (FLR) is a central component of an integrated programme of interventions to limit global warming to 2°C by growing trees in degraded landscapes . Under the Bonn Challenge and New York Declaration, nations have agreed to restore 350 million hectares by 2030 using FLR, a significant component of which will be via natural forest regeneration in the tropics. In combination, this offers great promise for conservation since secondary forests can recover significant amounts of carbon and biodiversity in relatively short time periods (Gilroy et al., 2014;Poorter et al., 2016;Lennox et al., 2018), including species of conservation concern (Gilroy et al., 2014;Basham et al., 2016). Promoting natural forest regeneration on marginal agricultural land offers a cost-effective opportunity to protect carbon and biodiversity through carbon-based payments for ecosystem services (e.g. REDD+), for instance in the tropical Andes (Gilroy et al., 2014).
Biodiversity loss in the tropics is often assessed using species richness-based measures of diversity and metrics of species composition and turnover. However, these metrics do not account for the differential role of individual species in an ecosystem and so may underestimate true biodiversity loss (Mouillot et al., 2013). Preserving a diversity of species' life histories and functional traits is important for maintaining ecosystem functioning and resilience (Cadotte et al., 2011). Changes in environmental conditions following landuse change can act as a filter, altering the composition and reducing the diversity of traits present in a community (Gray et al., 2007;Cardinale et al., 2012). Growing recognition of this problem has led to alternative measures of biodiversity being used to better assess the impacts of land-use change on functional composition.
Functional diversity (FD) quantifies the range of functional traits and ecological roles present in a community (Petchey & Gaston, 2002;Vill eger et al., 2008). The loss of FD is predicted to lead to ecosystem destabilization (Bregman et al., 2016) and declines in ecosystem service provision Cardinale et al., 2012). Therefore, in combination with understanding of how habitat change impacts the relative abundance of different functional traits (e.g. Edwards et al., 2013a;Cannon et al., 2019), FD is important in predicting the effects of future land-use management on ecosystem functioning. For example maintaining functionally diverse communities of ground beetles and bees is vital for pollination and natural pest control to safeguard future food production (Woodcock et al., 2013).
Conversion of natural habitats to agricultural land drives declines in FD. For example forest conversion to oil palm or pasture in Borneo and Colombia reduces the FD of dung beetles (Edwards et al., 2013a) and birds (Edwards et al., 2013b;Prescott et al., 2016;Chapman et al., 2018;Cannon et al., 2019). In turn, a pan-tropical analysis of avian responses to naturally regenerating secondary forests (tõ 100 years) revealed recovery of forest specialist species richness, functional divergence and functional dispersion over time to primary forest levels (Sayer et al., 2017). However, how secondary forest recovery impacts FD of other taxa is not well understood, which is especially critical in the context of invertebrates (Nichols et al., 2008;Manning et al., 2016). The FD of ground-foraging ants in lowland Brazilian Atlantic forest increased with time since abandonment of buffalo pastures (to~100 years; Bihn et al., 2010), but how this compares to true primary forest controls is unknown. Reforestation via active tree planting of pasture in Queensland, Australia, led to increases in dung beetle species and functional richness with concurrent return of ecosystem functioning (Derh e et al., 2016). Thus, how the recovery of naturally (passively) regenerating secondary forest impacts the FD of invertebrates and the abundance of their different functional traits remains a major unanswered question.
We fill this key knowledge gap by assessing the extent to which native dung beetle taxonomic (TD), functional diversity and abundance of associated individual functional traits, recover in naturally regenerating secondary forests of the Colombian Andes. Using a large-scale dataset spanning three regions, we compare dung beetle communities in secondary forests of different ages to that of cattle pasture and primary forest. The Colombian Andes are a threatened hotspot of global biodiversity (Myers et al., 2000) and represent a costeffective opportunity for gains in naturally regenerating secondary forests (Gilroy et al., 2014;Nanni et al., 2019). Dung beetles are an ideal taxon for assessing functional recovery, because they perform vital ecosystem functions, including seed dispersal, nutrient cycling and bioturbation (Nichols et al., 2008;Manning et al., 2016), are good indicators of change in other taxonomic groups, in particular mammals (Barlow et al., 2007;Nichols et al., 2009;Edwards et al., 2014), are sensitive to environmental change , and are taxonomically well-described (Spector, 2006).
Using this system, we predict three key hypotheses: (1) Dung beetle TD and FD will be greater in secondary forest and primary forest compared to cattle pasture; (2) dung beetle TD and FD will recover with increasing secondary forest age and (3) dung beetle communities in secondary forest will exhibit an assembly of functional traits more similar to that of primary forest than of cattle pasture.

Study area and dung beetle sampling
Three areas in the departments of Antioquia, Risaralda and Choc o in Colombia were sampled. The sites span an altitudinal range of 1290-2680 m above sea level, typified by subtropical and submontane cloud forest (Armenteras et al., 2003). Cattle farming is the dominant land-use in the region, with 95% of farmland devoted to cattle, a trend followed across the Colombian Andes (Etter et al., 2006). All sample sites were situated along the meeting point of agricultural land and large patches of contiguous forest (>1 000 000 ha), mostly comprising primary forest with some secondary forest cover (3-35 years old). All secondary forest points in this study were in relatively close proximity to, and connected with, contiguous primary forest. All sites are characterized by the same broad floristic habitat (Western Cordillera cloud forest).
Traps were placed within 400 9 400 m squares across the three sites, with squares allocated in proportion to habitat types, 38 in forest (23 primary, 15 secondary) and 20 in pasture, with some squares (n = 4 of 58; 7%) straddling habitat types. In Antioquia, we placed 9 pasture, 6 primary and 2 secondary squares, in Risaralda, 5 pasture, 9 primary and 7 secondary squares and in Choc o, 6 pasture, 8 primary and 6 secondary squares. Secondary forest ages were taken from Gilroy et al. (2014), and were obtained from a combination of records from local, land-owning NGOs and interviews with local people. A minimum of 300 m was left between squares in different habitats and 400 m for squares within the same habitat.
Clusters of five sampling traps were placed inside each square, with a minimum of 100 m left between sampling points to ensure community independence (Larsen & Forsyth, 2005). All secondary forest sample points were ≥60 m from the forest edge (45 were ≥100 m from the edge). Baited pitfall traps were used to sample dung beetles, with a total of 180 traps placed. Traps were baited with fresh human dung, which is known to attract the majority of species. Traps were collected every 24 h over a 4-day period, with dung replaced after 2 days. Plastic pint cups were used to create the pitfall traps. Traps were buried with the rim of the cup level to the ground, with cups partially filled with water and scent-free washing up liquid to immobilize trapped insects. Specimens were deposited in the Instituto Alexander von Humboldt, Colombia. Sampling was carried out in the relative dry period, from January to March and June to July 2012.

Functional traits
Six functional traits were analysed: body size, front leg area, front to rear leg ratio, behavioural guild, diel activity and diet range. ImageJ was used to take measurements of body size (length [base of head to elytra base] 9 width [of elytra]), front leg area (front femur area + front tibia area) and front to rear leg ratio ((front femur length + front tibia length)/(rear tibia length + rear femur length + rear spur length)) using photos of a subset of sampled individuals (n = 1-27). Measurements from multiple individuals were then used to calculate mean values for each species for the three traits. Trait information for each species' behavioural guild, diel activity and diet range were all obtained from the literature (Table S1). Where species-specific information was not available, we assumed that traits were common across a genus.

Taxonomic diversity
All analyses were performed in R version 3.5.3 (R Core Team, 2019). Species diversity was calculated using the Shannon-Weiner index, with evenness determined using Pielou's evenness index, with both calculated in the vegan package (Oksanen et al., 2011).

Functional diversity
For functional analysis within habitats, points were grouped together by habitat within each 400 9 400 m square. Calculation of functional indices requires the number of species (S) in a community to be greater than the number of axes (in this case S > 2). At the trap level, this condition was not met by 42 points and so grouping by habitat within a square was necessary. Squares that still did not have sufficient species to enable calculation of functional indices after grouping were dropped from the analysis (n = 3 of 58; 5%). Our functional analyses used two axes to enable the maximum amount of data to be used; increasing the number of axes would result in more habitats within squares not meeting the condition of more species than axes and thus being dropped from the analysis.
We assess three complimentary functional indices: functional richness (FRic) and functional evenness (FEve) based on the hypervolume concept (Vill eger et al., 2008); and dendrogram-based functional diversity (FD) (Petchey & Gaston, 2002). For the hypervolume indices, traits act as coordinates in functional space, identifying the species' functional niche (Vill eger et al., 2008). FRic is a measure of the volume of space occupied by constituent species and FEve describes the distribution of species' abundances within occupied functional space (Vill eger et al., 2008). All traits were equally weighted by abundance. This was carried out using the dbFD function in the FD package (Lalibert e et al., 2014).
The dendrogram-based functional index, FD, is the sum of all branch lengths of a functional dendrogram that connects all constituent species of a community (Petchey & Gaston, 2002). Analysis was carried out using the picante package (Kembel et al., 2010). We also considered additional indices of functional diversity (sesFD, which controls for the confounding impact of species richness on FD, FDiv and FDis; see Methods S1 for details). Moran's I was used to test for spatial autocorrelation in our results, implemented using the ape package in R (Paradis et al., 2004).

Comparing taxonomic and functional diversity between habitats
To compare taxonomic and functional diversity between habitats (cattle pasture, primary and secondary forest), linear mixed-effect models (LME) with maximum likelihood estimation (created using the lme4 package (Bates et al., 2014) were employed, with habitat and altitude as fixed effects, and site included as a random effect. Likelihood ratio tests (LRT) were performed to compare null models that excluded the fixed effect of habitat to the full models. For metrics for which the full model was the best fit (the full model had the lowest AIC value; Table S2), post hoc Tukey tests were performed using the multcomp package (Hothorn et al., 2008).

Comparing taxonomic and functional diversity across secondary forest age
To compare taxonomic and functional diversity over secondary forest age, LME with maximum likelihood estimations were created, with age and altitude as fixed effects, and site as a random effect. LRT were completed, comparing the full model to a null model with the fixed effect of age removed. Age was log transformed to normalize model residuals. Functional analysis was employed at the trap level, meaning points with too few species (S ≤ 2) were removed from analysis (n = 6).

Impact of habitat and altitude on functional traits
We used complementary methods to assess the association of environment and traits across our environmental gradient. First, we used RLQ ordination to identify the principal relationships between environmental variables (i.e. habitat type and altitude) and species' functional traits with reference to species adundances, using the ade4 package (Chessel et al., 2004). RLQ uses three matrices: species x trait matrix (Q), sample point x species abundance matrix (L) and a sample point x environmental variables matrix (R), to create a fourth matrix of traits x environmental variables (Dol edec et al., 1996). Second, we test for significant bivariate environment trait associations (i.e. relationships between an individual environmental variable and individual trait) using the fourthcorner method (Dray et al. 2014). For this we ran a permutation test (with 9999 permutations) using the model 6 combined approach that allows the simultaneous testing of both model 2 (permutation of sites) and model 4 (permutation of species), while also adjusting P values using the false discovery rate method.

Comparing taxonomic and functional diversity across habitats
A total of 17 686 individuals of 27 different species were recorded across all habitats. Primary forest had the greatest abundance of individuals (9750), followed by secondary forest (7351), with cattle pasture having the lowest number of individuals (585). Overall species richness was greatest in secondary (23 species) and primary (20 species) forest, with the lowest observed in cattle pasture (11 species). Given this, trap-level species richness (LME; X 2 = 72.509, d.f. = 2, P < 0.001) and abundance (X 2 = 40.076, d.f. = 2, P < 0.001) of primary and secondary forest was significantly greater than pasture, whereas there was no difference between primary and secondary forests (Fig. 1). Secondary forest taxonomic diversity and evenness did not differ from primary forest, but was significantly greater than in pasture ( Fig. 1; X 2 = 13.738, d.f. = 2, P < 0.01; evenness, X 2 = 15.236, d.f. = 2, P < 0.001).

Impact of habitat and altitude on functional trait composition
The RLQ ordination revealed RLQ axis 1 to be the principal driver of observed patterns, associated with primary habitats and a lack of species that have a dwelling nesting strategy (Fig. 3). The global RLQ permutation test showed that permutation Model 2 was significant (P = 0.0001), whereas Model 4 was not (P = 0.0766). This reveals an overall weak global relationship between species traits and environmental variables, indicating that associations were determined across individual traps (Model 2), but not across species (Model 4). Testing the direct links between RLQ axes and species traits revealed a positive association between RLQ axis 1 and mean leg ratio (Fig. 4). No environmental variables had significant associations (P > 0.07). RLQ and fourth corner analysis, focusing solely on secondary forest, found no evidence of a relationship between habitat variables (age, site and altitude) and species traits (Fig. S3).

Recovery of functional diversity and traits
Our results mirror those of a pan-tropical study on birds which found that forest specialist species richness, functional dispersion and functional divergence were similar in 100 year old secondary forests and primary forests (Sayer et al., 2017). In our study, FD increased with secondary forest age, but there was no effect of age on FRic, suggesting that secondary forests rapidly accumulate a greater diversity of functional traits than found in pasture. FD metrics are sensitive to trait selection (Petchey & Gaston, 2002). All traits we selected have established functional significance, relating to how beetles use resources and the amount and diversity of resources used (Table S1).
Both secondary and primary forests in our study had greater taxonomic and functional diversity than did pastures, supporting findings from previous studies on the impacts of forest loss on functional diversity Edwards et al., 2013a;Cannon et al., 2019) and biodiversity more generally (Barlow et al., 2007;Gibson et al., 2011;Edwards et al., 2014). Secondary forests could also facilitate the dispersal of functionally important species between previously isolated patches of forest (Kormann et al., 2016). Additionally, the regeneration of secondary forest may also increase ecosystem functioning within adjacent farmland through the spillover of pollinator and biological control species (Blanche et al., 2006;Karp et al., 2013). Our results thus reveal the importance of protecting forested habitats from conversion to pasture (Gibson et al., 2011;Laurance et al., 2014).
Functional evenness (FEve) did not differ between primary forest and the other habitat types, although pasture had greater FEve than secondary forest. Reduced FEve suggests lower resource-use efficiency by dung beetles within secondary forests compared to pasture (Mason et al., 2005). However, FEve is positively related to disturbance (Pakeman, 2011), and so highly disturbed sites (e.g. pasture) may have high FEve, whereas sites with less disturbance may have lower FEve as competition is more important in structuring communities.
Increased leg ratio is associated with roller and dweller nesting species. Rolling dung beetle species were positively associated with both forest habitats, especially primary forest, while they were almost entirely absent from pasture (Fig. 3). Dwellers were however more common in pasture than in either forest habitats (Fig. 3). Functional trait composition of primary and secondary forests was thus very similar (Fig. 3b). The recovery of rolling species in secondary forests is likely due to decreased soil temperature compared to pasture (Senior et al., 2017), owing to greater canopy cover, and an increase in the structure of the leaf litter layer, which in combination increases the survival rate of roller larvae (Larsen, 2012). The recovery of roller species in secondary forests is of particular functional importance as they play a vital role in distributing seeds and nutrients away from concentrated piles of dung (Nichols et al., 2008).

Environmental drivers of functional recovery and study caveats
Dung beetles are very sensitive to environmental changes , meaning community assemblage is strongly influenced by forest structure (Halffter & Arellano, 2002;Edwards et al., 2017). Secondary forest recovers microhabitats and favourable microclimates in the tropical Andes (Gonzalez del Pliego et al., 2016). This may explain the ability of sensitive, forest dung beetle species to recolonize secondary forests (Gilroy et al., 2014) and the associated recovery of functional diversity. More widely, dung beetles are good indicators of the presence of other taxonomic groups, particularly mammals given their reliance on dung as a nesting and feeding resource (Nichols et al., 2009). Therefore, the recovery of dung beetle functional diversity in secondary forest suggests a wider strengthening of ecosystem resilience and functioning in these habitats (Nichols et al., 2008).  Table S3.
Dung beetle community responses to land-use change vary with geographic location and altitude (Nichols et al., 2007). Higher-elevation dung beetle communities, such as those we have studied, tend to have a greater physiological tolerance to microclimatic changes than those from the lowlands (Escobar 2005;Ghalambor et al. 2006). Therefore, higher-elevation communities might be more able to recolonize secondary forests than lowland communities, possibly explaining the reduced species richness of dung beetles in lowland secondary versus primary forests in Brazil (Barlow et al., 2007). Such variability across elevation emphasizes the need for more geographically extensive studies of the value of secondary forests for dung beetle biodiversity and, more widely, for other taxonomic groups.
There are two key caveats in our study. First, all secondary forest habitats sampled were adjacent to primary forest, which presumably represent sources of individuals for recolonization . Therefore, more isolated patches of secondary forests may have reduced rates of taxonomic and functional diversity recovery, as many forest specialists are unable to cross the agricultural matrix (Feer & Hingrat, 2005;Larsen et al., 2008). Nonetheless, most secondary forest regeneration in the tropics occurs in close proximity to primary forest (Crk et al., 2009;Sloan et al., 2016), suggesting that our focus on secondary regrowth that is adjacent to contiguous primary forest yields broadly applicable results. Second, secondary forest populations may be sinks with below-replacement population growth, which are reliant on immigration from primary forest sources . Source-sink dynamics could thus erroneously enhance the perceived biological value of our secondary forests. Future research focusing on the effect of landscape configuration on the biological value of patches of secondary regrowth is a valuable next step.

Management recommendations
Our results demonstrate the strong potential for functional diversity recovery and associated conservation gains if pastures are abandoned and forests allowed to naturally regenerate. Secondary forests recover a diversity of functional traits comparable to those found within primary forest, strengthening ecosystem resilience, improving ecosystem functionality and ensuring the provision of ecosystem services. This offers great conservation promise given the expected increase in the extent of secondary forest cover via further land abandonment and FLR programmes. The low profitability of marginal agricultural land in the Tropical Andes (Gilroy et al., 2014) and elsewhere (Morton et al. in press), combined with high rates of land abandonment, suggest that these regions likely represent strong opportunities to promote low-cost forest regrowth. With additional carbon sequestration benefits (Gilroy et al., 2014;Poorter et al., 2016;Lennox et al., 2018), promotion of natural forest regrowth offers an attractive opportunity for conservation to recover and protect high levels of species and functional diversity.
Naturally regenerating forests tend, however, to be poorly protected. Laws, policies and socioeconomic conditions can frequently work against their long-term persistence (Reid et al., 2018). In Costa Rica, for example the laws that protect forests exclude young, regenerating sites; in fact, they are often targeted for clearing to prevent their being reclassified as forest and then legally protected (Sierra & Russman, 2006). In post-peace settlement Colombia, we may expect increased urban-to-rural migration as people reclaim land lost during the conflict and an expansion of alternative economic activities (e.g. mining), which in combination may lead to loss of secondary forests (Baptiste et al., 2017). We conclude therefore by highlighting that an urgent policy focus is needed on the legal underpinnings of forest regeneration and its subsequent longer term protection, supported by prioritization exercises to highlight particularly important areas of secondary forest for conservation action. we thank staff at the Instituto Alexander von Humboldt at Villa De Leyva, Colombia, in particular Claudia Alejandra Medina Uribe and F.Forero, for access to their extensive collection and logistical support, and A. Gonz alez and J. Stephens-Cardenas for dung beetle identification. Funding was provided to T.H. and D.P.E. by the Research Council of Norway (grant number 208836). This is publication #18 of the Biodiversity, Agriculture and Conservation in Colombia (Biodiversidad, Agricultura, y Conservaci on en Colombia [BACC]) project.
Authors' contributions R.W.D, D.P.E. and F.A.E conceived and designed the study questions and methodology. F.A.E collected and processed dung beetle sample data; R.W.D and F.A.E. produced the functional trait matrix; R.W.D led the writing of the paper. F.A.E and D.P.E. contributed critically to the drafts and gave final approval for publication.

Supporting information
Additional supporting information may be found online in the Supporting Information section at the end of the article.
Methods S1. Description of methods used to assess FDiv, FDis and sesFD of data. Figure S1. Map of the study area. Figure S2. FDis, FDiv and sesFD against habitat type with regression of secondary forest age. Figure S3. RLQ ordination of secondary forest habitat variables and traits. Table S1. List of functional traits with information on calculation and functional significance.