Temporal variation in selection on male and female traits in wild tree crickets

Abstract Understanding temporal variation in selection in natural populations is necessary to accurately estimate rates of divergence and macroevolutionary processes. Temporal variation in the strength and direction of selection on sex‐specific traits can also explain stasis in male and female phenotype and sexual dimorphism. I investigated changes in strength and form of viability selection (via predation by wasps) in a natural population of male and female tree crickets over 4 years. I found that although the source of viability stayed the same, viability selection affected males and females differently, and the strength, direction and form of selection varied considerably from year to year. In general, males experienced significant linear selection and significant selection differentials more frequently than females, and different male traits experienced significant linear selection each year. This yearly variation resulted in overall weak but significant convex selection on a composite male trait that mostly represented leg size and wing width. Significant selection on female phenotype was uncommon, but when it was detected, it was invariably nonlinear. Significant concave selection on traits representing female body size was observed in some years, as the largest and smallest females were preyed on less (the largest may have been too heavy for flying wasps to carry). Viability selection was significantly different between males and females in 2 of 4 years. Although viability selection via predation has the potential to drive phenotypic change and sexual dimorphism, temporal variation in selection may maintain stasis.


Introduction
Understanding the pace of phenotypic selection in natural populations is important to accurately calculate the rates of divergence and macroevolutionary processes, and to estimate how quickly populations can respond to change, man-made or otherwise (Siepielski et al. 2009). Selection analysis, the measurement of the relationship between phenotypic traits and relative fitness in a population (Lande and Arnold 1983), is an important tool for standardizing estimates of natural selection in order to detect such patterns. Since Lande and Arnold's (1983) landmark paper on quantifying selection using selection analyses, thousands of studies have used these methods, and several meta-analyses have compiled and compared selection gradients Kingsolver et al. 2001;Hereford et al. 2004;Kingsolver and Diamond 2011). However, compiled estimates of selection predict faster microevolution than is generally observed (Kinnison and Hendry 2001). This mismatch could be due to overestimating selection in the long term (over several generations) when using selection gradients measured in short term (over a single generation or single breeding season), as estimates of selection can be greatly affected by the duration of the episode of selection . Longer-term studies of selection on natural populations are valuable as they provide information on how the strength and direction of selection is affected by timescale (Kingsolver et al. 2001;Kingsolver and Pfennig 2007). In longer intervals of time, such as over several generations, temporal variation in the strength and direction of selection and variation in which traits are under selection can reduce the magnitude of net selection across years, which would result in slower evolution than predicted from a single episode of selection (Siepielski et al. 2009(Siepielski et al. , 2011; but also see Morrissey and Hadfield 2012).
Rates of evolution may be more difficult to predict when dealing with sexually dimorphic traits. Sexually dimorphic traits tend to evolve more slowly than monomorphic traits because males and females of the same species are highly genetically correlated (Lande 1980). Consistently different selection on males and females is necessary to produce sexual dimorphism; thus, temporal variation in the difference in the form and direction of selection between males and females can also slow the evolution of sexual dimorphism (Schulte-Hostedde et al. 2002;Reimchen and Nosil 2004). If such temporal variation occurs, then traits that produce optimal fitness at different sizes in males and females may remain at an intermediate size that is optimal for neither sex (i.e., intralocus sexual conflict, Cox and Calsbeek 2009).
However, some long-term studies can reveal consistent selection that can lead to rapid phenotypic change. For instance, rapid evolution has been observed in response to phenotypic selection in populations of stickleback (Bell et al. 2004;Aguirre and Bell 2012), Anolis lizards (Losos et al. 2006) and guppies (Reznick et al. 1997). And consistently different selection on male and female house finches (Carpodacus mexicanus) can lead to the rapid evolution of sexual dimorphism (Badyaev and Martin 2000;Badyaev 2005). These examples are exceptional although, and usually follow a major change, such as introduction of the species to a new environment or the introduction of a new predator to the current environment. Thus, successive generations would have been affected by new and likely strong selection on the same traits. For populations in more stable environments, this may not be the case.
Populations may experience temporal variation in viability selection as a result of several factors. Highly variable environmental factors can affect the relationship between traits and survival (Kalisz 1986;Grant and Grant 2002;Tarwater and Beissinger 2013). Different phenotypes may also confer survival advantages in one generation but not the next due to changes in major sources of mortality. For example, rotund shell shape in Physa snails is advantageous when fish are their dominant predator, but elongate shells are better protection against crayfish predators (DeWitt et al. 2000). Changes in relative predator abundance can also affect the evolution of sexually dimorphic traits, such as spine number in sticklebacks (Reimchen and Nosil 2004). However, it is not well known how selection from the same single source of viability may vary from generation to generation within the same population, or how this variation may affect the evolution of sexually dimorphic traits.
To investigate temporal variation in selection on phenotype and its potential to affect sexual dimorphism, I examined viability selection from a single source of mortality in a natural population of tree crickets, Oecanthus nigricornis Walker. I observed selection over 4 years by comparing the phenotypic distribution of surviving crickets to that of prey of a common cricket specialist wasp, Isodontia mexicana Saussure. Thus, I compared strength and shape of selection among years and between the sexes. Sex differences in viability selection are particularly of interest in this system because the predator preferentially hunts females O'Neill 2003, 2009;Ercit 2014), and this differential predation may affect male and female phenotypes differently, potentially producing sexual dimorphism.

