The P. australis samples were collected within an approximately 10 × 5-km area of grassland at the Pasture Ecology Research Station of Northeast Normal University, Changling, Jilin Province, in the Songnen Prairie of China (123°45′E, 44°38′N). The region is characteristic of a semi-arid and semi-wet continental monsoon climate. The mean annual precipitation is 313–581 mm, 70% of which falls between June and September. Annual evapotranspiration is 3- to 4-fold the amount of precipitation. Mean daily temperatures greater than 10°C can persist for 120–140 days per year and can reach 21.5–23.6°C (Li et al., 2014). Phragmites australis often grows as a companion species in the alkalized meadow but forms monodominant communities in local low-lying areas or saline-alkaline patches (Yang and Li, 2003).

The two sampled groups of P. australis occur widely because the two corresponding habitats lie in a mosaic of various contiguous habitats (Figure 1). Six sampling sites from each habitat were selected (Figure S1, Table S1). The two dryland habitats can be described as follows: (1) In a seasonal waterlogged low-lying meadow (H1), the monodominant reed community occurs with coverage greater than 85%. The soil is alkalized meadow soil with a pH of approximately 8–8.5. Water typically accumulates seasonally during the rainy season. The reed grows well in a reasonably moist soil, and its companion species are mainly Polygonum sibiricum, Scirpus planiculmis, and Echinochloa caudate. (2) In an alkaline patch without seasonally accumulated rainwater (H2), a monodominant reed community has also appeared. The soil is alkalized with a pH greater than 10 and the nutrients in the surface layer of soil have been lost. The soil is hardened and is poorly permeable. The reeds grow sparsely and in tufts with short stems. The community coverage is less than 20%. However, after the rainy season, annual salt-tolerant plants or halophytes flourish, such as Chloris virgata, Kochia scoparia, Suaeda heteroptera, Suaeda corniculata, Suaeda glauca, leading to a community landscape dominated by these species at the end of the growing season (Yang and Li, 2003).

On July 15, July 25, August 5, and August 15, 2012, which corresponded to 60, 70, 80, and 90 days of plant growth, respectively, samples were collected. After August 15, the reeds began to produce inflorescences and senesce. Because relatively more leaves are present at the corresponding positions during the mid-growth phase, this phase is suitable for determining the influence of leaf position. We counted positions from the shoot base to the apex irrespective of leaf death. Chlorophyll a (Chla), chlorophyll b (Chlb), total chlorophyll (Chl(a+b)), and total carotenoid (Car) concentrations were determined spectrophotometrically from 100-mg aliquots of acetone extracts of freshly collected leaf laminae by using 80% acetone according to Lichtenthaler (1987). The pigment concentrations were expressed in mg/g of leaf fresh weight, with three replicates per site. Fully expanded fresh leaves from six other individuals at each site were collected randomly, dried using silica gel, and stored at −20°C for genetic analysis. All of the sampled individuals were separated by a minimum of 30 m to decrease the likelihood of sampling clone-mates. As other authors have shown, in the case of species in which clonal propagation is the main form of reproduction, maximizing the number of sampled sites rather than the number of shoots at each site is more informative of genetic diversity (Kirk et al., 2011). Minimizing the impact of reed collection on local soil erosion was also taken into consideration. Three soil samples (15-cm length of side, 15-cm depth) from the individual site corresponding to the common garden in the field were measured for soil moisture, pH and EC, which were the main soil properties according to which habitats were classified (Nada et al., 2015; Table 1). Soil moisture was determined by oven drying a known weight of each soil sample at 105°C until a constant weight was reached. Soil pH and EC were determined using a PHS-3C pH meter and a DDS-307 EC meter (Shanghai precision and scientific instrument, Shanghai, China) in a 1:5 soil-water solution (Li et al., 2014).

