Introduction

At a global and local scale, changes in land use impose enormous pressures on biodiversity (^{Foley et al., 2005}; ^{Newbold, Hudson, Hill, Contu, & Lysenko, 2015}). Three variables linked with changes in land use are agricultural conversion, deforestation and human footprint (^{Meyer & Turner, 1992}). In the case of Mexico, the country has experienced a large change in land use (^{González-Abraham et al., 2015}, ^{Mas et al., 2004}; ^{Robson and Berkes, 2011}). For example, Mexican agricultural areas have dramatically increased over the last 20 years (^{Hansen et al., 2013}). At the same time, per capita forest area decreased by almost a third between 1980 and 2000 (^{Mas et al., 2004}; ^{Velázquez et al., 2002}). The human footprint in Mexico has been large and heterogeneous, and is explained by the complex physical geography, historic human settlements, technology use, and recent demographic explosion (^{González-Abraham et al., 2015}). Although it is intuitive that land use change, tree cover reduction and increased human footprint are driving Mexican organic diversity to extinction, precise measurements or even a rough approximation of this is lacking.

Land use change and deforestation have been shown to affect the biodiversity of many taxa (^{Barragán, Moreno, Escobar, Halffter, & Navarrete, 2011}; ^{Brown, 1997}; ^{Schulze et al., 2004}; ^{Vidal, López-García, & Rendón-Salinas, 2014}). However, such relationships are not necessarily straightforward and linear. For example, agricultural areas may drive reduced species richness but such effects can be more acute with the prevalence of cultivated species, intensive land management and short cultivated cycles (^{Scales & Marsdens, 2008}). In contrast, take the case of some coffee plantations which, despite altering vegetation structure of a community, still provide protection to a number of native species (^{Martínez et al., 2005}, ^{Perfecto and Vandermeer, 2015}; ^{Perfecto, Rice, Greenberg, & Van der Voort, 1996}). The impacts of land use change and deforestation on biodiversity can be taxa-dependent as some taxonomic groups are more resilient than others (^{Hansen et al., 2001}). For example, birds are commonly thought to be highly sensitive to land use change. However, some particular bird taxa are more likely than others to go extinct (^{Fischer et al., 2007}). In this case, the simplicity of landscape structure may render some species to become more likely to disappear if such species do not find their vital resources (^{Fischer et al., 2007}).

Conservation biologists have recently come to consensus that beyond species, there is a need to assess community level interactions (^{Rodewald, Rohr, Fortuna, & Bascompte, 2014}; ^{Tilianakis, Laliberté, Nielsen, & Bascompte, 2010}). Ecological interactions between species can affect community composition and dynamics (e.g. ^{Faust & Raes, 2012}; ^{Thébault & Fontaine, 2010}; ^{Verdú & Valiente-Banuet, 2008}). A classic example is that of microbial communities in which the presence or absence of certain species will affect the interactions of all species (^{Rodríguez-Martínez & Pascual, 2006}). This implies that communities occur non-randomly and that, in conservation terms, the community structure and dynamics should be assessed. There are several methods to systematically characterize and study ecological networks (^{Pascual & Dunne, 2005}). Among these methods, co-occurrence networks provide a valuable tool to infer statistical relationships among species or other taxonomic units in relatively large datasets, allowing one to go beyond composition and abundance and examine community structure (^{Deng et al., 2012}, ^{Faust and Raes, 2012}). Two of the structural parameters that have been proposed to reflect the relevance of particular nodes in an ecological network: (1) the number of connections a focal node has with other nodes, and (2) a focal node's betweenness centrality (the proportion of the shortest paths between all i and j pair of nodes that contains the focal node) (^{Jordan, 2009}; ^{Jordán, Liu, & Davis, 2006}). Many ecological networks exhibit modular organization, meaning that their elements tend to form semi-autonomous groups with more interactions within them than with other groups. Such a modular structure suggests that ecological interactions lead to the formation of well-defined groups of species, which could in turn be relevant with respect to conservation and ecosystem stability (^{Pascual and Dunne, 2005}, ^{Thébault and Fontaine, 2010}).