Study organisms
Oecanthus tree crickets (Gryllidae) are common in meadows of east and central North America (Capinera et al. 2005). In the area where this study took place, O. nigricornis are univoltine, and adults emerge in late July and persist approximately until frost arrives in October. They are typically found in open meadows with Solidago spp., Rubus spp., and Daucus spp., as males use these plants to call from, and females use stems as oviposition sites (Fulton 1915). Adult female tree crickets are larger than males, and males have greatly enlarged forewings (tegmina) which, in females, are less differentiated from the hind wings than in males. Male tegmina are used to produce a calling song that attracts receptive conspecific females (Walker 1957;Bell 1980;Toms 1993;Brown et al. 1996). In pair formation, male Oecanthus are mostly stationary, and females are mobile (Fulton 1915;Brown 1999).
Isodontia mexicana are common solitary wasps, found throughout southern Canada and the United States, and can be an important predator of O. nigricornis (Bohart and Menke 1976;Iwata 1976;O'Neill 2003, 2009;Ercit 2014). Female I. mexicana sting and paralyze their prey and carry them back to their nest to provision for their offspring (Iwata 1976). I. mexicana are often inhabitants of artificial trap nests (Krombein 1967), so their provisioning behavior is easy to observe. I. mexicana take significantly more adult female tree crickets than male O'Neill 2003, 2009;Ercit 2014), and this is due partly to female-biased sexual size dimorphism in prey (Ercit 2014), as well, adult female crickets with ovaries heavy with eggs may be easier to catch .

5119
and Ercit and Gwynne (2015), and thus collection methods of these studies are very similar. I sampled wasp prey that was provisioned in artificial trap nests, which were made based on the construction described in (Hallett 2001a,b). Nest blocks were grouped in five stacks of seven, within boxes covered with wooden lids and roofing shingles, and placed on wooden platforms 1 m off the ground.
I collected O. nigricornis prey of I. mexicana from trap nests and compared these to the overall distribution of crickets collected from the surrounding meadow. Samples were taken approximately weekly, starting every year when the first adult tree cricket appeared in a nest until all wasp provisioning activity had stopped. Prey samples were taken from the most recently provisioned cell of I. mexicana nest tunnels, and the entire contents of each recently provisioned cell were taken. Samples of the hunted cricket population were taken for comparison via sweep net on the same day or the previous day as prey samples, from within a 300-m radius of the wasp nest. These crickets will be hereafter referred to as "survivors." All samples were housed in small plastic containers and then fixed in 95% ethanol.