Genomic DNA was extracted from 72 individuals (36 from each habitat) using a modified cetyltrimethylammonium bromide (CTAB) method (Kidwell and Osborn, 1992). A standard AFLP analysis (Vos et al., 1995) was performed with minor modifications (Wang et al., 2005). Selective amplification primers fluorescently labeled at the 5′-end (Applied Biosystems Inc., Foster City, CA, USA) were used, along with an ABI-automated 3730XL DNA capillary sequencer. Restriction-ligation was conducted using an EcoRI/MseI endonuclease mixture (New England Biolabs, Massachusetts, USA) and double-stranded adaptors. The four EcoRI+3/MseI+3 primer pairs that provided the most reliable, consistently scorable bands were chosen after being screened for selective amplification (Table S2). A ROX-500-labeled internal size standard (Applied Biosystems) was added to each sample, and the software GeneMapper v.4.1 (Applied Biosystems) was utilized to collect and score the raw fluorescence data. We scored all well-resolved and reliable bands with a binary code, where zero represented an absent band, and one represented a band that was present. We excluded singleton observations from the dataset (i.e., markers with only one non-consensus band). One plant that continued to produce noisy sequencer electropherograms after the analyses were run several times was excluded. All scoring was performed “blindly” by the same person, who lacked any information about the samples. Non-overlapping peaks between 50 and 500 bp were included.

The rhizomes of P. australis were collected in May 2012 from three sites (genotypes) distributed in the western, central and eastern areas of the research station for each habitat; these three sites were among the six sites established for the field samples (Figures S1, S2). Three clones per site for each habitat were taken to the campus of Northeast Normal University, Changchun, Jilin Province (125°18′E, 43°51′N). This region experiences a continental monsoon, and the climate is similar to that of the region where the research station is located. All of the rhizomes were physically separated from one other and cultivated in an experimental plot in the sand for 6 weeks to produce at least 50 asexual individuals from the nodes of rhizomes per site. Then, healthy rhizomes were cut into 10-cm lengths, each containing one node, and only those with a uniform height aboveground were randomly transplanted into plastic pots (38-cm diameter, 30-cm height, three rhizome blocks per pot) containing peat soil. Each clone includes three pots.

In June 2015, when the plants had been grown for 3 years, showing comparable vegetative development, possible environmental maternal effect was removed fully. The ramets were 8 weeks old, and 54 plants (3 genotypes × 2 habitats × 3 replications × 3 treatments) were randomly assigned to three levels of alkali salt treatment. Two alkali (NaHCO₃ and Na₂CO₃) salts were selected based on the salt components of the salt-alkaline soils over much of northeast China and were mixed in a 9:1 molar ratio (NaHCO₃:Na₂CO₃) (Li et al., 2010). The alkaline stress groups were labeled AS1–AS3. The stress treatments were applied every 6 days, with the salt concentration gradually increasing for 7 weeks from 30 to 700 and 750 mmol L⁻¹ for treatment groups AS2 and AS3, respectively, thereby causing the soil pH to approach 10 and simulating H2. The AS1 group was continually treated with 30 mmol L⁻¹ salt until the pH value approached 8.0–8.5 to simulate H1. All of the pots were sheltered from the rain and watered when necessary. The inter-pot distances were sufficient to prevent competition among the plants for light. The photosynthetic pigment characteristics were measured from the fourth-youngest leaves of the three tallest shoots per pot. The properties of soil samples from 3 pots per habitat and per treatment were measured (Table 1).

Loci under selection present in the tails of the F_ST distribution generated across loci as a function of heterozygosity between habitats were identified based on F-statistics generated using Arlequin 3.5 (Excoffier and Lischer, 2010). A null distribution of F_ST close to the empirical distribution, which was originally developed by Beaumont and Nichols (1996), was acquired with 20,000 coalescent simulations and 100 simulated demes per group. A more rigorous identification of candidate loci was implemented in BayeScan 2.1, which directly calculates q-values for each locus (Foll and Gaggiotti, 2008). The sample size was set to 5000 and the thinning interval to 10 following 20 pilot runs. Then, the loci were split into two subsets, “selected” and “neutral.” Bayenv2 was run for 10⁶ iterations by four independent runs and testing for association with a single environmental variable using both a linear model and Spearman rank correlation while accounting for among-population structure through neutral parameterization by control data (Günther and Coop, 2013).