One relatively neglected functional group in conservation biology are complex life cycle (CLC) organisms. CLC animals are those that face an abrupt habitat shift as they go through different developmental stages and so they may be exposed to radically different ecological conditions and stressors (^{Werner, 1988}). For example, many animals with CLC spend their larval stage being aquatic while the adult is terrestrial. In terms of conservation strategies, CLC animals present a more complicated scenario for their protection as their different transitions (e.g. aquatic and terrestrial) need to be considered both for the stressors involved and their habitat protection. Odonate insects (dragonfly and damselflies) are an example of a CLC animal that has gained increased conservation attention (e.g. ^{Bried, Hassall, Simaika, Corser, & Ware, 2015}; ^{Clausnitzer et al., 2009}; ^{Kalkman et al., 2008}; ^{Samways, McGeoch, & New, 2013}). Calls for conservation efforts directed toward odonates focus on the transformation of the body waters they use as larvae and their aerial surroundings used by flying adults (^{Bried and Samways, 2015}, ^{Corbet, 1999}, ^{New, 2009}). The relationship between land use and odonate species richness shows mixed results. On one hand, it has been recorded that forest areas of agricultural areas may contribute to odonate richness as in the case of European farms (^{Sahlén, 2006}; ^{Ruggiero, Céréghino, Figuerola, Marty, & Angélibert, 2008}), but this is not necessarily so for tropical areas (e.g. Costa Rica; ^{Hofhansl & Schneeweihs, 2008}). Indeed, habitat changes affect populations of odonates very quickly and dramatically (^{Harabis & Dolny, 2011}). For example, ^{Junior et al. (2015)} observed that urbanization led to changes in odonate species composition and richness. Habitat quality (i.e. physical structure) affected odonate assemblages. ^{Oliveira-Junior et al. (2015)} found that moderately-degraded areas may host high species richness. ^{Kadoya, Akasaka, Aoki, and Takamura (2011)} conducted a heuristic analysis of contributing factors affecting odonate diversity and found that although water quality conditions largely explain species losses, land use change can also contribute (^{Kadoya et al., 2011}). Note, however, that many of these analyses have relied on simplified land use changes such as forest conversion to farms. Thus, further analyses are needed to examine how general is the relation between land use (at a broad scale level) and odonate diversity. Moreover, such analyses can be enriched if the component of community-level interactions is added. Some authors have begun such a community-level analyses: odonate communities have been related to habitat types (e.g. using canonical correspondence analysis; ^{Clausnitzer, 2003}, ^{Stewart and Samways, 1998}). Although valuable, these efforts do not provide information with regard to the role of species in maintaining community stability as the co-occurrence networks we explained above.

In this paper we examine how deforestation and land use affect odonate species richness in the Mexican territory. This was done using incidence of odonate geographical records according to the different land uses reported in Mexico from the year 2000 to 2014. As a second aim, we have analyzed the co-occurrence of species presence in the Mexican territory as a preliminary attempts at a community-level analysis.

Materials and methods

We searched for odonate species records in Mexico from 2000 to 2014 using ^{Conabio (2015)} databases, published literature and unpublished data. Records were depurated and georeferenced, if needed (note that 1 geographical species record consists of only 1 species, but 1 geographical coordinate may contain more than 1 species). We considered presence-absence species data as we did not take into account the abundance of total species records of each species in 1 geographical coordinate. This abundance data may have provided several useful resources to our study such as, for example, how human impact may have impacted species richness at the population level. Odonate species records were counted (even if there were several for the same species) per geographic coordinate and assigned to the corresponding land use entities as referred by the Mexican Land Use and Vegetation map, series V (Layer Union) (^{Inegi, 2011}) (Table 1). In each land use entity, we counted all species and geographic records. This was taken as species richness and presence of records respectively.

A grid of 15 × 15 km approximately for each pixel size was used. This grid size was used as the best logistical compromise between collection effort and the probability of finding a collected odonate. For example, a smaller grid size could have been used but given that there are few collections, richness could not have been calculated. Each pixel with at least 1 species was considered. All odonate species that fell inside each pixel were counted. The total species number observed per pixel was considered as a measure of species richness (herein, SOBS). Given that we were unaware of how complete richness was per pixel, we also calculated the Chao2 Index (^{Chao, 1987}; CHAO2). Both calculations were made with DIVA-GIS v.7.5 (^{Hijmans, Guarino, Cruz, & Rojas, 2002}).

To assess the relationship between odonate richness with tree cover presence, we obtained the tree coverage according to Global Forest Change (^{Hansen et al., 2013}). Trees were defined as vegetation taller than 5 meters, and so we used the series 2000 percent tree cover (^{Hansen et al., 2013}). Human footprint was assessed as the influence index (HII) from Socioeconomic Data and Applications Center (SEDAC). HII provides an updated view of anthropogenic environmental impact, considering human population density, human land and infrastructure, and human access (roads, railroads, navigable rivers or coastlines) (WCS, ^{CIESIN, 2005}). In each grid pixel, we calculated the average value of tree cover percentage (AVGTC) and the human influence index (AVGHII) with the ArcMap extension of ArcGIS v.9.3 (^{ESRI, 2008}). Then, average values of AVGTC and AVGHII were regressed between the SOBS and CHAO2 respectively. We searched for the best fit of a model regression in each case. Analyses were made using Minitab v.17 software (Minitab Inc.).

