# Molecular and quantitative genetic differentiation in
[i]Sitobion avenae[/i] populations from both sides of the
Qinling Mountains.

Xianliang Huang, Deguang Liu, Da Wang, Xiaoqin Shi, Jean-Christophe
Simon

► To cite this version:

Xianliang Huang, Deguang Liu, Da Wang, Xiaoqin Shi, Jean-Christophe Simon. Molecular and quantitative genetic differentiation in [i]Sitobion avenae[/i] populations from both sides of the Qinling Mountains.. PLoS ONE, 2015, 10 (3), pp.e0122343. 10.1371/journal.pone.0122343. hal-01274511v2

HAL Id: hal-01274511

https://hal.science/hal-01274511v2

Submitted on 1 Mar 2016

**HAL** is a multi-disciplinary open access
archive for the deposit and dissemination of sci-
entific research documents, whether they are pub-
lished or not. The documents may come from
teaching and research institutions in France or
abroad, or from public or private research centers.

L'archive ouverte pluridisciplinaire HAL, est
destinée au dépôt et à la diffusion de documents
scientifiques de niveau recherche, publiés ou non,
émanant des établissements d'enseignement et de
recherche français ou étrangers, des laboratoires
publics ou privés.

<figure><img src="image_1.png" /></figure>

Distributed under a Creative Commons CC BY 4.0 - Attribution - International License
RESEARCH ARTICLE

# Molecular and Quantitative Genetic Differentiation in *Sitobion avenae* Populations from Both Sides of the Qinling Mountains

Xianliang Huang<sup>1,2,3</sup>, Deguang Liu<sup>1,2,3*</sup>, Da Wang<sup>1,2,3</sup>, Xiaoqin Shi<sup>4</sup>, Jean-Christophe Simon<sup>5</sup>

1 State Key Laboratory of Crop Stress Biology for Arid Areas (Northwest A&F University), Yangling, Shaanxi Province, China, 2 Key Laboratory of Integrated Pest Management on Crops in Northwestern Loess Plateau, Ministry of Agriculture, Yangling, Shaanxi Province, China, 3 College of Plant Protection, Northwest A&F University, Yangling, Shaanxi Province, China, 4 Department of Foreign Languages, Northwest A&F University, Yangling, Shaanxi Province, China, 5 Institut National de la Recherche Agronomique (INRA), unité mixte de recherche (UMR) 1349, Institut de Génétique, Environnement et Protection des Plantes (IGEPP), Domaine de la Motte, Le Rheu, France

* dgliu@nwsuaf.edu.cn

<figure><img src="image_2.png" /></figure>

CrossMark
click for updates

<figure><img src="image_3.png" /></figure>

OPEN ACCESS

**Citation:** Huang X, Liu D, Wang D, Shi X, Simon J-C (2015) Molecular and Quantitative Genetic Differentiation in *Sitobion avenae* Populations from Both Sides of the Qinling Mountains. PLoS ONE 10 (3): e0122343. doi:10.1371/journal.pone.0122343

**Academic Editor:** Ulrich Melcher, Oklahoma State University, UNITED STATES

**Received:** December 5, 2014

**Accepted:** February 12, 2015

**Published:** March 30, 2015

**Copyright:** © 2015 Huang et al. This is an open access article distributed under the terms of the Creative Commons Attribution License, which permits unrestricted use, distribution, and reproduction in any medium, provided the original author and source are credited.

**Data Availability Statement:** All relevant data are within the paper.