The photosynthetic data were analyzed using SPSS Statistics 17.0 (IBM, New York). Student's t-test was used to test for differences between habitats at the plant level and the leaf level, and soil property differences. The data were tested for normal distributions and variance homogeneity with the Levene test, and data were log-transformed where necessary. Two one-way ANOVA tests were conducted with sampling date and treatment as the fixed factor and the pigment traits and soil properties as dependent variables, respectively. When the overall variation in an ANOVA was significant, post-hoc comparisons among means were performed using Tukey's HSD test at α = 0.05.

Correlations between pigment concentrations in each of AS1 and AS2 or AS3 were tested. The one-way ANOVA tests were conducted with habitat as an independent variable in each treatment. The patterns of variation in pigment traits between habitats and among treatments and genotypes were analyzed by a General Linear Model (GLM) implemented using the UNIANOVA procedure of PASW Statistics 18 (formerly SPSS) (IBM, New York), with treatment and habitat as fixed factors and genotype, nested within habitats, as a random factor. The change in variability between habitats was described by calculating SSpop/SStotal, where SSpop is the sum of squares between the two habitats and SStotal is the total variability identified from the one-way ANOVA.

To test whether the among-habitat divergence in pigment characters was greater than expected if they evolved neutrally, the Q_ST-F_ST approach was implemented (Defaveri and Merilä, 2014). We first decomposed the variation of each pigment trait, between habitats (σB2) and within habitats (σW2) in each treatment, by using the VARCOMP procedure (REML method) of PASW to calculate the observed Q_ST following eqn (1) (σB2 in place of σBE2). The σW2 is the variance component among clones within sites. Then, the expected among-habitat variance component (σBE2) and within-habitat variance component (σW2) values were used to calculate the expected Q_ST of a neutral trait as follows: (1)QST=σBE2/(σBE2+2σW2) The σBE2 of a neutrally evolving trait was obtained using F_ST, which was calculated from the above-mentioned neutral loci and theσW2 according to Equation (2). (2)σBE2=2FSTσW21-FST We simulated the distribution of the test statistic, Q_ST–F_ST, 10,000 times, which is the null hypothesis of evolution in a neutrally evolving trait. We then tested if the observed Q_ST–F_ST differed by the neutral expectations. The analyses were carried out in R using code from Lind et al. (2011).

The Chla, Chlb, Chl(a+b), and Car contents and the Chl(a+b)/Car ratio of the plants in H2, alkaline patch without seasonally accumulated rainwater were 25.4, 25.6, 25.5, 13.7, and 12.6% lower, respectively, than those of the plants in H1, the seasonally waterlogged low-lying meadow (P < 0.001) (Figure S3). However, no significant differences in the Chla/Chlb ratio were observed (P = 0.444). Those photosynthetic pigment traits decreased markedly with increasing time in each habitat, whereas only the Chla/Chlb ratio showed the opposite trend (Figure S4). The H2 plants consistently presented lower values than did the H1 plants. The values of Chla, Chlb, Chl(a+b), Car and the Chl(a+b)/Car ratio in H1 declined by 36.3, 45.9, 38.5, 21.5, and 22.7%, respectively, compared with the values at the first sampling date. Similarly, the values of those traits in H2 decreased by 45.0, 52.1, 46.6, 27.0, and 27.6%, respectively. Thus, the larger changes over time in H2 relative to H1 were indicative of habitat differences. The decrease in Car content and the Chl(a+b)/Car ratio over time varied across these two habitats, as shown in Figure S4.

Comparisons of pigment traits at the same leaf position between the habitats revealed significant differences between pairs of middle leaves but no significant differences in the basal or apical leaf pairs (Table S3, Figure S3). Most of the pairs showing significant differences shifted gradually to relatively higher positions with time. The pigment traits generally exhibited the same changes as shown in Figures S3, S4 regardless of leaf position except for Car, for which stable values were obtained from the middle leaves of position 9 in H1 and position 7 in H2 (data not shown).