Traits and measurements
Photographs were taken of all cricket samples using an AmScope (Irvine, CA, USA) 5MP microscope digital camera mounted on a Wild Heerbrugg M5A dissecting microscope. I then measured phenotypic traits in the digital photographs using ImageJ (National Institutes of Health, Bethesda, MD, USA) software.
I measured femur length, femur width, tibia length, pronotum length, tegmen width, and head width for all sampled crickets (Table 1). This suite of traits includes both sexually dimorphic (tegmen width and pronotum length) and traits that are monomorphic when accounting for allometry (leg measurements and head width). I included trait tegmen width because tegmina (forewings) are large sound-producing structures in males, and larger wings may attract more predator attention. I included pronotum length (as a proxy of body size, which is significantly larger in females) because previous results suggest that larger crickets are at higher risk of predation by wasps (Ercit 2014). I included leg measurements because leg size may be related to mobility rate (Kelly et al. 2008). Finally, I included head width because it is a sexually selected trait in male tree crickets (Ercit and Gwynne 2015). To reduce multicollinearity and to increase statistical power, I reduced leg measurements into a single principal component axis that explained 91% of variance. All three leg measurements loaded positively on this axis, and it was mostly influenced by tibia length (50%) and femur length (44%). After this reduction, variance inflation scores were all below 6. Body mass was not measured because all prey crickets necessarily weighed less than survivors as a consequence of paralysis and storage in wasp nests: Paralyzed crickets continued to metabolize their energy stores but could not eat. Instead, I used pronotum length as a proxy as it is the strongest measured predictor of body mass (Ercit 2014).

Statistical analysis
All statistics were carried out using R version 3.0.2 (R Development Core Team 2013). First, to see whether sampling years should be analyzed separately, I tested whether there were any significant interaction effects of year on the relationship between traits and fitness. I started with a saturated model that included all traits, quadratic trait terms, sampling date, and all interactions between traits and year. I then simplified the model using backwards stepwise model selection, and averaged the coefficients of terms where models had DAIC < 5.
To investigate the relationship between measured traits and fitness in male and female tree crickets in each year, I used several methods: Firstly, I conducted a cubic spline analysis (Schluter 1988;Schluter and Nychka 1994) to visualize selection. To do this, I conducted a principal component analysis to reduce the measured traits to a single PC axis. I then fit a cubic spline (in a generalized additive model) to the relationship between the PC trait (with the same original trait composition for all years and both sexes) and my estimate of fitness, and plotted this relationship for each sex in each year. Secondly, I calculated standardized selection differentials, which show, in standard deviations, how trait size has changed after selection (Arnold and Wade 1984). This was calculated as the covariance between fitness and standardized trait sizes. These values include phenotypic change as a result of both direct and indirect selection on that trait. I also calculated selection gradients (Lande and Arnold 1983) from multiple regression of standardized traits against my estimate of fitness. This term measures only the force of direct selection on that trait (Arnold and Wade 1984). Fitness was estimated as a score 0 if the cricket was prey of the predatory wasp and 1 if the cricket was a survivor sampled from the remaining population. I did not convert absolute fitness to relative fitness because converting to relative fitness gives the false impression that I sampled survivors and prey in proportion to the frequency in which they were hunted by wasps. For each year and sex, I found linear selection gradients (b) using multiple linear regression (which included only linear terms), and quadratic and correlational selection gradients (c) from separate multiple regression models that included quadratic and cross-product terms. Quadratic selection gradients were obtained by doubling the quadratic coefficients from nonlinear regression (Stinchcombe et al. 2008). As the estimate of fitness was binary, I used logistic regression to generate P-values of regression coefficients (Janzen and Stern 1998). Finally, I conducted canonical analyses to increase the ability to detect nonlinear and correlational selection (Phillips and Arnold 1989;Blows and Brooks 2003). This consisted of multiplying the matrix of standardized trait measurements by the matrix M (the diagonalization of the c-matrix) to obtain composite traits, and conducting a second round of linear and nonlinear regression on these composite traits. Significance of eigenvalues generated by canonical analysis were found using multiple permutation tests (Reynolds et al. 2009), and cross-product terms were added back into the model for the permutation test (Bisgaard and Ankenman 1996). If convex or concave selection was detected, I tested whether that selection was significantly stabilizing or disruptive (respectively) (Mitchell-Olds and Shaw 1987) using MOStest function in the R package "vegan" (Oksanen et al. 2012). To test whether directional, quadratic, and correlational selection was significantly different between males and females each year, I conducted partial F-tests (Chenoweth and Blows 2005). This consisted of conducting an analysis of variance on models of selection on all traits with and without sex as an interaction term. I also used similar partial F-tests to compare the difference in selection between pairs of years of this study. Significance values of partial F-tests were obtained by permutation tests.