**Funding:** This research was supported by the Specialized Research Fund for the Doctoral Program of Higher Education from Ministry of Education of China (No. 20110204120001) (http://www.moe.edu.cn/). DL received the funding. The funder had no role in study design, data collection and analysis, decision to publish, or preparation of the manuscript.

**Competing Interests:** The authors have declared that no competing interests exist.

## Abstract

Quantitative trait differences are often assumed to be correlated with molecular variation, but the relationship is not certain, and empirical evidence is still scarce. To address this issue, we sampled six populations of the cereal aphid *Sitobion avenae* from areas north and south of the Qinling Mountains, and characterized their molecular variation at seven microsatellite loci and quantitative variation at nine life-history traits. Our results demonstrated that southern populations had slightly longer developmental times of nymphs but much higher lifetime fecundity, compared to northern populations. Of the nine tested quantitative characters, eight differed significantly among populations within regions, as well as between northern and southern regions. Genetic differentiation in neutral markers was likely to have been caused by founder events and drift. Increased subdivision for quantitative characters was found in northern populations, but reduced in southern populations. This phenomenon was not found for molecular characters, suggesting the decoupling between molecular and quantitative variation. The pattern of relationships between F<sub>ST</sub> and Q<sub>ST</sub> indicated divergent selection and suggested that local adaptation play a role in the differentiation of life-history traits in tested *S. avenae* populations, particularly in those traits closely related to reproduction. The main role of natural selection over genetic drift was also supported by strong structural differences in G-matrices among *S. avenae* populations. However, cluster analyses did not result in two groups corresponding to northern and southern regions. Genetic differentiation between northern and southern populations in neutral markers was low, indicating considerable gene flow between them. The relationship between molecular and quantitative variation, as well as its implications for differentiation and evolution of *S. avenae* populations, was discussed.
## Introduction

Barriers to gene flow presented by mountain ranges can have significant consequences for ecol-
ogy and evolution of various organisms, and even lead to the isolation of populations. The iso-
lation of populations may in turn cause allopatric speciation in different environments where
reproductive isolation between populations could occur gradually and incidentally as a result
of mutation and drift, and indirect effects of local selections (causing local adaptations) [1-2].
The Qinling Mountains (with highest peak of 3, 767 m) in Shaanxi Province of China, which
extend for nearly 2,500 kilometers in the east-west direction, provide a unique scenario to
study such effects. The mountain ranges are located in the transition zone between two macro-
climatic regimes (i.e., subtropical and warm-temperate zones) [3]. The regions south of the
mountains have subtropical characteristics with wet summers and warm winters, while the re-
gions north of them belong to warm-temperate and temperate zones with relatively dry sum-
mers and cold winters [4].

Significant geographic variation in life-history traits (e.g., developmental time and fecundity) can occur among populations of a species separated by environmental barriers like the Qinling Mountains [4]. Determining the mechanisms that cause such variation is one of the basic objectives of ecological and evolutionary studies. It has long been recognized that geographic variation in life-history traits of an organism can include both genetic and environmental components [5]. The genetically-determined components of life-history trait variation are a significant factor which drives evolution by supplying raw material for evolutionary changes. Genetic variation can be affected not only by gene flow, but by natural selection, genetic drift, colonization history, mutation, and their interactions [6-7]. Among these causes, natural selection is generally considered as the major force generating changes in the level of genetic variation within and among populations (especially in heterogeneous environments), and thus driving adaptation and evolutionary changes [7-9]. A population's ability of responding to changing environments depends on additive genetic variation for quantitative traits that are ecologically relevant [10]. Due to technical and logistic difficulties in measuring such variation, many studies rely solely on molecular markers to evaluate adaptive genetic variation, however, the degree to which molecular variation reflects quantitative genetic variation between populations is still controversial [10].

The cosmopolitan aphid (*Sitobion avenae* (Fab.)), a serious pest on cereal crops [11], offers
a good model organism to study how differentiation can evolve among populations of small or-
ganisms with assumed high dispersal ability for several reasons. *Sitobion avenae* is widely dis-
tributed in areas around the huge Qinling Mountains, which can be a significant barrier for its
dispersal, since its occurrence on wild hosts is rare above 2, 000 m in the mountains (D. L., per-
sonal observation). This aphid reproduces by cyclical or obligate parthenogenesis, and can be
maintained clonally for an indefinite period under controlled laboratory conditions. So quanti-
tative traits of interest can be evaluated by using clonal individuals, and the partitioning of phe-
notypic variation into its genetic and environmental components can be achieved. Finally,
microsatellite markers have been developed for this species, allowing for evaluation of its mo-
lecular variation and population structure. In our previous study [4], considerable divergence
in life-history traits among *S. avenae* populations from north and south of the Qinling Moun-
tains was identified, but it is still not clear about the relative importance of different evolution-
ary forces in shaping and maintaining genetic variation in those populations. We hypothesize
that the separation of Qinling Mountains can cause increasing quantitative trait divergence
and local adaptation of *S. avenae* populations in northern and southern Shaanxi Province, thus
creating different population structuring patterns in those areas. To do so, we compared popu-
lations of *S. avenae* collected from areas north of the Qinling Mountains (hereafter referred to
as northern populations) to those from areas south of the mountains (hereafter referred to as southern populations) using both *Q*<sub>ST</sub> (an index of differentiation in quantitative traits [12]) and *F*<sub>ST</sub> (an index of differentiation in neutral genetic markers [12]). The aims of the present study were: (a) to characterize patterns of genetic differentiation among *S. avenae* populations at molecular traits (evaluated using seven microsatellite markers) and at quantitative life-history traits (measured under common laboratory conditions), (b) to assess the relationship between molecular and quantitative trait variation, and (c) to evaluate the relative importance of selection and drift in shaping and maintaining genetic variation between northern and southern populations.

## Materials and Methods

### Aphid sampling and colony establishment

From April to May in 2013, individuals of *S. avenae* were sampled from six locations including three each in areas north and south of the Qinling Mountains. The aphid clones came from wheat fields in the north (collected at Fuping, 34°46'46"N, 109°01'56"E; Tongchuan, 35°07'05"N, 109° 06'58"E; and Luochuan, 35°46'26"N, 109°24'29"E) and south (collected at Jinshui, 33°16'10"N, 107°47'13"E; Longting, 33°12'41"N, 107°38'31"E; and Mianxian, 33°11'36"N, 106°56'55"E) of the Qinling Mountains (no particular permissions were needed for sample collecting activities at sites mentioned above, and no endangered or protected species were involved in the collecting activities). Aphid samples were collected from the field, and insect colonies were maintained in the lab as described in detail previously in [9]. Briefly, over 15 wingless adults of *S. avenae* were collected in each of the five or more fields randomly selected at each location. Winter wheat seeds (*Triticum aestivum* cv. Aikang 58) were planted, and aphid colonies were kept in rearing rooms. For eliminating environmentally induced effects, the *S. avenae* populations were reared under common laboratory conditions for three generations before bioassays. Even without genetic changes among populations, maternal effects produce adaptive plastic responses that can be mistaken as local adaptation, but after three generations of cultivation under common laboratory conditions such effects can become negligible [13].

### Life-history data collection

Plants used in life-history tests were grown, and aphid life-history data were collected as described in detail previously in [9]. Briefly, wheat seedlings of one- to two-leaf stage (one per plant) received aphids transferred from rearing plants. Each plant with an aphid individual was enclosed with a transparent plastic cylinder (6 cm in diameter, 15 cm in height), and maintained in environmental growth chambers at 20±1°C, a relative humidity of 65±2%, and a photoperiod of 16:8 (L:D). Thirty clones were randomly selected from aphid colonies of each location for use in the tests, and three to four replicates were conducted for each clone. Test aphid individuals were observed twice daily from birth until the onset of reproduction, and molting occurrences and mortality counts were recorded. After reproduction started, numbers of offspring and dead aphids were recorded daily, and their offspring were then removed until the death of all test aphids.

### Quantitative trait analysis

Developmental time of first to fourth nymphal instars, the total developmental time of nymphs, adult lifespan, reproductive time, post-reproductive time, and lifetime fecundity were calculated as described in [9]. The abovementioned quantitative traits were analyzed with nested analyses of variance (nested ANOVAs) (sources of variance: ‘region’, ‘location’ nested
within ‘region’, and ‘clone’ nested within ‘location’) in SAS [14]. The variance explained by ‘clone’ (δ²<sub>clone</sub>) was determined. When overall variation in ANOVA was significant, separation of treatment means was carried out using Tukey tests at α < 0.05. Data were log-transformed when needed to meet the assumptions of normality and homoscedasticity (i.e., homogeneity of variance) in ANOVA.

Principal component analysis (PCA) was performed with all quantitative traits measured above after the raw data were log-transformed [14]. The purpose of PCA is to quantify the significance of the variables that explained the differences in the abovementioned quantitative traits, and identify similarity in group data structures [9]. The first two PCA components (PC1 and PC2) were plotted using the PROC GPLOT procedure. The factor weightings of each replicate for PC1 from the PCA were calculated, and they were used as a composite life-history factor (i.e., PC1) in subsequent analyses.

Our bioassays use clonal genotypes, and this design allows us to estimate the total variance of a particular quantitative trait (*V*<sub>P</sub>), which can be partitioned into among-clone genetic components *V*<sub>G</sub> (i.e., the broad-sense genetic variance) and within-clone components *V*<sub>E</sub> (i.e., environmental variance or residual variance) [9]. Broad-sense heritabilities were calculated as described previously in [9]. Genetic variance and covariance estimates for life-history traits were obtained with the restricted maximum likelihood (REML) method implemented in the software VCE 6.0.2 [15]. The genetic correlation between traits *x* and *y* was calculated from the genetic covariance estimate (cov(*x*, *y*)) and their additive variances as r = cov(*x*, *y*)/[(*v*<sub>*x*</sub>) × (*v*<sub>*y*</sub>)]<sup>0.5</sup>. The Flury hierarchical method was utilized to compare the resulting G matrices in the software CPCrand [16]. Based on maximum likelihood, this method can identify structural differences among G matrices by analyzing their eigenvectors and eigenvalues as described in [17]). Briefly, the method can test, in order, the models of unrelated structure, partial common principal components, common principal components, proportionality, and equality (see also in [17]). The statistical significance of genetic correlations and broad-sense heritabilities was assessed with likelihood-ratio tests (LRTs) (for more details, see also in [9, 17]).

