Composition and diversity of butterflies (Lepidoptera, Papilionoidea) along an atmospheric pollution gradient in the Monterrey Metropolitan Area, Mexico

Abstract This study compares the variation of richness, abundance and diversity of butterfly species along an atmospheric pollution gradient and during different seasons in the Monterrey Metropolitan Area, Mexico. Likewise, we analyse the influence of environmental variables on the abundance and richness of butterfly species and quantify the indicator species for each atmospheric pollution category. Based on spatial analysis of the main atmospheric pollutants and the vegetation cover conditions, four permanent sampling sites were delimited. The sampling was carried out monthly in each of the sites using aerial entomological nets and ten Van Someren-Rydon traps during May 2018 to April 2019. A total of 8,570 specimens belonging to six families and 209 species were collected. Both species richness and abundance were significantly different between all sites, except for the comparison between the moderate contamination site and the high contamination site; diversity decreased significantly with increasing levels of contamination. The seasonality effect was absent on species richness; however, for species abundance the differences between dry season and rainy season were significant in each site excepting the moderate contamination site. Regarding diversity, the seasonal effect showed different distribution patterns according to each order. Relative humidity, vegetation cover and three pollution variables were highly correlated with both abundance and species richness. From the total number of species found, only 47 had a significant indicator value. This study constitutes the first faunistic contribution of butterflies as indicators of the environmental quality of urban areas in Mexico, which will help in the development of strategies for the management, planning and conservation of urban biodiversity.


Introduction
The dynamics of the demographic growth that cities are facing represents a serious threat to the environment, to the health and quality of life of its inhabitants (Vlahov and Galea 2002). The excessive exploitation of the natural resources, land use changes, industrial and urban concentrations and the large quantity of pollutants being emitted to the atmosphere, damage the environment in a process that seems irreversible (García et al. 2013). These effects not only harm living beings, but also generate phenomena that affect the ecosystem (López et al. 2001). This unregulated urbanisation process and ecosystems degradation occurs more rapidly in countries located in regions classified as developing economies, particularly in Latin America, where it is estimated that 75% of the population lives in cities (UN-HABITAT 2010).
In Mexico, atmospheric pollution has deteriorated air quality in various cities, including the Valle de México Metropolitan Zone, the Guadalajara Metropolitan Zone and the Monterrey Metropolitan Zone (MMZ) (García et al. 2012;Cerón et al. 2014;Mancilla et al. 2015;Menchaca et al. 2015). It should be noted that the main problem is the perception of society which often does not realise the severity of the problem, mainly because there is no clear awareness of pollutant emissions, their concentrations and the damage they cause to health, urban infrastructure, ecosystems and biodiversity (Lezama 2010). The State of Nuevo León, in the northeast of Mexico, has an unregulated urban growth. Its main city, the Monterrey Metropolitan Zone, presents serious environmental problems: geological and hydrological risks, water scarcity, green areas loss, air pollution, amongst many others (Cantú et al. 2013;Badillo et al. 2015;Orta et al. 2016;Sanchez-Castillo et al. 2016;Sisto et al. 2016;Sanchez-Castillo et al. 2017).
Studies on species diversity in urban ecosystems with air pollution problems are necessary to understand the effect of anthropogenic development on the integrity and livelihood of the ecosystem (Mukherjee et al. 2015). However, arthropods in such environments are poorly studied despite being crucial components and indicators of urban ecosystems and biodiversity (McIntyre 2000;Magle et al. 2012;Bonebrake and Cooper 2014). Butterflies, in general, are very sensitive to changes in temperature, humidity and solar radiation produced by disturbances in their habitat, for which the inventory of their communities, through measures of diversity and richness, represents a valid tool for evaluating the state of alteration between an urbanised area with pollution problems and a natural environment (Kremen 1993;Wagner et al. 2003; García et al. 2007;Settele et al. 2008;Li et al. 2009). In this way, these insects act as biological indicators, which can reflect the state of the biota in terms of biodiversity, its variation along gradients, endemism or the degree of human intervention, including air pollution (Fagua 2001;Moreno et al. 2007;Butchart et al. 2010;Defra 2016). Various studies demonstrate that butterfly species richness decreases as urbanisation increases (Blair and Launer 1997;Blair 1999;Hardy and Dennis 1999;Brown and Freitas 2002;Di Mauro et al. 2007; Konvička and Kadlec 2011;Bonebrake and Cooper 2014;Ramírez and MacGregor 2016), not only because the construction of buildings and roads replaces or reduces the area of natural and semi-natural habitats, but because the quality of residual habitats is affected by various forms of pollution (Corke 1998;Hardy and Dennis 1999;Mulder et al. 2005;Jones and Leather 2012;Philips et al. 2017).
Some other studies have yielded results that support the intermediate disturbance hypothesis, in which species diversity peaked in areas with an intermediate level of habitat alteration (Dial and Roughgarden 1998;Niell 2001;Giuliano et al. 2004;Jones and Leather 2012). On the other hand, most insect studies investigating seasonality are based on relatively pristine ecosystems and few have examined the relationship between urban ecosystems seasonality and their butterfly assemblages. Only very few studies have explored the seasonality of urban butterflies, including Brown and Freitas (2002) in Sao Paulo, Brazil;Shapiro (2002) in California, USA;  in Singapore; Chowdhury et al. (2017) in Dhaka, Bangladesh; Gupta et al. (2019) in Delhi, India. The objectives presented here have the caveat that urban gradient studies are clearly a simplification of the complex patterns produced by urbanisation, such as pollution (Alberti et al. 2001;Hahs and McDonnell 2006;McKinney 2008). Therefore, the objectives of this study are: (1) Identify the butterfly species richness in the Monterrey Metropolitan Zone, Mexico; (2) Compare the variation in richness, abundance and diversity of butterfly species along an atmospheric pollution gradient and during the seasons of the year; (3) Analyse the influence of environmental variables (atmospheric pollutants, temperature, relative humidity, solar radiation and vegetation cover) on the abundance and richness of butterfly species; and (4) quantify the indicator value of the species per atmospheric pollution category.

