Contrasting patterns of neutral and functional genetic diversity in stable and disturbed environments

Abstract Genetic structure among and diversity within natural populations is influenced by a combination of ecological and evolutionary processes. These processes can differently influence neutral and functional genetic diversity and also vary according to environmental settings. To investigate the roles of interacting processes as drivers of population‐level genetic diversity in the wild, we compared neutral and functional structure and diversity between 20 Tetrix undulata pygmy grasshopper populations in disturbed and stable habitats. Genetic differentiation was evident among the different populations, but there was no genetic separation between stable and disturbed environments. The incidence of long‐winged phenotypes was higher in disturbed habitats, indicating that these populations were recently established by flight‐capable colonizers. Color morph diversity and dispersion of outlier genetic diversity, estimated using AFLP markers, were higher in disturbed than in stable environments, likely reflecting that color polymorphism and variation in other functionally important traits increase establishment success. Neutral genetic diversity estimated using AFLP markers was lower in disturbed habitats, indicating stronger eroding effects on neutral diversity of genetic drift associated with founding events in disturbed compared to stable habitats. Functional diversity and neutral diversity were negatively correlated across populations, highlighting the utility of outlier loci in genetics studies and reinforcing that estimates of genetic diversity based on neutral markers do not infer evolutionary potential and the ability of populations and species to cope with environmental change.