Selection over 4 years
Cubic spline analysis showed that the form of viability selection on principal component axis 1 was considerably different from year to year among both males and females (Fig. 1). PC axis 1 captured 69% of total variance and traits influenced this axis in the following proportion: tegmen width À0.35; pronotum length 0.53; head width 0.55; leg size 0.54. Selection on PC1 in males was especially variable, as it changed from negative linear to positive linear to concave to convex from 2009 to 2012. Partial F-tests support that selection on males varied from year to year: Linear selection on male traits was significantly different between 2009 and 2010 as well as between 2010 and 2011 (Table 2a). Cubic spline analysis suggests that the shape of selection on females changes considerably from year to year (Fig. 1), but these differences are not statistically significant (partial F-tests: Table 2a).
Viability selection analyzed across the entire 4-year period was weaker than in individual years for both males and females. The mean magnitude of directional selection gradients within years (calculated by averaging the absolute values of yearly b on individual traits in Table 3) was more than double the magnitude of that across the 4 years (the average of absolute values of b on individual traits over the entire 2009-2012 period) in both males ( x bÀwithinÀyear ¼ 0:190, x bÀacrossÀyear ¼ 0:092, two-tailed t-test, t 18 = 3.16, P < 0.01) and females ( x bÀwithinÀyear ¼ 0:098, x bÀacrossÀyear ¼ 0:025, two-tailed t-test, t 18 = 2.61, P = 0.02). Over the 4-year period, there was no significant linear selection on original traits in either sex (Table 3), but canonical analysis (Table 4)  linear selection (k = À0.385, P = 0.03) on male composite trait m 4 (influenced by leg size and tegmen width) resulting in males with moderately larger tegmina and smaller legs having a survival advantage (Fig. 2). Canonical analysis also found concave selection on a female composite trait m 1 (strongly influenced by head width, k = 0.230, P = 0.03), which resulted in a fitness trough for females with intermediate head widths (Fig. 3). Cubic spline analysis of selection of males and females over the same period shows weak linear selection in opposite directions (Fig. 1). This indicates that there is nonsignificant (P = 0.09), negative selection on male pronotum length, head width, and leg size, and positive selection on wing width (as wing width is loaded in opposition to the other traits on this axis). There is nonsignificant (P = 0.37), weak positive selection on the same trait axis in females. It is important to note that the relative trait loadings of the axis PC1 are different than those of both male axis m 4 and female axis m 1 , which is why they show different relationships between fitness and phenotype. The axis PC1 combines trait values across the 4 years and both sexes in a manner that captures the most variance, whereas axes from canonical rotation were calculated to show the combinations of traits under the strongest nonlinear selection within each period of selection for males and females, respectively. Thus the relative trait loadings of canonical axes under significant selection changed considerably from year to year.

Differences in selection on males and females
Cubic spline analysis indicated that the shape of phenotypic selection each year was quite different between males and females (Fig. 1). Selection on PC1 was strongly negative in males in 2009, but almost flat in females (Fig. 1), and partial F-tests confirm that linear selection was significantly different between males and females in 2009 (F = 2.86, P = 0.01, Table 2b). Partial F-tests also indicate that directional selection was significantly different (F = 2.42, P = 0.03) and quadratic selection was marginally different (F = 2.16, P = 0.08) in 2012. We can see from the cubic splines that selection on PC1 in 2012 was weakly negative in females, but was strongly convex and almost stabilizing in males.
In general, nonlinear selection was much more common than linear in females, and linear selection was more common than nonlinear in males (Tables 3 and 4).