Study area
The Monterrey Metropolitan Zone (MMZ) is the largest urban area in northeast Mexico and the third largest urban centre in the country, extending from 25°15' to 26°30' of north latitude and from 99°40' to 101°10' of west longitude ( Figure 1A, B). The area is limited by the coastal plain of the Gulf of Mexico and the Sierra Madre Oriental mountain range. The MMZ urban sprawl integrates the Municipality of Monterrey in the central portion, the Municipalities of Guadalupe, San Nicolas de los Garza and San Pedro Garza Garcia in the pericentral portion, Apodaca, Escobedo and Santa Catarina Municipalities in the periphery and El Carmen, Garcia, Santiago, Juarez, Cadereyta and Salinas Victoria in the surrounding area (Alanís 2005;González et al. 2011;Mancilla et al. 2015;Ybáñez and Barboza 2017). The MMZ has a vehicle fleet of 1.7 million (INEGI 2010) and 4.1 million of inhabitants (INEGI 2011), which is probably higher nowadays. Likewise, there is a variety of industrial complexes that include the production of glass, steel, cement and paper, amongst others (Menchaca et al. 2015). The centre of the city has an average altitude of 540 m a.s.l, its steppe climate is dry and warm with temperatures above 35 °C during the summer and below 8 °C during the winter (Alanís 2005;González et al. 2011;Menchaca et al. 2015).