& Hoffmann, 2006).
Researchers are increasingly aware that the consequences and benefits associated with high intraspecific diversity can in turn profoundly influence processes at higher levels of organization (Bolnick et al., 2011;Des Roches et al., 2018;Whitlock, 2014). A recent meta-analysis (Des Roches et al., 2018) indicates that, in aquatic systems, effects of intraspecific variation on ecological processes and ecosystem services are often comparable to, and sometimes stronger than, effects of species richness and species composition. This upward cascading includes effects of intraspecific diversity on community species composition, in part mediated by consequences of species traits, species interactions, and species sorting. An increased understanding of the causes and consequences of intraspecific diversity in a spatiotemporally changing world is thus essential for several areas in ecology and evolution.
In the present study, we investigate the roles of ecological and evolutionary processes for population genetic structure and as drivers of phenotypic, functional, and neutral genetic diversity, and examine whether and how the relative importance of different processes change depending on ecological settings and conditions.
To that end, we use data for populations of dispersal and color polymorphic pygmy grasshoppers (Tetrix undulata (Sowerby, 1806) Orthoptera, Tetrigidae) in disturbed and stable habitats (Figure 1).
Tetrix undulata is a small, diurnal, ground-dwelling, and widely distributed insect that inhabits biomes ranging from the Mediterranean to arctic regions of Europe (Holst, 1986). It usually occupies drier microhabitats in relatively open areas where it lives on the soil surface and feeds on microalgae growing on moist soils, mosses, and detritus (Holst, 1986). Adult and late instar nymphs hibernate during winter and emerge in April-May when reproduction ensues in our study area ( Figure 1). Females survive at most one reproductive season, produce multiple pods of egg (<35 eggs/clutch), and nymphs develop through five (males) or six (females) instars before eclosing. Because of its high reproductive capacity, T. undulata may rapidly become very numerous when conditions are favorable (Forsman, 2018).
Like many other insects, pygmy grasshoppers display dispersal polymorphism; a short-winged flight-incapable morph that lacks flight muscles coexists with a long-winged flight-capable morph (Berggren, Tinnert, & Forsman, 2012;Forsman, 2018;Harrison, 1980;Tinnert & Forsman, 2017) (Figure 1a). The discrepancy in dispersal capability among wing morphs generates a dynamic process, and a high initial incidence of long-winged individuals in recently founded populations is followed by a decline in frequency over time (Berggren et al., 2012;Forsman, 2018). Populations in disturbed environments are therefore expected to show a high proportion of long-winged phenotypes (Berggren et al., 2012;Denno et al., 1996) combined with low neutral genetic variation due to drift.
Here, we use amplified fragment length polymorphism (AFLP) data (Bensch & Åkesson, 2005;Vos et al., 1995) to investigate genetic structure and diversity in 20 T. undulata populations in southern Sweden (Figure 1). Next, we compare wing morph frequencies, functional phenotypic variability as estimated using data on color morph diversity, and neutral and functional genetic variability based on AFLP loci between populations in disturbed and stable habitats.
We hypothesized that disturbed environments inhabited by recently established populations should be characterized by a high incidence of dispersal phenotypes; high diversity in functional traits because of the enhancing effect of diversity on establishment and persistence (Forsman, 2014(Forsman, , 2016Forsman, Betzholtz, & Franzén, 2016;Forsman & Wennersten, 2016;Mills et al., 2018); and low neutral genetic diversity owing to the eroding effect of drift associated with small founder groups and fluctuating population sizes (Charlesworth & Charlesworth, 2017). The eroding effect of founder events is expected to leave a stronger signature on neutral than on functional genetic diversity in disturbed environments in part because establishment will typically fail if the group of colonizing individuals is not genetically diverse for functional traits . By contrast, we hypothesized that populations in stable environments should be characterized by a lower incidence of dispersal phenotypes, unless much influenced by recent immigration (Harrison, 1980;Roff & Fairbairn, 2007); lower variability in functionally important traits, owing to the variance-reducing effect of stabilizing selection (Arnold & Wade, 1984;Endler, 1986); and higher neutral genetic diversity, under the assumption that populations fluctuate less and that the eroding effect of founder events and bottlenecks therefore is weaker in stable environments (Charlesworth & Charlesworth, 2017).

| Sampling and study areas in contrasting environments
All applicable institutional and/or national guidelines for the care and use of animals were followed. We selected 20 sampling localities so that our dataset contained T. undulata populations with varying degrees of interpopulation distance and habitat stability. Sampling localities represented either undisturbed, relatively stable (n = 8) habitats or disturbed (n = 12) habitats ( Figure 1, Table 1). Habitats in stable environments included stream shorelines, pastures, meadows, and agricultural land, all of which are influenced by small-scale F I G U R E 1 Wing morphs, color morphs, and study area. (a) Tetrix undulata pygmy grasshopper female belonging to the macropterous morph with long and functional wings (left) and micropterous short-winged flightless morph (right). Photograph: J. Tinnert. (b) Tetrix undulata representing reddish brown, striped, striated brown, gray, and black color morphs. Photograph: A. Forsman. (c) Map of study area in the south of Sweden, showing location and Sample ID of 20 Tetrix undulata pygmy grasshopper populations in habitats representing stable (open circles) and disturbed (filled dots) environments. See Table 1 for a key to abbreviations of sampling locations substrate surface modifications through trampling by farm animals, wave action, grazing by cattle and wild birds (e.g., geese and ducks), and agricultural activities. The modifications of soil surface and ground vegetation that result from these activities favor pygmy grasshopper behaviors, such as feeding and egg-laying (Forsman, 2018). In contrast, sampling localities in disturbed environments TA B L E 1 Diversity in samples of Tetrix undulata pygmy grasshoppers collected from 20 sampling locations in stable and disturbed environments in the southeast of Sweden. N indicates number of individuals used for AFLP analyses, census population size is calculated based on information on number of collected individuals and capture probability, color morph unalikeability indicates color morph diversity estimated in a sample of (N) individuals, long-winged indicate the proportion phenotypes with functional wings in a sample of (N) adult individuals, PL indicates number of polymorphic loci, PPL indicates percentage polymorphic loci, and Hj, nHj neutral, and sHj outlier indicate genetic diversity estimated using all 1,419 AFLP loci, 1,208 neutral AFLP loci, and 28 outlier AFLP loci, respectively. Distance to centroid indicates within-group genetic dispersion as estimated using PERMDISP Overall, some species of pygmy grasshoppers seem well adapted to using a tracking strategy suitable for spatiotemporally changing environments and can thrive under these conditions (Forsman, 2018).
Disturbed sites were located, on average, at greater distances from the coast ( Figure 1) and at higher elevation (mean: 150, range: 36-238 m above sea level) compared with stable (mean: 36, range: 4-110) sites (Table 1). It might therefore be hypothesized that some other environmental factor(s) not related to disturbance have influenced the genetic differentiation and diversity of the populations under study. However, our analyses suggest that this is unlikely (see Results).
Grasshoppers were collected in spring and early summer (for details, see Forsman, Karlsson et al., 2011, 2012Tinnert, Hellgren, Lindberg, Koch-Schmidt, & Forsman, 2016;Tinnert & Forsman, 2017). Individuals were identified to species according to Holst (1986), classified according to their sex, wing morph (long-winged with functional wings or short-winged and flightless) (Berggren et al., 2012;Tinnert & Forsman, 2017), and color morph. All or a subset (depending on how many were collected) of the collected individuals were preserved in 90% ethanol until DNA extraction (see below). The number of individuals used for AFLP analyses is typically different from the sample sizes used to calculate the proportion of long-winged individuals and color morph diversity (Table 1).
These discrepancies arose because some of the collected individuals were nymphs, which can be used for AFLP analyses but not for classification of wing or color morph. For some sampling locations, only a subsample of the collected individuals was brought to the laboratory for classification of wing morph, color morph, and DNA extraction, whereas remaining individuals were released at the sampling location.

| Estimating population size and immigration rate
Census population size (Table 1) was estimated for each locality based on information on number of collected individuals adjusted for capture probability (Tinnert & Forsman, 2017), which was set to 0.4 based on previous extensive capture-mark-recapture studies (Berggren et al., 2012;Forsman & Appelqvist, 1999;Tinnert & Forsman, 2017;Tinnert, Hellgren et al., 2016). We searched for grasshoppers while walking slowly through the area during days with weather conditions suitable for grasshopper activity, that is, clear or overcast days with a temperature of at least 15°C (Forsman, Karlsson et al., 2011, 2012. The number of individuals collected at each site underestimates actual population size, but it is a reliable relative measure which is robust to factors such as differences in the area covered for sampling, time invested in sampling, number of people involved in each sampling event, habitat type, and weather conditions that could potentially influence total catch (Tinnert & Forsman, 2017;Tinnert, Hellgren et al., 2016).
The proportion of long-winged phenotypes ( Figure 1a) in each population was used as a proxy for recent colonization or immigration events (Table 1). Previous studies show that the proportion of pygmy grasshopper males and females that belong to the macropterous long-winged morph is highly correlated across samples from different populations and years, that the incidence of the long-winged morph does not differ consistently between males and females (Berggren et al., 2012;Forsman, 2018;Tinnert, Hellgren et al., 2016), and that capture probability is independent of both sex (Forsman & Appelqvist, 1999) and wing morph (Berggren et al., 2012). It is therefore unlikely that the estimates of the proportion of long-winged phenotypes reported in Table 1 were influenced to any important degree by sampling bias according to sex or wing morph, or by any differences in sex ratio among samples from different collection sites. Furthermore, only the macropterous phenotype is able to fly (Berggren et al., 2012). Collectively, this suggests that a high incidence of the long-winged flight-capable morph may be used as a proxy to identify T. undulata populations that have been recently established or that represent older populations that have been much influenced by recent immigration (and that might show signatures associated with admixture) (Tinnert, Berggren, & Forsman, 2016;Tinnert & Forsman, 2017;Tinnert, Hellgren et al., 2016). Available evidence indicates that wing morph in pygmy grasshoppers is heritable and not influenced to any important degree by developmental plasticity (Berggren et al., 2012;Forsman, 2018). Available evidence also indicates that wing morph is independent of color morph (Forsman, 2018;Forsman, Karlsson et al., 2011).

| Estimating color morph diversity
Pygmy grasshopper color morphs range from light gray via different shades of brown to black, some being uniform and others mottled or patterned with longitudinal stripes, vertical bars, or speckles (Ahnesjö & Forsman, 2003Forsman, 2018;Forsman et al., 2002Forsman et al., , 2007 (Figure 1b). Individuals vary also with regard to texture of the integument, the surface being either smooth, granular or consisting of longitudinal ridges and grooves. Split-brood experiments have shown that neither the patterning nor the overall darkness of pattern elements is influenced by substrate, temperature, or crowding (Ahnesjö & Forsman, 2003;Forsman, 2011;Karlsson, Johansson, Caesar, & Forsman, 2009). Further, color morph frequencies in samples of wild-caught individuals from different populations are highly correlated with those in captive-reared individuals, indicative of population-level heritability (Forsman, Karlsson et al., 2011), and there exists as yet no firm evidence that the color polymorphism in pygmy grasshoppers is affected by developmental plasticity (for a detailed discussion, see Karlsson et al., 2009).
To quantify color morph diversity in populations, we used the coefficient of unalikeability (Kader & Perry, 2007;Perry & Kader, 2005). Measures of population variability for categorical characters are different from those used for quantitative or continuous traits where variation about the mean is the norm. The coefficient of unalikeability (u 2 ) is a measure of the proportion of possible comparisons which are unalike and was calculated based on the equation: where p represents the proportion of the population in each of the i th categories of the categorical variable (i.e., color morph) (Kader & Perry, 2007;Perry & Kader, 2005). Unalikeability can take any value from 0.0 to 1.0, with higher values representing higher diversity. Estimates of color morph unalikeability for our study populations were independent of sample size (r = 0.35, n = 20, p = 0.14).

| DNA extraction and molecular genetics analyses
From each of the 20 sampling locations, we used 5-28 individuals for DNA extraction and genetic analysis (Table 1). DNA was extracted from the femur of each individual using phenol-chloroform method according to Sambrook (Sambrook, Fritch, & Maniatis, 2002) (see Supporting Information for details).
Analysis of AFLP was carried out as described previously (Bensch & Åkesson, 2005;Vos et al., 1995), using the restriction enzymes EcoRI and Tru1 pre-amplification primers M A X E A and combinations of four selective primers (pair 1-E TAG X M CGA , pair 2-E TAG X M CAG , pair 3-E TCG X M CAC , and pair 4-E TAG X M CAC ) (Tinnert & Forsman, 2017;Tinnert, Hellgren et al., 2016). Three negative and nine positive controls were included on each plate. The fragment analyses were performed by Uppsala Genome Center with an ABI3730XL DNA Analyzer (Applied Biosystems, USA). Resulting chromatograms were analyzed in GeneMapper 5.0 (Applied Biosystems) for allele calling (see Supporting Information for details). A total of 1,419 polymorphic sites were used for further analysis. Peak heights were normalized using AFLPscore (Whitlock, Hipperson, Mannarelli, Butlin, & Burke, 2008) and were scored as present if the peaks were >15% of the mean peak height of a particular locus. We obtained a resulting binary matrix, consisting of ones (presence) and zeros (absence of the fragment).

| Detection of loci under selection
We used a Bayesian approach implemented in BayeScan v.2.1 (Foll & Gaggiotti, 2008) to detect the outlier loci. The analysis was performed with 20 pilot runs and a 200,000 step burn-in followed by 200,000 iterations (a sample size of 5,000 and thinning interval of 20). Posterior odds (PO; the ratio of posterior probabilities of selection over neutrality) ratio was set to 100. Low polymorphic loci with dominant allele frequency <2% and >98% across the whole dataset were removed to minimize false discovery rate. The original dataset was divided into two subsets based on the results of the BayeScan: outlier loci under selection (sDATA) and neutral loci (nDATA). Loci with log 10 PO > 2 (decisive evidence for selection) were regarded as the outliers (loci under selection). To obtain a neutral dataset, loci with log 10 PO > 1 (strong evidence, Jeffrey, 1961) were removed from the dataset to exclude also potential outliers.

| Estimating genetic diversity within populations
Overall F ST and genetic diversity indices were calculated in AFLPsurv v1.0 (Vekemans, 2002) for all loci, all neutral loci (nDATA), and all outlier markers (sDATA). We used a Bayesian method with nonuniform prior distribution of allele frequencies (Zhivotovsky, 1999). We assumed Hardy-Weinberg equilibrium, as T. undulata is neither highly self-fertilizing nor haploid. Additionally, preliminary inbreeding coef-  (Nei, 1973) for both neutral (nHj) and outlier (oHj) loci. In addition to estimating within-population genetic diversity using Hj, PERMDISP tests performed with PRIMER v.7 (Clarke & Gorley, 2006) with PERMANOVA+ (Anderson, Gorley, & Clarke, 2008) were used to test the null hypothesis of homogeneity of within-group dispersions among populations and between stable and disturbed environments (Anderson, 2006). PERMDISP uses principal coordinates generated from a Jaccard distance matrix (Bonin, Ehrich, & Manel, 2007) to compare the average deviations from centroids (Anderson, 2006), and we did this separately for neutral and outlier loci. If the distances of individuals from one population to its centroid are larger than distances to the corresponding centroid in another population, this suggests that the diversity is greater in the former group (Parker, Anderson, Jenkins, & Brunton, 2012).
PERMDISP is normally used primarily for checking the assumption of homogeneity of variance and determine whether differences in within-group dispersion are likely to confound results and conclusions from PERMANOVA regarding variation in overall genetic structure among populations (Anderson, 2006;Anderson & Walsh, 2013). Here, in addition to that (see the next subsection), we used PERMDISP to test for a difference in the level of genetic diversity between populations in stable and disturbed environments.

| Analyses of genetic structure among populations
It is possible that environmental factors associated with geographic location (coastal versus inland) or with differences in elevation have resulted in genetic differentiation among populations. It is also possible that the existence of disturbed and stable environments has favored the evolution of genetically different ecotypes of pygmy grasshoppers. To evaluate these hypotheses, population genetic structure was investigated using estimates of pairwise F ST values (Weir & Cockerham, 1984) as implemented in ARLEQUIN v.3.5 (Excoffier & Lischer, 2010) using 10,000 permutations.
An analysis of molecular variance AMOVA (Excoffier, Smouse, & Quattro, 1992) was performed in ARLEQUIN to determine the hierarchical genetic structure based on F ST for the nDATA and sDATA datasets. Two separate calculations were performed for each dataset: one with structure and one without structure. For the former, to partition the total genetic variation and test for overall differentiation between grasshopper populations in disturbed and stable environments and according to geographic location or elevation, we performed a nested AMOVA where the populations were nested in either of two disturbance regimes (disturbed or stable). For the latter, all 20 populations were grouped together.
The statistical significance of the AMOVA was tested by 10,000 permutations.
In addition to AMOVA, a multivariate approach was used to investigate the effects of sampling location and environment (disturbed-stable) on the distribution of the genetic variation among T. undulata individuals. Jaccard distances between the individuals (Bonin et al., 2007), separately for neutral and outlier loci, were calculated in FAMD v1.3 (Schluter & Harris, 2006), and the resulting matrices were used for the further statistical analysis that was performed in PRIMER v7. The matrices were visualized using nonmetric multidimensional scaling ordination (MDS) (Kruskal, 1964).
Permutational multivariate analysis of variance (PERMANOVA) (Anderson, 2001;McArdle & Anderson, 2001) was used to test for differences in genetic structures among different sampling locations (populations) and between stable and disturbed environments.
The sampling locations were introduced as a random factor and disturbance regime as a fixed factor. All permutation tests used 9,999 unrestricted permutations and partial sums of squares type (Type III). A significance level of α ≤ 0.05 was accepted. Results from PERMDISP analyses indicated that there were differences in genetic dispersion between populations in disturbed and stable environments (see Results). Results from PERMANOVA pointing to any difference in overall genetic structure between environments should therefore be interpreted with caution (Anderson & Walsh, 2013).
To examine whether genetic structure and pairwise divergence among populations in neutral AFLP loci were associated with the divergence in outlier AFLP loci putatively influenced by selection, we used a Mantel test. Mantel tests were performed in the zt software (Bonnet & Van de Peer, 2002) with 9,999 permutations.

| Statistical analyses
Analyses of genetic structure uncovered significant differentiation among most study populations and no signature of isolation by distance (see Results), justifying that the different study populations can be considered as independent observations. We tested for differences between populations inhabiting disturbed and stable environments in proportion of long-winged phenotypes, color morph diversity (coefficient of unalikeability), neutral and outlier genetic diversity, and census population size using either t tests or nonparametric Kruskal-Wallis tests, depending on whether data violated assumptions of homogeneity of variances or normality. We investigated whether genetic diversity as estimated based on neutral and outlier AFLP data was associated with the level of within-population diversity in a functionally important phenotypic trait (color pattern) using correlation analysis. Prior to statistical analyses, the proportion of long-winged individuals was arcsine-square root-transformed, and census population size was log10-transformed.

| Outlier loci identified using selection test
Of the 1,419 AFLP loci, 1,277 had dominant allele frequencies >2% and <98%, and they were included into the outlier analyses.
According to BayeScan analysis, 28 (of the 1,277 = 2.2%) loci exceeded the threshold for "decisive" evidence of selection with a FDR of ≤0.001043 as the outliers, and they were used for the outlier dataset (sDATA). Alpha values were positive for all 28 loci, suggesting diversifying selection. In addition to these 28 loci, 41 loci that exceeded the threshold of "strong" evidence of selection were excluded from the neutral dataset (nDATA) that included 1,208 loci. These 41 loci also had positive alpha values, except one with negative value suggesting balancing or purifying selection.

| Population genetic variance within populations
Analyses of the 1,419 polymorphic AFLP markers (the dataset ALL) revealed high within-population genetic diversity estimated with proportion polymorphic loci (PPL) and average heterozygosity (Nei's genetic diversity, Hj) across all sampling localities (PPL = 60.5 to 82.5; Hj = 0.23 to 0.35, Table 1). Hj and PPL were correlated across sampling locations (r = 0.84, n = 20, p < 0.0001). Estimates of Nei's genetic diversity and of mean deviations from centroid for neutral (DevCneu) and for outlier loci (DevCout) are reported for each population in Table 1 Estimates of neutral and functional genetic diversity were not strongly or consistently associated with sample sizes (the number of individuals from each site used for AFLP analyses, nHj: r = −0.19, p = 0.42 sHj: r = −0.34, p = 0.14; DevCneu r = 0.13, p = 0.57; DevCout r = 0.28, p = 0.24, all n = 20).

| Population genetic structure
Results of the nested AMOVA did not reveal any significant overall genetic differences depending on disturbance regime (Supporting Information Table S1), indicating that disturbed and stable environments were not populated by grasshoppers that represented different reproductively isolated evolutionary lineages. This lack of difference also demonstrates that populations that inhabited sites that were closer to the east coast and at lower elevation were not genetically separated or different overall from inland populations ( Figure 1c). Sampling location accounted for 9% of the total genetic variation, and 91% was explained by variation among individuals within sampling locations (Supporting Information Table   S1). Comparisons of results of outlier and neutral AFLP loci showed that sampling location accounted for a greater portion of the total variance in functional (32%) than in neutral (8%) genetic variation (Supporting Information Table S1).
Results from PERMANOVA were similar to those reported above for AMOVA, for both neutral and outlier loci. When populations were nested within environments (stable or disturbed), there was no significant signature of environment (neutral: F 1,315 = 1.33, p = 0.065; outlier: F 1,315 = 1.60, p = 0.14), whereas the location effects were strong (neutral: F 18,315 = 3.00, p = 0.0001; outlier:

| Associations of functional phenotypic diversity with genetic diversity, immigration, and census population size
Color morph diversity decreased with increasing neutral genetic diversity across all populations (deviation from centroid: r = −0.51, n = 20, p = 0.022) but was not associated with outlier genetic diversity (r = 0.22, n = 20, p = 0.36).
Color morph diversity increased with increasing immigration rate (as estimated using proportion of long-winged phenotypes, r = 0.46, n = 20, p = 0.04).
F I G U R E 3 Results from contrasting environment comparisons. Comparisons of (a) immigration rate (as estimated using proportion long-winged phenotypes), (b) functional color morph diversity (as estimated using unalikeability), (c) neutral genetic diversity (Hj, as estimated from AFLP data), (d) outlier genetic diversity (Hj, as estimated from AFLP data), (e) genetic dispersion of neutral ALFP loci (estimated using PERMDISP as deviations (mean ± SE) from centroid), and (f) genetic dispersion of outlier AFLP loci between 20 populations of Tetrix undulata pygmy grasshoppers collected from disturbed or from stable environments (see text and Table 1  ics, the well-studied pygmy grasshoppers lend themselves admirably for investigating drivers of population-level genetic diversity (Forsman, 2018;Tinnert & Forsman, 2017). Overall, our present results illustrate how the genetic structure and diversity of natural populations are influenced by a combination of ecological and evolutionary processes, the relative importance of which differ for neutral and functional genetic diversity and vary according to ecological conditions, such as environmental disturbance regime.
Specifically, analyses based on both neutral and outlier AFLP markers revealed significant differentiation among most study populations (174 and 175 of 190 for the neutral and outlier loci, respectively). Results from comparisons between contrasting environments uncovered that populations in disturbed habitats had a higher incidence of long-winged phenotypes, as expected if these populations were recently established by flight-capable colonizers.
Intrapopulation functional color morph diversity was also higher in disturbed than in stable environments. This outcome was expected given that color polymorphism increases establishment success (Forsman, 2014;Forsman et al., 2012;. Neutral genetic diversity varied among populations, and, unlike functional diversity, it was lower (not higher) overall in disturbed than in stable habitats. This pattern indicates that the eroding effects on neutral diversity of genetic drift associated with recent founding events or population bottlenecks were stronger in disturbed compared to stable habitats. Finally, there was a statistically significant negative overall correlation between functional phenotypic and neutral genetic diversity across populations.

| Contrasting environments-diversity in disturbed and stable habitats
The pattern of genetic structure and intrapopulation diversity seen in these grasshoppers can be accounted for by combined effects of ecological and evolutionary processes associated with environmental change, as summarized in Figure 4 and discussed below.
According to literature, T. undulata is an almost exclusively shortwinged species (Holst, 1986). This is in concordance with the incidence of the long-winged morph (0-17%) in populations in stable habitats. The higher proportion of the long-winged dispersive phenotype in disturbed habitats (59% on average) suggests that these populations were founded by individuals carrying genes for longwinged phenotypes shortly after the disturbance event (Berggren et al., 2012). The lack of association across populations of neutral genetic diversity with long-winged phenotypes is also in agreement with the interpretation that the populations in disturbed environments represented new establishments, rather than older populations having been influenced by gene flow and admixture. If the founders of recently established populations had originated from different source populations, one would expect both functional and neutral genetic diversity to be high in disturbed environments (Hedrick, 2006;Simberloff, 2009). However, pygmy grasshopper populations can be established by very few individuals , indicating that founder groups of divergent origins are not necessarily very common. Older populations that are much influenced by immigration might also display high frequencies of long-winged individuals, but this should not be accompanied by low neutral diversity.
Color morph diversity and within-group dispersion of outlier genetic diversity were higher in disturbed than in stable habitats. This is most likely a consequence of the higher survival and establishment success of more color morph diverse pygmy grasshopper founder groups . The positive effects of functional diversity on colonization success, stability, and persistence of populations (Forsman, 2014(Forsman, , 2016Forsman, Betzholtz, & Franzén, 2015;Forsman & Wennersten, 2016;Hughes et al., 2008;Karpestam et al., 2016) are probably sufficiently strong to outbalance the eroding effect of drift. Additionally, opposing selection in males and females (Forsman, 2018;Forsman & Appelqvist, 1999;Karpestam, Merilaita, & Forsman, 2014) together with the promiscuous mating behavior of pygmy grasshoppers Johansson, Caesar, & Forsman, 2013) may prevent the erosion of functional genetic diversity.
The higher color morph diversity and the higher genetic diversity of outlier AFLP loci in disturbed compared with stable F I G U R E 4 A schematic illustration and summary of how random and deterministic processes differently influence neutral and functional genetic diversity of populations in stable and disturbed environments. Green (up) and red (down) thick arrows refer to positive and negative effects, respectively, and indicate that the magnitude of the effect is greater than in the other environmental disturbance regime. Black thin arrows indicate that the effects might be present albeit of relatively small magnitude compared with those indicated by red and green thick arrows. na indicates not applicable environments may be attributable in part to different selection regimes. In disturbed environments, oscillating and temporally changing selection associated with succession and habitat modifications resulting in rapid temporal shifts in morph frequencies (Forsman, Karlsson et al., 2011;Karpestam et al., 2013) may promote and preserve genetic variation. In stable environments, by contrast, stabilizing selection likely reduces genetic variance (Arnold & Wade, 1984;Endler, 1986) and favors the adaptation to optimum phenotypes (de Vladar & Barton, 2014).
The difference in neutral genetic diversity between populations in stable and disturbed environments was opposite to that seen in functional phenotypic diversity and in outlier AFLP loci. That neutral genetic diversity was lower in populations that occupied disturbed compared with stable habitats is in agreement with the interpretation that populations in disturbed habitats were the results of recent colonization events and affected by genetic bottlenecks that lowered primarily neutral genetic diversity. There is experimental evidence that pygmy grasshopper populations can be established by founder groups consisting of no more than six individuals , thus leaving ample opportunity for eroding effects of founder events on neutral diversity. These eroding neutral processes were countered by other, deterministic and selective, processes that influenced and were influenced by functional (but not neutral) diversity, as discussed above.
Across populations, there was a negative correlation between color morph diversity and neutral AFLP diversity. However, color morph diversity was not associated with diversity in outlier AFLP loci. These findings are in agreement with predictions from theory and add to the body of evidence (Holderegger et al., 2006;Leinonen, O'Hara, Cano, & Merilä, 2008;Reed & Frankham, 2001;Whitlock, 2014;Willi et al., 2006) that estimates of neutral genetic diversity cannot be used as substitutes for adaptive or functional genetic diversity to infer evolutionary potential and the ability of populations and species to cope with environmental change.
One interpretation of the negative association between functional and neutral diversity would be that selection has been more efficient in removing unfit color morphs in larger populations, under the assumption that neutral diversity can be considered a proxy of effective population size. However, field studies that have compared quantitative genetic variation in populations of different size have not revealed any clear pattern (Willi et al., 2006). In our study, the correlation between functional color morph diversity and census population size was positive (not negative). Moreover, both functional diversity and census population size were larger (not smaller) in disturbed than in stable environments. These patterns cannot be accounted for by the potential response to selection or evolvability being greater in larger populations.
Our present analyses of genetic structure advance our understanding of the relative roles of selection, ecological events, and neutral processes as drivers of evolutionary divergence between populations of pygmy grasshoppers (Tinnert & Forsman, 2017;Tinnert, Hellgren et al., 2016) and other ecologically similar organisms. That neutral genetic differentiation between populations was associated with functional genetic differentiation (as revealed by Mantel test based on pairwise F ST values for neutral and outlier AFLP loci) might indicate that divergence in functional traits has not only been influenced by selection, but also by gene flow and stochastic processes associated with founder events and drift in small populations (Dudaniec et al., 2018;Leinonen, McCairns, O'Hara, & Merila, 2013;Tibblin et al., 2015). However, pairwise population genetic differentiation was more pronounced on average for outlier loci than for neutral loci, and sampling location accounted for a greater proportion of the total variance in outlier (32%) than in neutral (8%) AFLP loci, thus implicating also natural selection and local adaptation as drivers of population divergence (Dudaniec et al., 2018;Johansson et al., 2016;Landguth & Balkenhol, 2012;Quintela et al., 2014).
That there was no genetic differentiation between populations in the two contrasting environments is an important finding in that it shows that stable and disturbed habitats were not populated by reproductively isolated evolutionary lineages or ecotypes. In our study, the disturbed sites were located on average at greater distances from the coast and at higher elevation than stable sites.
However, two lines of evidence indicate that this geographic clustering did not confound the results and conclusions regarding the effects of disturbance on genetic diversity. First, there was no overall genetic structure or separation between populations in stable and disturbed environments. It is therefore unlikely that environmental conditions or ecological and evolutionary processes related to distance from the coast, but independent of the degree of disturbance, contributed in any important way to the genetic and phenotypic differences between stable and disturbed sites. Second, not only were populations that were closer to the coast genetically similar to inland populations, our estimates of neutral and functional genetic and phenotypic diversity within populations were independent of elevation. A plethora of factors can impact phenotypic, functional, and neutral genetic diversity of natural populations. The correlational approach of our study means that we cannot with certainty identify all the contributing and nonimportant mechanisms. However, our findings are in agreement with the notion that the relative importance of ecological and evolutionary processes as drivers of neutral and functional genetic diversity changes along the stability-disturbance gradient.

| Extensions and caveats
Depending on the questions and hypothesis under investigation, contrasting environments can consist of habitats that represent different environmental conditions, disturbance regimes (as in the present study), populations inhabiting large (source) as opposed to small (sink) habitat patches, populations in ephemeral or temporal as opposed to permanent habitat patches, or populations at the core of the distribution range as opposed to populations in marginal areas that represent invasion fronts of expanding species (Johansson et al., 2016;Noguerales et al., 2016;Quintela et al., 2014;Shine, Brown, & Phillips, 2011;Sunde, Tamario, Tibblin, Larsson, & Forsman, 2018). To broaden utility, environments could be characterized not as binary states but along continuous gradients, such as age of habitat or time since disturbance event (Forsman, Karlsson et al., 2011).
Wing dimorphism is widespread among insects (Harrison, 1980;Schwander & Leimar, 2011). Our study provided more evidence that the incidence of long-winged morph could thus be used as a proxy for recent establishment and immigration in many Orthoptera, Heteroptera, and Coleoptera species. In other animal groups, genetic expression of locomotion and dispersal-enhancing quantitative phenotypic traits can be used as surrogates to infer recent establishment, recolonization, or immigration events (Arnold & Bennett, 1988;Berthouly-Salazar, van Rensburg, Le Roux, van Vuuren, & Hui, 2012;Forsman, Merilä, & Ebenhard, 2011;Shine et al., 2011).
A benefit of using color morph diversity as a proxy for functionally important genetic and phenotypic variation in the present study is that color pattern of pygmy grasshoppers is genetically and developmentally associated with morphology, physiology, life history, and behaviors (Ahnesjö & Forsman, 2003Forsman, 1999Forsman, , 2018Forsman et al., 2002Forsman et al., , 2007Forsman et al., , 2012 and references therein).
Differences in color morph diversity between populations therefore likely reflect imprints of recent establishment and effects of selection acting directly on color pattern and on traits that are associated with color pattern, possibly resulting in a stronger signature.
Developmental plasticity and phenotypic flexibility can severely impact variation in quantitative traits (Bradshaw, 1965;Forsman, 2015;West-Eberhard, 2003) and potentially affect both the magnitude and direction of the difference in comparisons between neutral molecular and functional genetic diversity. This is not an issue in the present study because pygmy grasshopper color patterns and dispersal phenotypes are not affected to any important degree by developmental plasticity (Berggren et al., 2012;Forsman, 2018;Forsman, Karlsson et al., 2011;Karlsson et al., 2009;Nabours, 1929 and references therein). Available evidence also indicates that color morph and wing morph are independent in pygmy grasshoppers (Forsman, 2018).
Interpretation of results from comparisons based on estimates of quantitative genetic variation versus neutral variation is potentially complicated by the conversion of nonadditive (epistatic and dominance) to additive genetic variation and by the physical inhibition of recombination associated with chromosomal rearrangements in small populations (Neiman & Linksvayer, 2005;Walsh & Lynch, 2012;Willi et al., 2006).
The increase in additive genetic variance resulting from these phenomena offers an alternative explanation to why functional phenotypic diversity may be higher than neutral diversity in recently established and bottlenecked populations, even in cases when functional diversity has not contributed to enhanced establishment or stabilized population dynamics. We cannot completely discard the possibility that such conversion of genetic variance has influenced the outcome of the comparison of color morph diversity in the present study. However, our finding that the difference between stable and disturbed environments in neutral AFLP diversity was opposite the diversity in outlier AFLP loci cannot be explained by conversion of nonadditive to additive genetic variation.
The color polymorphism in pygmy grasshoppers likely has a polygenic inheritance (Fisher, 1939;Nabours, 1929). A transcriptomic analysis on Tetrix japonica-a closely related pygmy grasshopper species (see Figure S1 in Forsman, 2018)-has identified several putative pigment-related genes and differential expression patterns of them between adults and larval phases in this species (Qiu, Liu, Lu, & Huang, 2017). Moreover, putative genes involved in juvenile hormone metabolism and signaling pathways that regulate the body color and flight activity in insects are available for the same species (Qiu et al., 2017). A genomic approach taking advantage of advances in sequences technologies such as restriction site-associated DNA sequencing (e.g., RADseq) and of the available genome-wide sequencing data and protein expression profiles from related species can generate a more comprehensive understanding of the genes/ loci under selection. More specifically, quantitative trait analyses can reveal the interplay between the loci that regulate polygenic traits (such as color patterns and wing morphs) to gain insights into the rapid evolution of relevant phenotypic traits by comparing DNA profiles from contrasting environments (Dudaniec et al., 2018;Johansson et al., 2016;Noguerales et al., 2016;Quintela et al., 2014;Zhi-Xiang et al., 2018).
When interpreting results from comparisons of neutral and functional genetic diversity in populations between contrasting environments, it is necessary to take into consideration any systematic differences in effective population size. It is not only the eroding effect of genetic drift that depends on effective population size.
Evolvability, or the efficiency by which selection can remove deleterious alleles and drive advantageous alleles to fixation, increases with increasing effective population size (Willi et al., 2006). If populations in disturbed, marginal, or temporal habitats are generally smaller than populations in stable, central, or continuous habitats, then this may contribute to a pattern of greater functional than neutral genetic diversity in the former type of environments. It is however unlikely that such a confounding effect of population size biased the results of the grasshopper population comparisons in the present study. If anything, study populations in disturbed environments were larger (not smaller) compared with populations in stable environments and yet they harbored higher functional compared to neutral diversity.

| CON CLUS IONS
Our results for T. undulata grasshoppers illustrate that genetic diversity within natural populations is determined by a combination of ecological and evolutionary processes, the importance of which are different for neutral and functional genetic diversity and vary between environmental disturbance regimes (Figure 4).
A higher proportion of the long-winged phenotype in disturbed habitats suggested that these populations were founded by flightcapable individuals shortly after the disturbance event. Both color morph diversity and outlier genetic dispersion were higher in disturbed than in stable habitats, a likely consequence of higher survival and establishment success of more color morph diverse founder groups , combined with stabilizing selection reducing functional genetic variance in stable environments .
Neutral genetic diversity was lower in populations that occupied disturbed habitats, as a result of recent colonization events and bottlenecks that reduce neutral genetic diversity. We also found that functional phenotypic color morph diversity was negatively correlated with neutral genetic diversity across populations, underscoring that diversity estimates based on neutral markers should not be used as substitutes for adaptive or functional genetic diversity to infer evolutionary potential and the ability of populations and species to cope with environmental change. Our findings highlight the utility of combining morphological data with outlier and neutral loci in analyzing the consequences of environmental change for population genetic structure and diversity, thus informing about key processes in ecology, evolution, and protection of biodiversity.