$Q_{ST}$ values were evaluated using genetic variances and defined as $Q_{ST} = V_{GB} / (V_{GB} + 2V_{GW})$, where $V_{GB}$ is the genetic component of the variance between population means and $V_{GW}$ is the average genetic variance within populations [18]. $Q_{ST}$ was estimated between each pair of populations, as well as between northern and southern regions. Differentiation among populations within a region was evaluated similarly as $Q_{SR}$. Pairwise population $Q_{ST}$ values (a global index of differentiation) were also computed for the composite life-history trait, which was the first component obtained from the PCA performed on abovementioned quantitative traits. Following Chapuis *et al.* [7], $Q_{ST}$ was considered different from corresponding $F_{ST}$ if their respective 95% confidence intervals did not overlap. Another method of evaluating the variation in genetic variances for comparison among characters and populations is the coefficient of genetic variance ($CV_G = V_G^{1/2}/X$, where X is the mean phenotype) [19].

## Microsatellite data collection and analysis

Over 40 clones per population were randomly selected and genotyped at seven microsatellite loci (Sm10, Sm12, S17b, Sm17, S16b, Sa4Σ and S5L) following the procedures described in Simon et al. [20].

Using identified microsatellite loci, population differentiation was assessed by calculating pairwise $F_{ST}$ (θ) [21] between populations and between both regions. Deviations from Hardy-Weinberg equilibrium and unbiased estimates of $F_{IS}$ (f) [21] were estimated. Microsatellite diversity was assessed by calculating unbiased estimates of gene diversity ($H_E$) [22] for each sample. Region-specific $F_{ST}$ (hereafter referred to as $F_{SR}$), $F_{IS}$, $H_E$, $H_o$ (observed heterozygosity),
<table><caption>Table 1. Mean phenotypic values (SE) of all life-history traits for *Sitobion avenae* populations from both sides of Qinling Mountains (i.e., northern and southern regions).</caption><thead><tr><th rowspan="2">Traits</th><th rowspan="2">Northern regions</th><th rowspan="2">Southern regions</th><th colspan="2">Difference between regions</th><th colspan="2">Difference among populations nested in regions</th></tr><tr><th>F</th><th>P</th><th>F</th><th>P</th></tr></thead><tbody><tr><td>DT1</td><td>1.66 (0.02)</td><td>1.75 (0.04)</td><td>F<sub>1,215</sub> = 5.26</td><td>0.023</td><td>F<sub>4,215</sub> = 7.48</td><td>&lt; 0.001</td></tr><tr><td>DT2</td><td>1.51 (0.02)</td><td>1.65 (0.03)</td><td>F<sub>1,215</sub> = 30.09</td><td>&lt; 0.001</td><td>F<sub>4,215</sub> = 5.26</td><td>&lt; 0.001</td></tr><tr><td>DT3</td><td>1.54 (0.02)</td><td>1.69 (0.03)</td><td>F<sub>1,213</sub> = 25.60</td><td>&lt; 0.001</td><td>F<sub>4,213</sub> = 8.30</td><td>&lt; 0.001</td></tr><tr><td>DT4</td><td>2.02 (0.02)</td><td>2.00 (0.03)</td><td>F<sub>1,212</sub> = 0.20</td><td>0.658</td><td>F<sub>4,212</sub> = 1.66</td><td>0.161</td></tr><tr><td>DT5</td><td>6.73 (0.03)</td><td>7.09 (0.06)</td><td>F<sub>1,212</sub> = 53.77</td><td>&lt; 0.001</td><td>F<sub>4,212</sub> = 9.96</td><td>&lt; 0.001</td></tr><tr><td>Lifetime fecundity</td><td>56.54 (0.87)</td><td>72.11 (1.57)</td><td>F<sub>1,199</sub> = 216.05</td><td>&lt; 0.001</td><td>F<sub>4,199</sub> = 131.23</td><td>&lt; 0.001</td></tr><tr><td>Post-reproductive time</td><td>8.11 (0.30)</td><td>15.20 (0.74)</td><td>F<sub>1,199</sub> = 259.10</td><td>&lt; 0.001</td><td>F<sub>4,199</sub> = 66.77</td><td>&lt; 0.001</td></tr><tr><td>Adult lifespan</td><td>28.67 (0.47)</td><td>37.75 (0.58)</td><td>F<sub>1,199</sub> = 353.29</td><td>&lt; 0.001</td><td>F<sub>4,215</sub> = 268.05</td><td>&lt; 0.001</td></tr><tr><td>Reproductive time</td><td>20.56 (0.30)</td><td>22.55 (0.67)</td><td>F<sub>1,199</sub> = 20.86</td><td>&lt; 0.001</td><td>F<sub>4,199</sub> = 111.25</td><td>&lt; 0.001</td></tr></tbody></table>

Note: DT1-DT4, the developmental time of 1st to 4th instar nymphs; DT5, the total developmental time of nymphs; all traits but lifetime fecundity measured in days; lifetime fecundity, the number of offspring produced in a female's life.

doi:10.1371/journal.pone.0122343.t001

and $R_S$ (allele richness) were also estimated and their differences between regions were tested by performing 1000 permutations of alleles among individuals. The abovementioned F-statistics and related tests were conducted in FSTAT version 2.9.3 [23].

We also analyzed the microsatellite data with an individual based principal component analysis (PCA) using the program PCA-GEN version 1.2.1 [24]. This analysis uses allele frequencies to define new variables (components) that can characterize the neutral genetic variation among populations. PCA is a preferable method for visualizing microsatellite data in situations of potential gene flow because traditional bifurcating trees constructed from genetic distances can be difficult to interpret [25].

## Results

### Differences in quantitative traits

Northern populations showed significantly reduced developmental time for the first, second and third instar nymphs, but not for the fourth instar nymphs (Table 1). The total developmental time of nymphs for northern populations was also significantly lower than that for southern populations. Southern populations had much higher lifetime fecundity than northern populations. Post-reproductive time and adult lifespan of southern populations were much higher than those of northern populations, whereas the reproductive time of southern populations was only slightly higher than that of northern populations. Similarly, all the abovementioned quantitative traits showed significant variation among populations within a region except the developmental time of fourth instar nymphs. Of the nine quantitative characters, eight differed significantly among populations within regions, as well as between the two regions.

### Quantitative genetic differentiation

Genetic difference between northern and southern populations was shown in Table 2. There were no significant differences between northern and southern populations in genetic variance ($\delta^2_{clone}$). However, the coefficient of genetic variance ($CV_G$) and mean broad sense heritability
<table><caption>Table 2. Comparison of within-region molecular and quantitative genetic differentiation between northern and southern regions.</caption><thead><tr><th>Parameters</th><th>Northern region</th><th>Southern region</th><th>P value</th></tr></thead><tbody><tr><td>δ<sup>2</sup><sub>clone</sub></td><td>18.9913</td><td>7.4387</td><td>0.056</td></tr><tr><td>CV<sub>G</sub></td><td>0.1417</td><td>0.0741</td><td><strong>0.048</strong></td></tr><tr><td>H<sup>2</sup></td><td>0.4736</td><td>0.3579</td><td><strong>0.034</strong></td></tr><tr><td>Ho</td><td>0.744</td><td>0.630</td><td>0.602</td></tr><tr><td>R<sub>S</sub></td><td>2.3</td><td>3.7</td><td>1.000</td></tr><tr><td>F<sub>IS</sub></td><td>-0.679</td><td>-0.550</td><td>1.000</td></tr></tbody></table>