From the odonate species records, we generated a file containing a list of all species in the records and their count in each of the land-use types. This file was processed using CoNet (http://psbweb05.psb.ugent.be/conet/), a Cytoscape 3.2.1 application that has been used to infer statistical relationships among species or other taxonomical units within communities (^{Faust and Raes, 2012}, ^{Faust et al., 2012}). In order to infer the significant co-occurrence links among Odonata species in the records, we set CoNet to use the count mode and the Pearson correlation index with a threshold of 0.9. Co-occurrence links were graphed using Cytoscape 3.2.1. This rendered a network in which nodes represent odonate species and edges stand for significant co-occurrence relationships. Connectivity and betweenness centrality of the nodes were calculated using the Cytoscape 3.2.1 built-in network analyzer. Finally, the network modules (potential communities) were identified with the Cytoscape 3.2.1 application ModuLand (^{Szalay-Bekö et al., 2012}).

Results

We found 4,071 species records from 275 odonate species (information available upon request). Species richness and records were mainly from forest/rainforest and agricultural areas which accounted for approximately 60%, with around 30% of all species records for each area. Forest in general encompassed 23.4% of total richness, followed by agricultural areas with 15.9% of the whole species total, and urban areas with 10.5%.

Our results suggest that odonate species and records occurring in the Mexican territory are scarce in areas without vegetation and with grassland-rainforest (Table 1, Figure 1, Figure 2). Odonates collected in urban areas represented 7.7% of records and 10.5% of species richness. Desert areas contributed with 7.6% of the species richness (6.1% of the records). There was a significant correlation between odonate records and species richness per land use (Pearson r = 0.987, p < 0.00001). The highest SOBS and CHAO2 values were found mostly for pixels located in Michoacán, Nayarit, San Luis Potosí and Veracruz states, while the lowest values were found in Sonora, Chihuahua, Coahuila and Zacatecas states with only 1 species per pixel (Figure 3, Figure 4).

Figure 1 Land use, vegetation entity, odonate species richness (black bars) and geographical records (gray bars) in the Mexican territory 2. 

Figure 2 Odonate species records for the Mexican territory. Collection range spanned from 2000 to 2014. 

Figure 3 Observed odonate species (SOBS) in a 15 × 15 km (approximately) pixel size for the Mexican territory. Green pixels indicate from 1 to 23 species, yellow pixels indicate 21–26 species, orange pixels indicate 27–40 species, while red pixels indicate 40–65 species. 

Figure 4 Number of odonate observed spe cies (SOBS) (a) and Chao 2 (CHAO2) (b) regressed against tree cover average (AVGTC). 

A positive linear relation between SOBS and AVGTC (n = 406, linear R² = 0.05, p < 0.005) was found. CHAO2 showed a similar tendency (n = 406, R² = 0.0139, p = 0.017) (Fig. 4a, b). The regression model with the best fit was quadratic for both SOBS (n = 406, quadratic R² = 0.0437, p = 0.001; Fig. 5a) and CHAO2 (n = 406, quadratic R² = 0.0139, p < 0.05; Fig. 5b) versus AVGHII.

Figure 5 Number of observed species (SOBS) (a) and Chao 2 (CHAO2) (b) regressed against human influence index (AVGHII). 

Regarding the analysis at the community level, the co-occurrence network generated from the records suggest that Odonata species indeed tended to associate in well-defined modules or communities (Fig. 6a). Moreover, this analysis uncovers some species as nodes with a particular high connectivity or betweenness centrality, pointing at potentially relevant species in terms of the conservation of whole communities (Fig. 6b; Table 2).

Figure 6 Structure of the co-occurrence network inferred from the odonate species counts in the different land-use types. Nodes represent species and edges correspond to the statistically significant co-occurrence links. Colors in (a) represent the different modules that were identified and that could correspond to different communities. In (b), the size of each node is proportional to its betweenness centrality, while its color corresponds to how densely connected a node's neighborhood is. 

Discussion

The relatively low number of odonate records in this analysis may affect our results and thus our discussion below has to be framed within this situation. However, waiting to increase sampling effort in such a big country may be meaningless for Mexico given that the number of specialists has not augmented in the last 3 decades (^{González-Soriano & Novelo-Gutiérrez, 2014}). This suggests that sampling effort will not increase.

In relation to the land use, our results indicated (not surprisingly) a prominent role of forest or tropical forest areas in preserving odonate presence (i.e. over one third of records were collected in these areas). It is well-known that oligiotrophic waters (with macrophyte presence) and good riparian vegetation coverage are essential components for the odonate life cycle (e.g. ^{Remsburg and Turner, 2009}, ^{Samways and Steytler, 1996}, ^{Stewart and Samways, 1998}). These components can be present in forested and even agricultural water bodies. Coherent with this, a substantial portion of the records (30%) and 20% of all species were located in agricultural areas. This indicates that agricultural areas are also important reservoirs for Mexican odonate species. This finding echoes previous odonate studies in other countries. For example, in lowland Pacific Costa Rica, species’ accumulation curves showed highest species richness in forest margin habitats followed by agricultural areas (^{Hofhansl & Schneeweihs, 2008}). Not all agricultural practices impact species richness negatively (e.g. ^{Hofmann and Mason, 2005}, ^{Lenat, 1984}). The precise mechanism(s) that explains why cultivated areas are key for odonate species richness and communities in the Mexican scenario is unclear. We speculate that conditions that affect larval and adult survival in aquatic and terrestrial places play a key role in explaining this pattern. With such a large number of odonate species in our study, describing such conditions is an impossible task. Agricultural area sampling was more intense in our database than forested area sampling: 60% of records are from agricultural areas while only 6% of records are from forested areas. Rather than odonates being less common in these areas, an alternative explanation is that forested areas with combined land use have been less collected. A bias in collecting effort is actually one problem in our analysis. This may explain, for example, why we found a positive relationship between odonate records and species richness. As a matter of fact, odonate collection in the country of Mexico has rarely been systematically planned with the exception of research dissertations (^{González-Soriano & Novelo-Gutiérrez, 2014}). This may have resulted in biased sampling such as collecting near roads or populated areas (^{Cuevas-Yáñez, Rivas, Muñoz, & Córdoba Aguilar, 2015}), a phenomenon that has been also described for other insects such as butterflies (^{New, 2009}). Surprisingly, richness was fairly high in urban areas. Urban aquatic areas are stressed sites where odonate richness is usually impaired. At the general level, 3 sources of stress can negatively affect odonate diversity: landscape fragmentation (e.g. ^{Chovanec et al., 2000}), vegetation loss/or modification (^{Jeanmougin, Leprieur, Loís, & Clergeau, 2014}) and sewage discharge and/or stormwater runoff (^{Henriques-Oliveira, Baptista, & Nessimian, 2007}). Despite these stressors, a prominent number of odonate species have been found to occur in urban areas (reviewed by ^{Villalobos-Jiménez, Dunn, & Hassall, 2016}). These odonates include generalists, specialists and even threatened species. This indicates: (a) a prominent role for urban aquatic ecosystems, and (b) a large level of tolerance to urbanization by odonates. For our study, one case is that of Mexico City which despite being one of the largest cities in the world, it has a fairly high diversity (^{González-Soriano & Novelo-Gutiérrez, 2007}). One reason is that some areas have remained safe to the stressors above despite being within the urban boundary.

Our results also indicated some geographical areas of particular importance for odonate conservation: including regions in Michoacán, Nayarit, San Luis Potosí and Veracruz states. ^{González-Soriano, Noguera, and Oñate (2011)} suggested the Huasteca region of San Luis Potosí as an odonate hot spot which concurs with our findings. Clearly, these areas deserve immediate attention in terms of conservation. Sadly, studies that assess the odonate conservation status in the country of Mexico have been rare (e.g. ^{Cuevas-Yáñez et al., 2015}). This, along with current Mexican legislation which does not include any odonate as threatened (NOM-059-Semarnat-2010; ^{Semarnat, 2010}), implies that efforts should be directed at working with local and federal governments to outline an action plan toward odonate conservation assessment.

A weak but positive relationship between tree cover and odonate species richness was also observed. With regards to how novel these results are, scarce information exists on tree cover and dragonfly richness in tropical areas. In fact, ^{Clausnitzer et al. (2009)} highlighted that the highest odonate richness are in running waters in tropical forests partly arguing in favor of the role of the tree cover effect. Mexican neotropical forest regions have highly complex habitats. Forest complexity provides odonates habitats where they may shelter at night, or safely secure against adverse weather conditions (e.g. rain or strong wind) or predation. If farm sites are near forest, the latter can function as a refuge or buffer that allows establishment of odonate species (^{Paulson, 2006}). Although it is true that shaded conditions have a negative effect on odonate presence (^{Remsburg & Turner, 2009}), forest cover still plays an important role in almost every part of the life cycle of the average odonate. Mechanisms that mediate odonate survival in general can be manifold – thus there is not a sole role for landscape biotic integrity. For example, in Japan, more than half of the odonate species depend on its forest, although at a broad-scale, climate largely explains species distribution (^{Tsubaki & Tsuji, 2006}). Despite this, the Mexican scenario for odonate conservation should emphasize the maintenance of tree cover even after native forests are transformed. One example is the case of coffee plantations which retain canopy integrity and, in combination with native forest, may become key refuges for several local plant and animal species (^{Martínez et al., 2005}, ^{Perfecto and Vandermeer, 2015}, ^{Perfecto et al., 1996}). Actually, odonate species that are currently in danger, e.g. Paraphlebia zoe, are likely to be protected thanks to such combination of coffee plantation-original cloud forests in Veracruz and Hidalgo states in Mexico (^{Cuevas-Yáñez et al., 2015}).

Our analysis showed a quadratic relationship between human influence and odonate richness and richness completeness. Nonlinear regressions have already described for odonate community structure (e.g. Southern High Plains in the USA; ^{Hernández, Reece, & McIntyre, 2006}). Notice, however, that our regression coefficients were low. This suggests that other variables not considered in our study are involved in predicting odonate richness. The effects of human activity on odonate populations have been documented in lentic (standing) waters in temperate zones (^{Findlay & Houlahan, 1997}; ^{Hamasaki et al., 2009}; ^{Hall, McCauley, & Fortin, 2015}; ^{Raebel et al., 2012}), but rarely in tropical zones and lotic systems (^{Monteiro-Junior, Juen, & Hamada, 2015}). Anthropogenic disturbance in the form of, for example, road construction, may have significant effects on wetlands biodiversity (^{Findlay & Houlahan, 1997}). In our analysis, the highest richness was in sites with an intermediate human footprint. This tendency has been also observed in odonates of the Amazon basin (^{Junior et al., 2015}). ^{Stewart and Samways (1998)} also noticed that high dragonfly species richness occurs in streams with moderate disturbance, although the species that live in such places are not endangered or threatened. This phenomenon may be explained by a high metapopulation succession dynamics in sites altered by anthropogenic factors or that such high dynamics are weighted by less complex adaptations in habitat selection or may function as refuges for species survival (^{Harabis & Dolny, 2011}). Another explanation is that collection and sampling may be biased by the collectors’ preferences for places close to urban areas (e.g. aquatic places close to roads) to avoid the challenge implied by areas of difficult access.

The co-occurrence network suggested that odonate species are organized in non-random communities and that some species may reinforce each other's presence. While this network offers a first look at the potential community organization of Odonata species, it is important to bear in mind that it represents statistical relationships and that more in-depth studies will be necessary to uncover the direct or indirect interactions underlying co-occurrences. Related to this, among the future challenges is the functional or taxonomical characterization of the modules, or proposed communities, which could help uncover some of the principles behind their assembly. Species with high connectivity or betweenness centrality seem to show very diverse biological properties and/or ecological needs: modules typically (1) include both odonate suborders: Anisoptera (e.g. Macrothemis imitans) and Zygoptera (e.g. Ischnura capreola); (2) include small (e.g. Ischnura ramburi) and big (e.g. Remartinia luteipennis); (3) engage in territorial (e.g. Hetaerina titia) and non-territorial (e.g. Enallagma novahispaniae) male mating behavior; (4) use lotic (e.g. Brechmorhoga tepeaca), and lentic (e.g. Orthemis discolor) waters, and, (5) use open (e.g. Erythemis vesiculosa) but also closed (e.g. Brechmorhoga praecox) habitats. Nevertheless, if anything, these species seem to share a widespread distribution. How this property (or others) may affect community stability deserves further studies.

Finally, what is to be learned in terms of the functional property of odonates being CFC organisms and our results? Given that odonates make use of aquatic and terrestrial areas, we need to pay attention to the different stressors that occur in such distinct places and the preservation of both habitats (^{Stoks & Córdoba-Aguilar, 2010}). On one hand, stressors in the aquatic environment are usually agrochemicals which in many cases increase eutrophication or decrease macrophyte presence (e.g. ^{Kadoya et al., 2011}, ^{Raebel et al., 2012}, ^{Weisner et al., 2007}). On the other hand, stressors such as the reduction of tree cover imply a number of negative actions for odonates such as: fewer perching refuges for resting and predation protection, less available food, increased sun exposure and less humidity. Although the aim of our paper is not to detail a conservation plan, in general the above stressors need to be reduced to preserve odonate presence. Indeed, the challenge is huge.