A total of 1132 AFLP fragments from 71 sampled individuals were unambiguously scored from primer combinations, with a relatively low scoring error rate of 0.99%. All of the six outlier loci identified (0.53% of the total) exhibited F_ST values that were markedly greater than expected (Figure 2, mean ± SE = 0.8000 ± 0.0140) and are therefore likely to be subject to directional selection or linked to loci under selection. The remaining 1126 loci exhibited lower levels of differentiation (mean ± SE = 0.0364 ± 0.0029). The outliers identified by Arlequin 3.5 and BayeScan2.1 were consistent. Bayenv was used to serve as a cross-check for consistency with those of different model approaches. The six outliers showed high Bayes factors in the top 1% of the BF-values (BF >3) and the high absolute values of ρ (|ρ|>0.4). Therefore, the signal can be considered to be robust. The six putative outliers also showed high absolute values of ρ (|ρ|>0.4) when the environmental variables without normalization were used (Table S4).

The ANOVA tests revealed significant effects of habitat, alkali salt treatments, habitat-by-treatment and treatment-by-genotype interaction on most of the traits (Table 2) except Chl(a+b)/Car, for which no significant effect of habitat was found. No significant variation in any trait was found among genotypes. One-way ANOVA also revealed significant between-habitat variation in each treatment for all traits except Chl(a+b)/Car in the group AS1 (P = 0.872). Comparison of reaction norms showed that the local genotypes had higher values of Chla, Chlb, Chl(a+b), and Car in their simulated local habitat than did genotypes from the other habitat (Figure 3). These between-habitat pigment concentration differences were sustained not only during the common garden growth of group AS1 but also in group AS2 and AS3 under increasing pH, with larger percentages of variability between habitats in the latter two groups (Table 3). Soil moisture, pH and EC differed significantly between the two habitats (all P ≤ 0.001) (Table 1). The fitness trade-off (negative correlation) expected under a pattern of local adaptation was found in the groups AS1 and AS3 for Chla, Chlb and Car but not Chl(a+b)/Car (Figure 4). Chla, Chlb and Car within each treatment showed evidence of adaptive divergence between habitats, with Q_ST–F_ST values significantly larger than the Q_ST–F_ST values expected for a neutrally evolving trait (Figure 5). The data for Chl(a+b) are not shown as the patterns of changes in Chla and Chlb were the same.

The flexibility of photosynthetic pigment traits is an important adaptive mechanism through which plant increase mean fitness to meet the challenges inherent in a variable environment (Aubin-Horth and Renn, 2009). In this study, photosynthetic pigment traits in two habitats with different severities of soil alkaline conditions and moisture were investigated. The Chla, Chlb, Chl(a+b), and Car contents and the Chl(a+b)/Car ratio were shown to be generally lower in Habitat 2 and decreased over time. These results were expected because Chl loss is often a common effect of stress and occurs during senescence (Dinç et al., 2012). The age-related variation of the Chla/Chlb ratio and TChl content appears to be an adaptive response to stress conditions (Ebrahimiyan et al., 2013; Zhang et al., 2013), as recently inferred for juvenile mastic trees and P. australis (Bragato et al., 2009; Juvany et al., 2013). Water stress, as one of the main causes of alkali stress, commonly causes a decrease in leaf Chl concentration in grasses, which helps reduce light absorption. Because only Chla/Chlb ratio did not differ significantly between habitats, the relative high value of Chla/Chlb ratio in Habitat 2 is an index of drought adaptation due to a smaller light-harvesting complex, making photosystem II less susceptible during stress (Caudle et al., 2014). The consistent patterns of Chl(a+b), Car and the Chl(a+b)/Car ratio were likely due to the lower rate of loss of Car than of Chl(a+b), which reflects an adaptive strategy because carotenoid plays an important role in photoprotection and antioxidation (Ebrahimiyan et al., 2013; Zhang et al., 2013). The information about increased antioxidant capacity under different stresses has been reported largely (Naudts et al., 2014; Xu et al., 2014, 2015). Decrease in antenna pigments such as Chl b suggested that plants experienced a process of antenna size readjustment. All of these mechanisms may act in coordination to allow successful survival for P. australis.