Note: The genetic variance ($\delta^2_{clone}$), the coefficient of genetic variance ($CV_G$) and the broad sense heritability ($H^2$) were measured for each region and over all quantitative traits; allelic richness ($R_S$), observed heterozygosity ($Ho$), and inbreeding coefficient ($F_{IS}$) were estimated for each region and over all loci, and significance tests were performed between northern and southern regions by randomization procedures using FSTAT software.

doi:10.1371/journal.pone.0122343.t002

$(H^2)$ for the populations in the northern region were both significantly higher than those in the southern region.

Significant genetic correlations were found between life-history traits for *S. avenae* populations from both regions (Table 3). Significantly negative correlation (i.e., trade-off) between DT1 and DT2 was found for southern populations, but not for northern populations. DT1, DT2 and DT3 were all shown to be significantly and negatively correlated with FEC, POS, SPA and RET for northern populations, whereas all the genetic correlations were non-significant for southern populations except that DT2 was significantly correlated with FEC and SPA. DT5 was significantly and negatively correlated with FEC, POS, SPA and RET for northern populations, whereas it was only significantly correlated with FEC among the characters for southern populations. The correlations of POS-FEC and POS-RET were significantly positive for northern populations, but they were significantly negative for southern populations. The correlations of FEC-SPA, FEC-RET and SPA-RET were all significantly positive for populations from both regions.

G-matrix comparisons by Flury's method and jump-up approach (that is, at each step in the hierarchy, the hypothesis is tested against the hypothesis of unrelated structure) showed

<table><caption>Table 3. Genetic correlations between life-history traits for northern (above the diagonal) and southern (below the diagonal) populations of *Sitobion avenae*.</caption><thead><tr><th></th><th>DT1</th><th>DT2</th><th>DT3</th><th>DT4</th><th>DT5</th><th>FEC</th><th>POS</th><th>SPA</th><th>RET</th></tr></thead><tbody><tr><th>DT1</th><td></td><td>-0.143</td><td>0.076</td><td>0.011</td><td>0.489***</td><td>-0.388*</td><td>-0.163*</td><td>-0.279*</td><td>-0.282*</td></tr><tr><th>DT2</th><td>-0.168*</td><td></td><td>-0.214*</td><td>-0.016</td><td>0.377*</td><td>-0.300*</td><td>-0.310*</td><td>-0.327*</td><td>-0.219*</td></tr><tr><th>DT3</th><td>-0.043</td><td>-0.226*</td><td></td><td>0.042</td><td>0.506***</td><td>-0.460***</td><td>-0.202*</td><td>-0.364*</td><td>-0.350*</td></tr><tr><th>DT4</th><td>-0.124</td><td>0.299*</td><td>-0.334*</td><td></td><td>0.468***</td><td>-0.150</td><td>-0.020</td><td>0.030</td><td>0.062</td></tr><tr><th>DT5</th><td>0.478**</td><td>0.438**</td><td>0.317*</td><td>0.437**</td><td></td><td>-0.471***</td><td>-0.281*</td><td>-0.339*</td><td>-0.299*</td></tr><tr><th>FEC</th><td>-0.018</td><td>-0.371**</td><td>-0.038</td><td>-0.285*</td><td>-0.263*</td><td></td><td>0.619***</td><td>0.852***</td><td>0.857***</td></tr><tr><th>POS</th><td>-0.090</td><td>-0.135</td><td>0.045</td><td>-0.133</td><td>-0.117</td><td>-0.404***</td><td></td><td>0.848***</td><td>0.527***</td></tr><tr><th>SPA</th><td>0.021</td><td>-0.191*</td><td>0.027</td><td>-0.082</td><td>-0.072</td><td>0.253*</td><td>0.494***</td><td></td><td>0.875***</td></tr><tr><th>RET</th><td>0.068</td><td>-0.063</td><td>-0.086</td><td>0.017</td><td>0.035</td><td>0.675***</td><td>-0.670***</td><td>0.305*</td><td></td></tr></tbody></table>

Note: DT1-DT4, the developmental time of 1st to 4th instar nymphs; DT5, the total developmental time of nymphs; FEC, lifetime fecundity; POS, post-reproductive time; SPA, adult lifespan; RET, reproductive time; statistical significance of genetic correlations evaluated using likelihood-ratio tests; *, P < 0.05; **, P < 0.01; ***, P < 0.001.

doi:10.1371/journal.pone.0122343.t003
<table><caption>Table 4. Comparisons of G-matrices for quantitative traits of *Sitobion avenae* populations from north and south of the Qinling Mountains.</caption><thead><tr><th rowspan="2">Population source</th><th rowspan="2">G matrices</th><th colspan="3">Flury hierarchy</th></tr><tr><th>LRT</th><th>P-value</th><th>Verdict</th></tr></thead><tbody><tr><td rowspan="3">Northern region</td><td>Fuping vs. Luochuan</td><td>199.6</td><td>&lt; 0.001</td><td>Unrelated</td></tr><tr><td>Fuping vs. Tongchuan</td><td>110.9</td><td>&lt; 0.001</td><td>Unrelated</td></tr><tr><td>Luochuan vs. Tongchuan</td><td>35.6</td><td>&lt; 0.01</td><td>CPC(1)</td></tr><tr><td rowspan="3">Southern region</td><td>Jinshui vs. Longting</td><td>90.3</td><td>&lt; 0.001</td><td>Unrelated</td></tr><tr><td>Jinshui vs. Mianxian</td><td>74.0</td><td>&lt; 0.001</td><td>Unrelated</td></tr><tr><td>Longting vs. Mianxian</td><td>75.0</td><td>&lt; 0.001</td><td>CPC(1)</td></tr><tr><td>Between the two regions</td><td></td><td>88.2</td><td>&lt; 0.001</td><td>Unrelated</td></tr></tbody></table>

Note: The verdict of the Flury hierarchy is the model shown to be the best in explaining the difference between paired matrices; the *P*-values are for the test of equality of two matrices; CPC(1), one of the nine possible components shared in common; unrelated, no relations between the matrices.

doi:10.1371/journal.pone.0122343.t004

significant differences between paired matrices of life-history traits (Table 4). The CPC(1) model best explained the differences between matrices for populations Luochuan and Tongchuan from the northern region (LRT = 35.6, *P* < 0.01), and between those for Longting and Mianxian from the southern region (LRT = 75.0, *P* < 0.001), in other words, matrices shared only one of the nine possible principal components. However, unrelated structures were found for all other comparisons (LRT = 74.0–199.6, *P* < 0.001). So, G-matrices comparison showed significant structural differences among populations within regions, as well as between regions. The strong structural differences identified can be interpreted as a result of divergent selection among populations.

## Microsatellite differentiation

Microsatellite variation among populations (Table 5) was evident, with overall gene diversity for southern populations ranging from 0.285 in population Mianxian to 0.766 in population Jinshui. Overall gene diversity for single northern populations varied from 0.393 in population Luochuan to 0.538 in population Fuping. So, variation in gene diversity among populations within a particular region tended to be higher in the south than in the north. Levels of gene diversity tended to be lower in the north (0.466 ± 0.049) than in the south (0.530 ± 0.086), although the difference was not significant (*P* > 0.05). The differences in allelic richness (*R*<sub>s</sub>), observed heterozygosity (*Ho*), and inbreeding coefficient (*F*<sub>IS</sub>) between the two regions were not significant (all *P* > 0.05, Table 2).