Yearly viability selection on males
There were significant interactions between year and the male traits of head width, leg size, wing width, pronotum size (Supplementary Table S1a), so each year was analyzed separately.
Among males, linear selection was frequently detected, but the traits under significant selection changed every year. In 2009 males, significant directional selection for smaller legs (b = À0.298, P = 0.01) and wider tegmina (b = À0.298, P = 0.01) was detected. Significant selection differentials indicate a reduction in the pronotum length (S = À0.202, P = 0.01) and leg size (S = À0.264, P < 0.01) in 2009. After canonical analysis, I found significant linear selection on composite traits m 3 (for larger tegmen and smaller legs, h = À0.252, P < 0.01) and m 4 (for wider heads and smaller legs, h = 0.249, P = 0.02, Table 4). In 2010, there was no significant selection gradients before or after canonical rotation, but a significant selection differential indicated that tegmen width increased after predation (S = 0.203, P = 0.05, Table 3).
Directional viability selection on males did not consistently predict changes in trait size in the next generation (Tables 1 and 2). In 2011, significant selection for narrower heads was detected, and 2012 males had significantly narrower heads ( x 2011 = 1.61 mm x 2012 = 1.54 mm, two-tailed t-test, t 66 = 2.49, P = 0.02). However, in 2009, I saw significant selection for smaller pronotum length and leg size, yet in 2010, male pronota and legs were significantly larger ( x 2009 = 2.20 mm x 2010 = 2.37 mm, two-tailed t-test, t 77 = À5.39, P < 0.01).

Yearly viability selection on females
As with males, I found significant interactions between a trait (tegmen width) and sampling year in selection on females, so I analyzed selection separately for each year (Supplementary Table S1b). There were no significant selection differentials or directional selection gradients on any original female traits in any year (Table 3), but canonical analysis did reveal nonlinear selection on composite female traits (Table 4). Although the relative trait loadings for canonical traits changed considerably from year to year (see M-matrices, Table 4), composite traits that mostly represented pronotum length were frequently under significant selection. In 2009, there was marginally significant nonlinear convex selection on composite trait m 4 (k = À0.542, P = 0.06), which represented pronotum length and leg size equally. In 2010, there was significant concave (k = 1.615, P = 0.01) selection on composite trait m 1 (which is strongly influenced by pronotum length) and significantly convex (k = À0.341, P < 0.01) and stabilizing selection (Mitchell-Olds Shaw test, P = 0.03) on m 4 (mostly representing head width). In 2011, I could not conduct a canonical analysis on female traits because only eight adult females were found in I. mexicana nests. In 2012, I found significantly concave (k = 0.810, P = 0.01) and disruptive (P = 0.02) selection on composite trait m 1 , which mostly represents pronotum length.

Discussion
The magnitude of viability selection gradients on both male and female O. nigricornis within years was larger than across years (Table 3). Over the 4-year period, there were no significant directional selection gradients on any original traits in both sexes. Males were subject to significant directional selection within years, but which traits were under selection changed significantly from year to year, resulting in no significant selection gradients on any one trait over the 4-year period. However, a significant selection differential shows that male head width became slightly narrower over the 4 years. Similar temporal variation in selection is commonly observed in long-term selection studies (e.g., Kalisz 1986;Gibbs and Grant 1987;Milner et al. 1999;Punzalan et al. 2010;Siepielski et al. 2011), and this variation may dampen the strength of directional selection (Chaine and Lyon 2008;Siepielski et al. 2009); but also see Morrissey and Hadfield 2012). Yearly variation in selection instead resulted in significant nonlinear selection on both males and females over the 4-year period (Table 4). Males experienced significant convex selection on composite trait m 4 that resulted in males with relatively larger tegmina and smaller legs having a survival advantage against predatory wasps, but this advantage diminishes as the trait value increases (Fig. 2). These results contrast with the findings in Ercit and Gwynne (2015) that males with smaller tegmina and larger legs had a survival advantage in 2012. However, as the form of total selection over 4 years on the tegmen/ legs size trait is convex, male crickets from 2012 may have had larger-than-average traits on axis m 4 and may represent the downward slope seen on the right side of the graph in Figure 2. Females experienced concave selection on composite trait m 1 that indicates a fitness trough for females with intermediate head widths, and females with large and small head widths are more likely to survive wasp predation. This selection appears to be disruptive (Fig. 3), but it is not significantly so. It is not clear why female head width is important in viability selection over the 4 years, especially as I did not detect significant selection on it within any single year. Selection on female head width may be a statistical artifact, or it may be that head width in females is more strongly correlated to body mass than estimated, and concave selection on head width may result from disruptive selection on body mass. Cubic Table 3. Vectors of standardized directional selection gradients (b) (and their associated standard errors) and selection differentials (S) for viability selection on male and female Oecanthus nigricornis over 4 years.