Delimitation of the pollution and vegetation cover gradients
Since November 1992, the MMZ operates a network of air quality monitoring stations known as the Integral Environmental Monitoring System (SIMA). The SIMA network currently consists of 13 registration stations, located following the criteria of meteorological, epidemiological, land use and population density studies. The concentrations registered in the monitoring stations are: PM 10 (particulate matter of less than 10 µm), PM 2.5 (particulate matter of less than 2.5 µm), carbon monoxide (CO), ozone (O 3 ), nitrogen dioxide (NO 2 ), nitrogen oxides (NO x ) and sulphur dioxide (SO 2 ). In addition, some meteorological variables are reported, such as barometric pressure (Bp), rainfall (R), relative humidity (Rh), solar radiation (Sr), temperature (T) and the direction (Wd) and magnitude of wind (Ws) (Arreola and González 1999;González et al. 2011;Mancilla et al. 2015). The data recorded by SIMA stations for air quality and meteorological variables during the period from 2008 to 2017 were obtained from the website of the National System of Air Quality Information (SINAICA). Descriptive measures for each of the months and each registered year were obtained in the Statistica 13.3 software (TIBCO Software Inc. 2017).
To identify the main air quality descriptor pollutants in the MMZ during the period 2008-2017, a Principal Component Analysis (PCA) was carried out. Subsequently, to differentiate the changes in the spatial distribution of the air quality indicator pollutants in the MMZ, maps were prepared using the annual average information per monitoring station. The creation of maps was carried out using Inverse Distance Weighting Interpolation (IDW), with a value of 2 as Coefficient of Distance and the pixel size of the output raster re-defined to 10 metres. As reference of the extension for each interpolation, the minimum and maximum distances were taken from the vector sections corresponding to the urban areas that form the MMZ; such vectors were obtained from the National Land Use and Vegetation Series 6 layer (INEGI 2016). The procedures described above were performed using the QGis 3.2 software (QGIS Development Team 2018). As a result, three categories of air pollution were generated: low (0.19 to 26.51 µg/m 3 of NO x ; 3.22 to 10.56 µg/m 3 of PM 2.5 ), moderate (26.51 to 52.83 µg/m 3 of NO x ; 10.56 to 17.92 µg/m 3 of PM 2.5 ) and elevated (52.83 to 79.19 µg/m 3 of NO x ; 17.92 to 25.3 µg/m 3 of PM 2.5 ) ( Figure 1C). Percentage of vegetation cover was determined through an analysis of MODIS images for the period 2008-2017, obtained from GIOVANNI online server. Consequently, three categories of vegetation cover were designated: low (23 to 40%), moderate (40 to 57%) and high (57 to 74%) ( Figure 1D).