<table><caption>Table 5. Gene diversity (SE) for each locus in each population.</caption><thead><tr><th>Population</th><th>Sm10</th><th>Sm12</th><th>S17b</th><th>Sm17</th><th>S16b</th><th>Sa4Σ</th><th>S5L</th><th>All loci</th></tr></thead><tbody><tr><td>Fuping</td><td>0.500</td><td>0.500</td><td>0.607</td><td>0.619</td><td>0.567</td><td>0.476</td><td>0.500</td><td>0.538 (0.022)</td></tr><tr><td>Tongchuan</td><td>0.500</td><td>0.650</td><td>0.350</td><td>0.650</td><td>0.500</td><td>0.250</td><td>0.350</td><td>0.464 (0.058)</td></tr><tr><td>Luochuan</td><td>0.500</td><td>0.500</td><td>0.350</td><td>0.600</td><td>0.600</td><td>0.200</td><td>0.000</td><td>0.393 (0.085)</td></tr><tr><td>Jinshui</td><td>0.833</td><td>0.808</td><td>0.802</td><td>0.771</td><td>0.825</td><td>0.731</td><td>0.594</td><td>0.766 (0.032)</td></tr><tr><td>Longting</td><td>0.500</td><td>0.500</td><td>0.000</td><td>0.500</td><td>0.833</td><td>0.000</td><td>0.000</td><td>0.333 (0.126)</td></tr><tr><td>Mianxian</td><td>0.500</td><td>0.500</td><td>0.000</td><td>0.500</td><td>0.500</td><td>0.000</td><td>0.000</td><td>0.285 (0.101)</td></tr><tr><td>All populations</td><td>0.576</td><td>0.622</td><td>0.377</td><td>0.631</td><td>0.654</td><td>0.305</td><td>0.260</td><td>0.489 (0.064)</td></tr><tr><td>Northern populations</td><td>0.501</td><td>0.550</td><td>0.436</td><td>0.623</td><td>0.561</td><td>0.309</td><td>0.283</td><td>0.466 (0.049)</td></tr><tr><td>Southern populations</td><td>0.671</td><td>0.760</td><td>0.318</td><td>0.645</td><td>0.754</td><td>0.313</td><td>0.248</td><td>0.530 (0.086)</td></tr></tbody></table>

doi:10.1371/journal.pone.0122343.t005
## Q<sub>ST</sub> versus F<sub>ST</sub>

Pairwise Q<sub>ST</sub> and F<sub>ST</sub> values for all tested *S. avenae* populations were shown in Table 6. Pairwise F<sub>ST</sub> values ranged from 0 for population Longting vs. population Mianxian to 0.1228 for population Jinshui vs. population Mianxian (Table 6). Molecular variation was found to be very low between populations Luochuan, Tongchuan, Longting and Mianxian, with all of the pairwise F<sub>ST</sub> values smaller than 0.021. Substantial molecular variation was found between population Fuping and population Longting (F<sub>ST</sub> = 0.1068) or Mianxian (F<sub>ST</sub> = 0.1083). Substantial molecular variation was also found between population Jinshui and population Longting (F<sub>ST</sub> = 0.1197) or Mianxian (F<sub>ST</sub> = 0.1228).

Pairwise Q<sub>ST</sub> values varied from 0.0094 for population Luochuan vs. population Tongchuan to 0.6215 for population Fuping vs. population Luochuan. In terms of quantitative traits, population Fuping was significantly different from population Luochuan (Q<sub>ST</sub> = 0.6215) or Tongchuan (Q<sub>ST</sub> = 0.5554). Population Luochuan from north of the mountains differed substantially from all of the populations from south with pairwise Q<sub>ST</sub> values ranging from 0.3181 to 0.5936, and similar substantial variation was found between population Tongchuan and any population from south (pairwise Q<sub>ST</sub> varied from 0.2624 to 0.524). However, differentiation in quantitative traits was low between population Fuping of the northern region and any population of the southern region (pairwise Q<sub>ST</sub> varied from 0.0102 to 0.0808). Quantitative trait differentiation between populations from south of the mountains was also found to be low with Q<sub>ST</sub> values varying from 0.0377 to 0.0612. Of all the pairwise comparisons, eight of them showed significantly higher Q<sub>ST</sub> than corresponding F<sub>ST</sub>.

The Q<sub>ST</sub> values of the developmental times of first to fourth instar nymphs and the total developmental time of nymphs were not significantly higher than corresponding F<sub>ST</sub> values (Fig. 1). However, lifetime fecundity, post-reproductive time, adult lifespan and reproductive time all showed significant differentiation among populations and significantly exceeded the neutral expectation set by F<sub>ST</sub> (i.e., Q<sub>ST</sub> > F<sub>ST</sub>).

The magnitude of subdivision in quantitative traits (i.e., Q<sub>SR</sub>) was significantly higher for northern populations than that for southern population, but there were no significant differences in molecular subdivision (i.e., F<sub>SR</sub>) between northern and southern regions (P > 0.05, Fig. 2). Quantitative trait differentiation among populations (Q<sub>SR</sub>) in the north was significantly higher than molecular subdivision among populations (F<sub>SR</sub>) in the same region. F<sub>SR</sub> appeared to be higher than Q<sub>SR</sub> in the south, but no significant differences were found between the two parameters.

## Principal component analyses (PCA)

The first two components from PCA of life-history traits explained 91.4% of the total data variability (78.4% and 13.0% respectively for PC1 and PC2; Fig. 3, upper panel). Post-reproductive

<table><caption>Table 6. Pairwise Q<sub>ST</sub> calculated from a composite life-history trait (below the diagonal) and F<sub>ST</sub> (above the diagonal) estimates based on seven microsatellite loci for *Sitobion avenae* populations (the composite life-history trait obtained from PC1 of PCA analysis of all tested life-history traits; *, significant differences between Q<sub>ST</sub> and corresponding F<sub>ST</sub>).</caption><thead><tr><th></th><th>Fuping</th><th>Luochuan</th><th>Tongchuan</th><th>Jinshui</th><th>Longting</th><th>Mianxian</th></tr></thead><tbody><tr><td>Fuping</td><td></td><td>0.0879</td><td>0.0648</td><td>0.0419</td><td>0.1068</td><td>0.1083</td></tr><tr><td>Luochuan</td><td>0.6215*</td><td></td><td>0.0102</td><td>0.0978</td><td>0.0035</td><td>0.0036</td></tr><tr><td>Tongchuan</td><td>0.5554*</td><td>0.0094</td><td></td><td>0.0700</td><td>0.0206</td><td>0.0209</td></tr><tr><td>Jinshui</td><td>0.0808</td><td>0.5936*</td><td>0.5240*</td><td></td><td>0.1197</td><td>0.1228</td></tr><tr><td>Longting</td><td>0.0278</td><td>0.3181*</td><td>0.2624*</td><td>0.0612</td><td></td><td>0.0000</td></tr><tr><td>Mianxian</td><td>0.0102</td><td>0.4519*</td><td>0.3845*</td><td>0.0377</td><td>0.0591</td><td></td></tr></tbody></table>