Males
Females Bolded values are significant at a = 0.05. None of the selection differential or gradients were significant, thus none are bolded.
spline analysis and partial F-tests showed that the shape and direction of viability selection were different between males and females within years. However, total selection over the 4-year period was not significantly different between males and females. Such temporal variation in the difference in selection between males and females may slow the evolution of sexual dimorphism. It is interesting that males experienced mostly linear selection, whereas for females it was mainly nonlinear. Linear viability selection on males may be connected to The M-matrix of relative loadings of the original traits on the new canonical axes is also included. Bolded values indicate significance at a = 0.05).  sexual selection on males: The results of Ercit and Gwynne (2015) show that traits that made males successful at mating also made them more likely to be killed by wasps in 2012. If mating per se is risky for males, male traits that attract females will also be subject to viability selection by I. mexicana. If this sexual selection is mostly linear (as it predominantly was in 2012 [Ercit and Gwynne 2015;]), and if predation risk in males increases linearly with mating success, we may expect opposing viability selection on males to also be linear. If sexually attractive males attract more predators, this may also explain why male traits under viability selection change from year to year. In other animals, the male traits that are related to mating success can vary between years (e.g., Hughes et al. 1999;Chaine and Lyon 2008). If this is the case in tree crickets, the traits of successfully mating males would vary between years, and so might the traits of males killed in risky mating behaviors. The observed nonlinear selection on females may be due to biases and limitations of the predator. In 2010 and 2012 females, I saw significant concave selection on composite traits that represent body size (pronotum length), which can be expected based on results from (Ercit 2014): I. mexicana take large prey, but the largest females may be too heavy for the wasp to transport (Marden 1987;Coelho and Ladage 1999). Thus, the observed temporal variation in viability selection may be caused partly by variation in predator size, possibly exacerbated by sampling error, as I only sampled prey from the few dozen wasps that nested in the trap nests each year.
One limitation of this study was that I was not able to estimate what proportion of the population was killed by wasp predation each year. Variation in relative predator and prey populations could significantly affect the intensity of viability selection (Benkman 2013). In the presence of abundant nesting habitat, some solitary sphecids can have significant impacts on prey density (Dukas 2005). If wasps overhunt one prey population, some solitary wasps will switch prey species (Polidori et al. 2007), and indeed, I observed such prey-switching during the course of this study. In 2011, part of the reason why so few O. nigricornis were collected was because most wasps were provisioning other Oecanthus species. Thus, the intensity of predation by wasps on this study population of O. nigricornis likely varied greatly from year to year, and in turn, affected the intensity of selection.
Although males in our population experienced significant directional selection on traits, I did not reliably see significant change in that trait in the next generation. This result is not surprising, as I have only measured selection from one component of fitnessviability selection from a single predator. Tree crickets also experience viability selection from other predators such as spiders and birds, as well as from environmental factors. Furthermore, within a generation, viability selection may be counteracted by fecundity or sexual selection (as in Ercit and Gwynne 2015). Even if total natural selection (selection from every component of fitness) was significantly directional, this selection may be acting on phenotypic variance caused by environmental rather than heritable variables, as was found in a study of collared flycatchers (Alatalo et al. 1990). In another multigenerational study, Milner et al. (1999) found repeated selection for larger body weight in Soay sheep did not result in any change in population mean weight, and this was likely due to selection acting on phenotypic variance caused by the environment.
The results presented here are consistent with other multigenerational studies of selection that show that the direction of viability selection is variable between generations (summarized in Siepielski et al. 2011). The results of this study underline the importance of temporal scale in selection studies. If short-term studies are extrapolated to long-term selection , this may overestimate the rates of evolution and sexual dimorphism (Kinnison and Hendry 2001). Thus, this study adds to our knowledge of how selection acts in different timescales on a natural population and may help in future studies to estimate rates of evolution in nature.