Selection of sampling sites
Four permanent sampling sites were delimited, based on the spatial overlapping of four geographic elements: (1) the IDW analysis of the main atmospheric pollutants ( Figure 1C), (2) the vegetation cover ( Figure 1D

Sampling and processing of specimens
Monthly samplings were carried out for each of the sites, during the period from May 2018 to April 2019, resulting in a total of six samplings per season: dry season (November, December, January, February, March and April) and rainy season (May, June, July, August, September and October). The seasons were defined based on historical data of monthly total values of temperature and rain (average from 2008 to 2017), which were obtained from the SIMA stations located within the study area ( Figure 2). Therefore, a total of 48 samplings were considered (six samplings for two seasons and four sites).
The sampling of individuals was carried out using an entomological aerial net. In each of the sites, tours were made inside a pre-established quadrant of 150 m × 150 m, following the techniques recommended by Villarreal et al. (2006). The sampling time at each site was nine hours in the period from 08:00 h to 17:00 h. Furthermore, along with the use of aerial nets, Van Someren-Rydon traps were used (Rydon 1964). Ten traps were placed, five at one end of the quadrant and five at the opposite end, at a distance of 30 m from each other and between 1 and 2.5 m from the ground, with an exposure time of nine hours (08:00 h to 17:00 h). The bait used for the traps consisted of a fermented mixture of seasonal fruits: banana (Musa paradisiaca), pineapple (Ananas comosus), mango (Mangifera indica) and guava (Psidium guajava).   The collected specimens were mounted according to the described procedure of Andrade et al. (2013). For taxonomic identification of specimens, the works of Scott (1986), Llorente et al. (1997), Luis et al. (2003), Garwood and Lehman (2005), Vargas et al. (2008), Luis et al. (2010) and Glassberg (2018) were consulted. The interactive list and the phylogenetic ordering of Warren et al. (2012) were taken as reference. All the specimens were labelled and deposited in the entomological collection of the Conservation Department of the Forestry Sciences Faculty at the Autonomous University of Nuevo León, Linares, Nuevo León, Mexico.

Data analysis
The observed species richness was measured as the total number of species in the study area, as well as at each of the sites. The constancy index was determined and the species were classified as: constant (species found more than 50% of the time during sampling), accessory (species present between 25 to 50%) and accidental (species in less than 25%) (Sackis and Morais 2008). Significant differences in the number of species between sites were determined using a generalised linear model (GLM) (they were modelled as Poisson-distributed variables with a log link) and one-way ANOVA, with Statistica software 13.3. The sampling efficiency was calculated for the entire study area and for each study site using the interpolation and extrapolation methodology proposed by Chao and Jost (2012), available in the package iNEXT .
Five species categories were considered according to the total registered abundance: rare (species with one individual), scarce (from 2 to 5), frequent (from 6 to 21), common (from 22 to 81) and abundant (with 82 or more individuals) (Luna et al. 2010). Differences in the abundance of butterfly communities at the sites were calculated with a GLM and one-way ANOVA. For the analysis of alpha diversity, we adopted the analytical method of Chao and Jost (2015) to obtain diversity profiles in which diversity is evaluated in terms of "effective numbers of species" (qD), an approach that is equivalent to Hill's numbers (Hill 1973). The analysis was made for all the study area and for each study site using the 3.5.3 version of R (R Development Core Team 2019), with the package SpadeR . To examine the differences in species composition amongst the four sites, we performed a non-metric multidimensional scaling analysis (NMDS), based on the similarity matrix using the Bray-Curtis Index. A one-way PERMANOVA was also performed to test for differences in species composition between sites. Both analyses were performed using R 3.5.3 with the package Vegan (Oksanen et al. 2019).
The seasonal effect was measured separately, comparing the species richness, abundance and diversity observed per study site during the rainy season (May to October 2018) and dry season (November of 2018 to April 2019). The indexes and statistical tests mentioned above were used for such comparisons: GLM and nested ANOVA tests for differences in species richness and abundance, estimation of species richness and alpha diversity index, which were performed in Statistica 13.3 and R 3.5.3. In addition, a two-way PERMANOVA and NMDS analyses were carried out, to include the seasonal effect in the species composition, with the aim of grouping sites and seasons. These analyses were performed in Statistica 13.3 and R 3.5.3. A Canonical Correlation Test was applied between the community parameters (number of species and abundance) and the different environmental variables: monthly averages of the main variables of atmospheric pollution (with the highest loading values in the PCA previously obtained) (NO 2 , NO x , CO and PM 2.5 ), climatological variables (temperature, relative humidity and solar radiation) and vegetation cover variables extracted from the SIMA stations and from MODIS images nearest to the sampling sites, using Statistica 13.3.
Finally, to calculate the association value of each butterfly species with the habitat type, the Indicator Value Index, or IndVal, was used (Dufrêne and Legendre 1997). This index is based on the degree of specificity (exclusivity of the species to a particular site, based on its abundance) and the degree of fidelity (frequency of occurrence within the same habitat) (Tejeda et al. 2008), expressed in a percentage value. The analysis was performed with the labsdv package on platform R 3.5.3, using 1,000 random permutations to define the level of significance. The indicator species with an index equal to or greater than 75% were classified as "characteristics", which are defined by their high specificity for a given habitat; species with a value less than 75%, but equal to or greater than 50% are considered as "detectors", which show different degrees of preference for diverse habitats (McGeoch et al. 2002).
The sample coverage estimator indicated that our inventory for the MMZ is 99.7% complete. In Figure 3, we plotted the proposed diversity with the method of Chao and Jost (2015) and the confidence intervals for q = 0, 1 and 2 for each pollution category and station.
Both species richness and abundance were significantly different (p < 0.05) between all sites, except for the comparison between Site 3 (moderate pollution) and Site 4 (high pollution) (Table 2, Figure 4A, B). Both parameters (species richness and abundance) decreased with increasing levels of pollution. In Site 1 (pollution free) 2,683 individuals and 199 species were registered, representing a sampling coverage of 99.7%. In Site 2 (low pollution), the values were reduced to 2,334 individuals and   For 0 D, 1 D and 2 D, Site 1 (pollution free) had the highest diversity. All comparisons between sites were significantly different (with 95% confidence intervals) ( Figure 3A, B and C). The one-way PERMANOVA test detected significant differences in species composition amongst Sites 1, 2 and 3 (free, low and moderate contamination) (SS total = 6.63; SS within-group = 5.05; F = 4.575, p < 0.001). Butterflies collected during each month for Sites 1, 2 and 3 (free, low and moderate contamination) formed separate groups in the NMDS diagram (Stress = 0.23) ( Figure 5A).

Effect of seasonality on butterfly changes per pollution category
In the MMZ, the highest abundance was registered during the rainy season, while the highest species richness was shown during the dry season. The highest completeness of the inventory was recorded during the rainy season ( Table 2).
The differences in the abundance of the dry and rainy season were significant (p < 0.05) in each site with the exception of Site 3 (moderate contamination) (Table 3, Figure 4D). Regarding diversity, the seasonal effect was absent in Site 1 (pollution free) for 0 D and 1 D, while for 2 D, Site 1 (pollution free), 2 (low contamination) and 4 (high contamination) did not show significant differences between seasons ( Figure 3D, F). Two-way PERMANOVA allowed us to identify a significant effect of season (F = 7.702, df = 1, p < 0.001) and pollution (site) (F = 5.682, df = 3, p < 0.001) on species composition, as well as an interaction effect between the two factors (F = 2.315, df = 3, p < 0.001). Butterflies collected each month in each of the study sites formed separated groups by seasons in the NMDS ordination diagram (Stress = 0.23) ( Figure 5B).
The variables NO 2 , NO x , PM 2.5 , relative humidity and vegetation cover were highly correlated, both with abundance and species richness. The individual interaction of NO 2 , NO x , PM 2.5 and relative humidity with abundance and species richness resulted negative, while the interaction of vegetation cover with both ecological parameters resulted positive (Table 4).    (Clench, 1946) 18  (Table 5).

Discussion
This study constitutes the first faunistic contribution of butterflies as indicators of the environmental quality of an urban area in Mexico and the first inventory of butterflies systematically carried out in the State of Nuevo León. The 209 registered species in the MMZ constitute 78.6% of the described richness so far in Nuevo Leon, according to Luz and Madero (2011) in collaboration with the North American Butterfly Association (NABA) and 10.2% to what was recorded for Mexico (Warren 2000;Llorente et al. 2006). However, it should be noted that urban gradient studies are clearly a simplification of the complex patterns produced by urbanisation (Alberti et al. 2001;Hahs and McDonnell 2006;McKinney 2008). The specific impacts of urbanisation on species richness vary, depending on variables, such as geographic location and many historical and economic factors that are unique to each city (McKinney 2008). Urbanisation intensity is correlated with increased disturbance and structural simplification of the remaining vegetation through landscaping practices that remove woody plants, leaf litter and other microhabitats from natural communities (Marzluff and Ewing 2001). Combination of all these factors reduces the area and quality of animal habitat and tends to increase with the intensity of urbanisation (Alberti et al. 2001;Hahs and McDonnell 2006). Studies using spatial gradients in urban areas have shown that the development of these areas can strongly and negatively affect many sensitive butterfly species (Blair and Launer 1997;Blair 1999;Clark et al. 2007). The close relationship between the abundance of host plants and the persistence status of butterflies suggests that the decline of plants may cause the co-extinction of some associated butterflies  or the host plant's own rarity (rather than decrease) could be associated with another trait of the butterfly that makes it vulnerable to extirpation (Harrison 1991). Corke (1999), Mulder et al. (2005) and Öckinger et al. (2006) suggest that the butterfly decline in urban areas could be a secondary effect of heavy metal stress presence on local plants, not resulting in a decrease in the number of host-plants, but in a selective pressure of pollutants on the plant vigour, subsequently affecting their associated fauna. Either way, our results corroborate similar studies of declining butterfly populations, suggesting that habitat degradation may be a devastating threat to the persistence of certain sensitive taxa, such as butterflies characteristic of unpolluted and low-pollution sites (Schultz and Dlugosch 1999;Weiss 1999;Wagner and Van Driesche 2010;Bonebrake and Cooper 2014).
A study of butterfly communities in fragments of urban forests in Brazil (Brown and Freitas 2002) similarly found that the most important factors affecting diversity and composition, excluding site size and sampling time, were connectivity, vegetation, flowers and negative human impact, such as pollution (indirect effect of urbanisation).
Observations of butterfly diversity provide information on variations in species richness and abundance formed by vegetation throughout the landscape and the interaction between species (Öckinger et al. 2009). Although local determinants of diversity, such as competition and predation, remained undetermined in these studies, to a large extent, landscape characteristics influence butterfly richness and abundance in different geographic areas (Öckinger et al. 2006; Öckinger et al. 2009). Spatial scale differ-ences in butterfly diversity can be attributed to heterogeneity at the landscape level, while timescale differences can be attributed to changes in climatic conditions at both local and regional scales (Mukherjee et al. 2015). In the current context, it can be assumed that butterfly diversity varies in the four sampling sites as a matter of differences in pollutants concentration and composition of vegetation.
In general, the differences in species distribution in the four areas were prominent, although the abundance of the different species was not profound (indirect effect of the high degree of urbanisation), possibly due to the high concentrations of the contaminants, as well as the corresponding abundance of host plants in the affected areas. The observed variations in species richness in areas without apparent pollution provide an impression of differences in food plants abundance and landscape characteristics in the region. Previous studies on the diversity of butterflies in landscapes with high pollution in contrast to the regions of moderate and low pollution show that the richness increased with the availability of green space and the heterogeneity of habitats in terms of the available plant species and dominant microenvironmental conditions (Kuussaari et al. 2007). According to these studies, our observations register a greater diversity in the areas with no apparent pollution and low urbanisation, followed by the areas of low, moderate and high pollution and urbanisation (Blair and Launer 1997;Kitahara and Sei 2001;Hogsden and Hutchinson 2004).
Regardless of variations between different landscapes, observations of butterfly diversity in the study area suggest that conservation management is necessary to ensure the livelihood of the different ecosystem services derived from butterflies. The abundance of butterflies in urban landscapes will promote the pollination and hence the propagation of different plant species that can reduce the decrease in vegetation, consequently diminishing other variables, such as noise and mainly pollution levels (Mukherjee et al. 2015;Selmi et al. 2016;Alfie and Salinas 2017). To understand the ways these insects respond to urbanisation, different authors have suggested three ecological patterns that are related to our results: (1) there are fewer butterfly species in highly urbanised areas (Knapp et al. 2008;Soga et al. 2014); (2) the number of specialised butterflies decreases with increasing urbanisation, a case demonstrated by the number of indicator species in each gradient category (Bergerot et al. 2011;Lizée et al. 2011;Koike 2012, 2013); and (3) urbanisation can lead to local disappearance of rare and not abundant, specialised butterfly species, as shown in this study (Fattorini 2011;Soga and Koike 2012).
We found that the variables associated with the increase in urbanisation (NO 2 , NO x and PM 2.5 ) were negatively correlated with the richness of butterflies, while the measures associated with less developed areas (green space) were positively correlated. These results are consistent with those of Ruszczyk (1986), Ruszczyk and DeAraujo (1992) and Stefanescu et al. (2004), who found lower diversity of species with higher urbanisation. Growing more trees and shrubs that provide nectar for adult butterflies and using a greater variety of larval food plants species in planting schemes, may be more effective in maintaining populations of these insects than simply increasing the amount of plant cover cultivated .
The richness and distribution of butterfly species fluctuates according to their life cycle, which is linked to seasonal changes. However, compared to butterflies in temperate climates, seasonal variation generally does not have a great impact on tropical butterflies, which are reported as well distributed throughout the year, the case corresponding to the present study, as there is no seasonal differentiation for most of the comparisons (Hamer et al. 2005). In the MMZ, the butterflies showed a highest species richness during the dry season with 88% of the total species observed during the evaluation period. This finding is contrary to other studies, which report higher numbers in the rainy season (Devries et al. 1997;Hamer et al. 2005;Hernández-Mejía et al. 2008;Meléndez-Jaramillo et al. 2019). This can be attributed to human intervention, when irrigation provides higher food resources for butterflies and also attracts a higher number of species than normal during the dry season.
Biodiversity inventories provide crucial reference information for future ecological and conservation studies. The existence of species lists at various stages of the urbanisation process allows documentation of changes in species composition over time. However, few lists of butterfly species have been published in cities, most of which are restricted to few countries, for example, Brasil, Argentina or India (Núñez 2008;Chowdhury and Soren 2011;Silva et al. 2012). Until now, significant efforts have been made in and around cities to conserve endangered butterfly species (Daniels 2009;Ramírez and MacGregor 2016). Butterfly conservation in urban areas is a feasible task, since many species can thrive in these areas. Hopefully, creative urban planning and management, such as habitat design and planting of native, nectar-rich plants, could improve urban habitats for butterflies. However, all actions must be monitored and must build on prior knowledge about the biology and ecology of the target species to be successful (Snep et al. 2006;Kadlec et al. 2008).

Conclusions
For the first time in Mexico, butterflies were systematically sampled in order to monitor the environmental quality in an urban area. A total of 8,570 specimens belonging to six families, 19 subfamilies, 31 tribes, 138 genera and 209 species of butterflies were collected for the study area. The highest species abundance and richness, as well as alpha diversity, are recorded at the site free from air pollution, that is associated with a less impacted landscape. Both species richness and abundance were significantly different between all sites, except for the comparison between the moderate contamination site and the high contamination site, while diversity decreased significantly with increasing levels of contamination. The overall trend of distribution of butterflies to the levels of air pollution shown in the Monterrey Metropolitan Area is a decrease, this being in agreement with the general disturbance hypothesis.
The seasonality effect was absent on species richness; however, for species abundance, the differences between dry season and rainy season were significant in each site, excepting the moderate contamination site. Regarding diversity, the seasonal effect showed different distribution patterns according to each order. The variables NO 2 , NO x , PM 2.5 , relative humidity and vegetation cover, were highly correlated, both with species abundance and richness, so they could be the main reasons for the variation of butterfly communities in this study.
This work is one of the first studies of butterflies in a specific area of northeast Mexico, in which the environmental quality and seasonality in an urban area were analysed. The information presented here provides benchmarks that allow the comparison of the diversity and richness of Papilionoidea species at regional and national levels. This information can be used as an initial step to analyse the possible use of butterflies as an indicator group of the biodiversity in Mexico.