doi:10.1371/journal.pone.0122343.t006
<figure><img src="image_4.png" /><figcaption>Fig 1. Differentiation ($Q_{ST}$ ±SE) between northern and southern populations for different quantitative traits [the dotted line represents $F_{ST}$ (±SE) between the two populations based on allelic variation in microsatellite loci; *, $Q_{ST}$ significantly higher than $F_{ST}$; DT1-DT4, the developmental time of 1st to 4th instar nymphs; DT5, the total developmental time of nymphs; $Q_{ST}$, the index of quantitative variation based on life-history traits; $F_{ST}$, the index of molecular variation derived from microsatellite markers].</figcaption></figure>

doi:10.1371/journal.pone.0122343.g001

time (loading: 0.868) and adult lifespan (loading: 0.367) contributed the most to the first principal component (PC1). PC2 was associated mainly with reproductive time (loading: 0.596), lifetime fecundity (loading: 0.571), adult lifespan (loading: 0.316), and post reproductive time (loading: -0.451), the last one with negative correlation. The three southern populations (i.e., Jinshui, Longting and Mianxian) clustered together in the lower right of the plot, whereas the three northern populations (i.e., Fuping, Luochuan, and Tongchuan) were separated in the plot.

PCA on microsatellite data showed that the first two axes accounted for 95.6% of the total inertia (77.7% and 17.9% respectively for PC1 and PC2; Fig. 3, lower panel). This analysis clustered four populations (i.e., Longting, Mianxian, Luochuan and Tongchuan) together in the middle left of the plot, whereas the Jinshui and Fuping population fell in the upper and lower right of the plot, respectively.

## Discussion

### Relationship between molecular and quantitative variation

Many studies use molecular markers as a surrogate for adaptive genetic variation [10], and very few studies have addressed the relationship between molecular and quantitative subdivision for different populations of a species. *Sitobion avenae* provides a good model to explore
<figure><img src="image_5.png" /><figcaption>**Fig 2.** Comparison between molecular (F<sub>SR</sub>) and quantitative genetic (Q<sub>SR</sub>) divergence within each region (*, significant differences between Q<sub>SR</sub> and F<sub>SR</sub> based on non-overlapping confidence intervals; NS, non-significant differences between Q<sub>SR</sub> and F<sub>SR</sub>; F<sub>SR</sub> of the northern region was not significantly different from that of the southern region, P = 0.458).</figcaption></figure>

doi:10.1371/journal.pone.0122343.g002

this relationship because its clonal individuals can be reared in common laboratory conditions for several generations (to minimize confounding environmental effects) and then phenotyped [4], allowing for better estimates of heritable variation and comparisons between Q<sub>ST</sub> and F<sub>ST</sub>.

In our study, northern and southern populations showed significantly different quantitative traits and G-matrices, so high quantitative variation was evident between populations from the two regions. However, microsatellite data did not show the expected variation between both types of populations. In addition, populations within the northern region had increased subdivision in quantitative traits (shown by higher Q<sub>SR</sub> and CV<sub>G</sub>), but not in molecular characters. The decoupling of molecular and quantitative variation was also evident in PCA analyses where clustering patterns from molecular and quantitative trait data were distinct. The conversion of non-additive variation underlying quantitative traits to additive variation could occur due to isolation and drift following bottlenecks, and quantitative variation can be increased as a result [10, 26]. So, frequent bottlenecks may contribute to the decoupling of molecular and quantitative variation in our case, since *S. avenae* lives on highly seasonal and sometimes unpredictable habitats [27]. In our analyses of microsatellite markers, overall molecular differences between northern and southern populations appeared to be very low (F<sub>ST</sub> = 0.018), indicating populations from the two regions shared the majority of identified genotypes. Therefore, high dispersal ability of *S. avenae* and the resulting gene flow can be an important factor in retarding the molecular differentiation between populations from both sides of the Qinling Mountains. The different pattern for molecular and quantitative traits might also be simply explained by differences in mutation rates, because mutation rates for quantitative
<figure><img src="image_6.png" /><figcaption>Fig 3. Plot of PC1 versus PC2 from principal component analysis of life-history traits (upper panel) and microsatellite loci (lower panel) as variables (northern populations: Fuping, FP; Luochuan, LC; Tongchuan, TC; southern populations: Jinshui, JS; Longting, LT; Mianxian, MX).</figcaption></figure>

doi:10.1371/journal.pone.0122343.g003

characters are typically several orders of magnitude higher than those for molecular characters
[28].

Indeed, the extent to which molecular markers reflect quantitative genetic subdivision be-
tween populations is still controversial [10]. Our results suggested that molecular variation
might not be a reliable indicator for quantitative variation, and could only provide a conservative estimate of adaptive divergence. The decoupling of molecular and quantitative variation in our study is consistent with the findings of McKay and Latta [29], Palo et al. [30] and Chapuis et al. [7] that no significant correlations existed between molecular and quantitative variation, although positive correlations have been reported in both plants and animals [19, 31–32]. So, our results provided data cautioning against the use of molecular metrics of variation as a surrogate for the ability of *S. avenae* populations to respond to environmental changes. Further studies are still needed to clarify the relationship between molecular and quantitative variation in insects, as well as in other organisms.

## Differentiation and evolution of *S. avenae* populations

In our current study, all tested life-history characters but one (i.e., DT4) differed significantly between northern and southern regions, as well as among populations within a region. However, northern populations showed a higher DT4 than southern populations in our previous study [4]. Compared to those from northern populations, individuals from southern populations generally had longer developmental times, but higher fecundity, adult lifespan and reproductive time in this study. In contrast, no significant differences in DT2 and reproductive time were found between northern and southern populations in our previous study [4]. The seemingly contradictory results in some life-history traits may be attributed to different population sampling time of both studies since environmental conditions can vary significantly at the same location in the same month of different years. Despite the abovementioned inconsistency, consistent differences between northern and southern populations were found in DT1, post-reproductive time, adult lifespan, and lifetime fecundity in both studies. Northern and southern populations also differed significantly in genetic correlations between quantitative traits and in G-matrices. So, genetic variation among *S. avenae* populations from both sides of the Qinling Mountains was evident.