At the low and high relative leaf positions, the differences between habitats were generally small, whereas the differences at the middle leaf position were much larger (Table S3). This pattern may be due to the different physiological status of these leaves (Lippert et al., 2001). The upper leaves can neither accumulate mineral elements nor inhibit ion toxicity efficiently and show pigment contents to be enhanced, whereas the lower, older leaves undergo senescence, with gradual declines in their ion-transportability and photosynthetic capacities. In contrast, middle leaves are mature and may be functionally important for adaptation to stressful environments, considering their susceptibility to various stresses including N supply, herbicide, light intensity, tropospheric ozone and so on (Zhang et al., 2010; Yoon et al., 2011; Wang et al., 2014).

Speciation is often initiated by the adaptive divergence of populations occupying ecologically distinct environments (Endler, 1977; Hendry et al., 2002). Local adaptation to diverse environments is an important evolutionary step before divergence. Although local adaptation at a small scale has been considered rare due to high gene flow, it is a common result of the driving force of divergent selection at a small scale, as noted by Kawecki and Ebert (2004). Examples of this phenomenon have recently accumulated in species with high gene flow (Defaveri and Merilä, 2014), which are relevant for specific ecological processes involving elevated pH values in the soil (Bone and Farres, 2001). Because plants are generally sessile, and if environmental gradients are strong, gene flow might be non-random, all of these could thus generate adaptive genetic differentiation in quantitative characters at small spatial scales (Shapiro et al., 2006; Ross et al., 2009). It is therefore not surprising that this study found evidence for local adaptation.

Levels of genetic differentiation can be highly variable across a genome (Nosil et al., 2009). In the case of model organisms, it is easy to track the “outlier locus” under selection based on whole-genome information, whereas for non-model organisms, such as P. australis, it is difficult to identify candidate genes due to restricted genomic resources (Wang et al., 2012). AFLPs have frequently been used to scan genomes in non-model organisms. They can provide a large number of random loci scattered throughout the genome and detect candidate genomic regions without identifying loci of interest a priori (Herrera and Bazaga, 2008; Feng et al., 2016). The most commonplace approaches applied for AFLP genome scans have been Dfdist, Bayescan analyses until now (Yang et al., 2016). The six outlier loci were identified consistently by these two methods and examined by Bayenv method which provide a useful statistic for ranking markers by their allele frequency correlation with an environmental variable. Although the sampled habitat patches were not adjacent, the various patches among which the common reeds were distributed were parapatric. Thus, it appears reasonable to find that the proportion of the genome exhibiting outlier behavior (0.53%) at a small scale was relatively low relative to the range (0.4–9.0%) reported for outlier loci in previous studies (Nosil et al., 2009).

Although the field experiment provided confirmatory evidence of the differences in Chla, Chlb, Chl(a+b), Car, and the Chl(a+b)/Car ratio, the results of the common garden experiment did not provide robust evidence concerning Chl(a+b)/Car ratio. In the common garden experiment, we used population replication to rule out other forces (e.g., drift or migration) as the cause of differentiation as different samples from similar habitats should show similar responses to given environmental conditions if the differentiation is driven by those conditions (Defaveri and Merilä, 2014). Because the climatic conditions were very similar at the small scale considered here, soil properties should be the key factors responsible for selection regimes. It seems likely that soil conditions, especially soil moisture, pH and EC, determined the habitat-specific variation and drove local adaptation in photosynthetic pigment contents.

In conclusion, photosynthetic pigment concentrations and the Chl(a+b)/Car ratio in the field were shown to significantly decrease in the habitat with a higher soil pH value and less moisture. Only the Chla/Chlb ratio showed no alkali salt dependent trend in P. australis. The common garden experiment indicated that local adaptation can be maintained over small spatial scales. Additional support for this inference was provided by the population-genomic approach that detected outlier loci, and Q_ST–F_ST comparisons. Soil moisture, pH and EC differentiated among the two habitats, which may have acted as selective agents to produce such an adaptive divergence.

The authors declare that the research was conducted in the absence of any commercial or financial relationships that could be construed as a potential conflict of interest.