Between-population differentiation for neutral genetic markers—seven microsatellite loci— varied with *F*<sub>ST</sub> ranging from 0 to 0.1228 (*F*<sub>ST</sub> > 0.1 implies a very high degree of differentiation among populations, see more in Simon *et al.* [20] and Vialatte *et al.* [33]. The high genetic differentiation in neutral markers found between particular populations can be partially explained as a founder effect, due to fast growth and large population sizes after a bottleneck in population foundation by limited number of individuals. In several cases of pairwise comparisons (e.g., Luochuan vs. Tongchuan), no significant differences were found between *F*<sub>ST</sub> and *Q*<sub>ST</sub>, so it can not be ruled out that genetic drift alone could account for the differentiation of *S. avenae* populations in those cases.

The role genetic drift plays in genetic differentiation of *S. avenae* populations is well documented [33–35], but the impact of natural selection has received little attention. *Sitobion avenae* showed a high potential to respond to natural selection, and rapid changes in ecologically relevant traits (e.g., phenotypic plasticity) have been reported in response to changes in selective pressure [27]. However, empirical data on the actual occurrence of local adaptation are scarce. In our study, differentiation of four quantitative traits closely related to reproduction (i.e., lifetime fecundity, post-reproductive time, adult lifespan and reproductive time) exceeded the neutral expectation set by *F*<sub>ST</sub> (i.e., *Q*<sub>ST</sub> > *F*<sub>ST</sub>), indicating divergent selection among *S. avenae* populations. There were over half pairwise population comparisons in the *Q*<sub>ST</sub>–*F*<sub>ST</sub> analysis for the composite life-history trait that yielded higher *Q*<sub>ST</sub> values than their corresponding *F*<sub>ST</sub> estimates (i.e. more genetic divergence in the trait than expected by neutral genetic drift only), providing additional evidence of divergent selection. Populations of *S. avenae* that are better adapted to seasonal and ephemeral habitats are expected to have higher fitness under
our experimental conditions (e.g. higher lifetime fecundity). In this study and [4], southern populations showed much higher fitness than northern populations, suggesting the occurrence of local adaptation. So, our results provide substantial evidence for patterns that suggest local adaptation in *S. avenae*. Another way to test for local adaptation is to perform reciprocal transfer of individuals in wild populations, and it is logistically demanding and extremely difficult to carry out such experiments [7,36]. However, such experiments may further substantiate the occurrence of local adaptation in *S. avenae*.

The occurrence of local adaptation in *S. avenae* is not unexpected in our study, because the separation of Qinling Mountains might cause the isolation of its populations. In addition, the pattern of natural selection can vary greatly between areas north and south to the Qinling Mountains because of different environmental conditions (esp., weather). The mountains cause harsh winter weather conditions in the north and mild ones in the south. The differences in winter weather conditions can have significant consequences for *S. avenae*'s life-history, for example, the reproductive mode of *S. avenae* clones in Romania (harsh winter) was different from that in France (mild winter) [34]. Thus, local adaptation, as a result of natural selection of local conditions, appears to be a common phenomenon for this insect. In this study and [4], natural selection seemed to favor individuals with smaller DT1, shorter post-reproductive time and shorter adult lifespan in the north of the mountains; on the contrary, those with larger DT1, longer post-reproductive time and longer adult lifespan were selected for in the south. We observed significant heritabilities, which show that populations have the potential to respond to local selection. As a result, it is likely that some individuals have adapted to local conditions.

Despite the highly likely occurrence of local adaptation in northern and southern populations, *S. avenae* individuals seemed to be able to disperse across the Qinling Mountains because cluster analyses using both quantitative trait (in this study and [4]) and microsatellite data didn't result in two groups corresponding to northern and southern regions. The amount of dispersal of *S. avenae* individuals across the huge mountains was thought to be very low, because the Tianshan Mountains in Xinjiang of China led to clear geographic separation between northern and southern populations of a similar wheat aphid, *D. noxia* [37]. On the contrary, a substantial amount of *S. avenae* dispersal across the Qinling mountains could occur, because the overall low genetic differentiation between northern and southern populations (FST = 0.018) indicated a significant amount of gene flow between them. Therefore, the Qinling Mountains did not result in the isolation of *S. avenae* populations in our study area, and allopatric speciation for *S. avenae* appears to be highly unlikely.

*Sitobion avenae* presents a good model to study how differentiation can evolve among populations of small organisms with high dispersal ability, which is crucial to understand speciation and the match between environment and phenotype [12]. Further experiments are needed to determine what specific factors are important in differentiation and evolution of *S. avenae* populations. Temperature is a particularly likely candidate for the cause of differences between populations identified in our study. Because development is closely tied to temperature, natural selection might cause individuals from colder, northern climates to develop more rapidly than southern ones when raised under common laboratory conditions. The structure of G-matrix for *S. avenae*'s quantitative traits also has implications for its population differentiation and evolution. Interestingly, quite a few negative covariances (i.e., trade-offs) were found between quantitative traits measured in our study. These trade-offs (also shown by negative genetic correlations) may play a role in slowing the evolution of habitat specific *S. avenae* genotypes and thus slowing the differentiation of *S. avenae* populations in our study. The G-matrix structure of different quantitative traits for aphid genotypes have complex relationships with factors like genotype specialization, trade-offs, and genotype-by-environment interactions [38]. The
identified strong structural differences in G-matrices for *S. aveane* populations represented a rare phenomenon, because empirical studies have generally supported the stability of G-matrices for populations of the same species [39]. Further studies are needed to explore the implications of G-matrix structure and stability in maintaining and shaping genetic differentiation in *S. avenae*, as well as in its evolutionary potential.

## Acknowledgments

We want to thank P. Dai and Z.-H. Yi (Northwest A&F University, China) for their laboratory and field assistance. We are indebted to Y.-G. Hu and X.-S. Hu (Northwest A&F University, China) who provided plant seeds in this study.

## Author Contributions

Conceived and designed the experiments: DL XH. Performed the experiments: XH JCS DW. Analyzed the data: DL XH XS. Wrote the paper: DL XH JCS.

## References

1. Hoskin CJ, Higgie M, McDonald KR, Moritz C. Reinforcement drives rapid allopatric speciation. Nature. 2005; 437: 1353-1356. PMID: 16251964

2. Via S. Natural selection in action during speciation. PNAS. 2009; 106 Suppl. 1: 9939-9946. doi: 10.1073/pnas.0901397106 PMID: 19528641

3. Bai J, Yan J-P, Su K-H. Differential analysis of abrupt climate change between southern and northern Qinling Mountains in the past 50 years. J Shaanxi Norm Univ (Nat. Sci. Edit.). 2010; 38: 98-105.

4. Huang X-L, Liu D-G, Gao S-X, Chen H. Differential performance of *Sitobion avenae* populations from both sides of the Qinling Mountains under common garden conditions. Environ Entomol. 2013; 42: 1174-1183. doi: 10.1603/EN13132 PMID: 24246739

5. Endler JA. Geographic variation, speciation, and clines. Princeton: Princeton University Press; 1977.

6. Jolivet C, Bernasconi G. Molecular and quantitative genetic differentiation in European populations of *Silene latifolia* (Caryophyllaceae). Ann Bot. 2007; 100: 119-127. PMID: 17565967

7. Chapuis E, Trouve S, Façon B, Degen L, Goudet J. High quantitative and no molecular differentiation of a freshwater snail (*Galba truncatula*) between temporary and permanent water habitats. Mol Ecol. 2007; 16: 3484-3496. PMID: 17688547

8. Byers DL. Evolution in heterogeneous environments and the potential of maintenance of genetic variation in traits of adaptive significance. Genetica. 2005; 123: 107-124. PMID: 15881684

9. Gao S-X, Liu D-G, Chen H, Meng X-X. Fitness traits and underlying genetic variation related to host plant specialization in the aphid *Sitobion avenae*. Insect Sci. 2014; 21: 352-362. doi: 10.1111/1744-7917.12085 PMID: 24243773

10. Edmonds S, Harrison JS. Molecular and quantitative trait variation within and among populations of the intertidal copepod *Tigriopus californicus*. Evolution. 2003; 57: 2277-2285. PMID: 14628915

11. Gao S-X, Liu D-G. Differential performance of *Sitobion avenae* clones from wheat and barley with implications for its management through alternative cultural practices. J Econ Entomol. 2013; 106: 1294-1301. PMID: 23865194

12. Campillo S, García-Roger EM, Carmona MJ, Gómez A, Serra M. Selection on life-history traits and genetic population divergence in rotifers. J Evol Biol. 2009; 22: 2542-2553. doi: 10.1111/j.1420-9101.2009.01871.x PMID: 19878499

13. Pitchers WR, Brooks R, Jennions MD, Tregenza T, Dworkin I, Hunt J. Limited plasticity in the phenotypic variance-covariance matrix for male advertisement calls in the black field cricket, *Teleogryllus commodus*. J Evol Biol. 2013; 26: 1060-1078. doi: 10.1111/jeb.12120 PMID: 23530814

14. SAS. SAS System, version 6.12, SAS Institute Inc., Cary, NC, USA; 1998.

15. Neumaier A, Groeneveld E. Restricted maximum likelihood estimation of covariances in sparse linear models. Genet Selection Evol. 1998; 30: 3-26.

16. Phillips PC, Arnold DE. Hierarchical comparison of genetic variance-covariance matrices. I. Using the Flury hierarchy. Evolution. 1999; 53: 1506-1515.
17. Carter MJ, Simon J-C, Nespolo RF. The effects of reproductive specialization on energy costs and fitness genetic variances in cyclical and obligate parthenogenetic aphids. Ecol Evol. 2012; 2: 1414–1425.
doi: 10.1002/ece3.247 PMID: 22957150

18. Spitze K. Population structure in *Daphnia obtusa*: quantitative genetic and allozymic variation. Genetics. 1993; 135: 367–374. PMID: 8244001

19. Morgan KK, Hicks J, Spitze K, Latta L, Pfrender ME, Weaver CS, et al. Patterns of genetic architecture for life-history traits and molecular markers in a subdivided species. Evolution. 2001; 55: 1753–1761.
PMID: 11681731

20. Simon J-C, Baumann S, Sunnucks P, Hebert PDN, Pierre J-S, Le Gallic J-F, et al. Reproductive mode and population genetic structure of the cereal aphid *Sitobion avenae* studied using phenotypic and microsatellite markers. Molecular Ecology. 1999; 8: 531–545. PMID: 10327655

21. Weir BS, Cockerham CC. Estimating *F*-statistics for the analysis of population structure. Evolution. 1984; 38: 1358–1370.

22. Nei M. Molecular evolutionary genetics. New York: Columbia University Press; 1987.

23. Goudet J. FSTAT: a computer program to calculate *F*-statistics. J Hered. 1995; 86: 485–486.

24. Goudet, J. PCA-GEN. for Windows v. 1.2 [www document], 1999. Available: http://www2.unil.ch/popgen/softwares/pcagen.htm.

25. Lozier JD, Roderick GK, Mills NJ. Genetic evidence from mitochondrial, nuclear, and endosymbiont markers for the evolution of host plant associated species in the aphid genus *Hyalopterus* (Hemiptera: Aphididae). Evolution. 2007; 61: 1353–1367. PMID: 17542845

26. Lopez-Fanjul C, Fernandez A, Toro M. The effect of epistasis on the excess of the additive and nonadditive variances after population bottlenecks. Evolution. 2002; 56: 865–876. PMID: 12093023

27. Dai X-J, Gao S-X, Liu D-G. Genetic basis and selection for life-history trait plasticity on alternative host plants for the cereal aphid *Sitobion avenae*. PloS ONE. 2014; 9(9): e106179. doi: 10.1371/journal.pone.0106179 PMID: 25181493

28. Pfrender ME, Spitze K, Hicks J, Morgan K, Latta L, Lynch M. Lack of concordance between genetic diversity estimates at the molecular and quantitative-trait levels. Conserv Genet. 2000; 1: 263–269.

29. McKay JK, Latta RG. Adaptive population divergence: markers, QTL and traits. Trends Ecol Evol. 2002; 17: 285–291.

30. Palo JU, O'Hara RB, Laugen AT, Laurila A, Primmer CR, Merila J. Latitudinal divergence of common frog (*Rana temporaria*) life history traits by natural selection: evidence from a comparison of molecular and quantitative genetic data. Mol Ecol. 2003; 12: 1963–1978. PMID: 12803645

31. Lynch M, Pfrender M, Spitze K, Lehman N, Hicks J, Allen D, et al. The quantitative and molecular genetic architecture of a subdivided species. Evolution. 1999; 53: 100–110.

32. Steinger T, Haldimann P, Leiss A, Müller- Schärer H. Does natural selection promote population divergence? A comparative analysis of population structure using amplified fragment length polymorphism markers and quantitative traits. Mol Ecol. 2002; 11: 2583–2590. PMID: 12453241

33. Vialatte A, Plantegenest M, Simon J-C, Dedryver C-A. Farm-scale assessment of movement patterns and colonization dynamics of the grain aphid in arable crops and hedgerows. Agric For Entomol. 2007; 9: 337–346.

34. Papura D, Simon J-C, Halkett F, Delmotte F, Le Gallic JF, Dedryver CA. Predominance of sexual reproduction in Romanian populations of the aphid *Sitobion avenae* inferred from phenotypic and genetic structure. Heredity. 2003; 90: 397–404. PMID: 12714986

35. Dedryver C-A, Le Gallic J-F, Haack L, Halketty F, Outreman Y, Simon J-C. Seasonal and annual genotypic variation and the effect of climate on population genetic structure of the cereal aphid *Sitobion avenae* in northern France. Bull Entomol Res. 2008; 98: 159–168. PMID: 18076784

36. Kawecki TJ, Ebert D. Conceptual issues in local adaptation. Ecol Lett. 2004; 7: 1225–1241.

37. Zhang B, Edwards OR, Kang L, Fuller SJ. Russian wheat aphid (*Diuraphis noxia*) in China: native range expansion or recent introduction? Mol Ecol. 2012; 21: 2130–2144. doi: 10.1111/j.1365-294X.2012.05517.x PMID: 22417053

38. Nespolo RF, Halkett F, Figueroa CC, Plantegenest M, Simon J-C. Evolution of trade-offs between sexual and asexual phases and the role of reproductive plasticity in the genetic architecture of aphid life histories. Evolution. 2009; 63: 2402–2412. doi: 10.1111/j.1558-5646.2009.00706.x PMID: 19473379

39. Cano JM, Laurila A, Palo J, Merila J. Population differentiation in *G* matrix structure due to natural selection in *Rana temporaria*. Evolution. 2004; 58: 2013–2020. PMID: 15521458