Introduction
The central tenet of life history theory is that organisms have evolved to make strategic allocations of time and energy to essential components of maintenance, growth, and reproduction. These include, among others, the trade-off between growth and reproduction and the trade-off between quality and quantity of offspring (Roff, Reference Roff2002; Stearns, Reference Stearns1992). Although life history theory is ostensibly a theory about the evolution of between-species differences in life history characteristics, such as maturational timing, age at first reproduction, and mating patterns, models of within-species adaptation to environmental exposures during development have been extrapolated from life history theory (see e.g., Belsky et al., Reference Belsky, Steinberg and Draper1991; Ellis et al., Reference Ellis, Figueredo, Brumbach and Schlomer2009). This extrapolation is not without controversy (e.g., Frankenhuis & Netttle, Reference Frankenhuis and Nettle2020; Nettle et al., Reference Nettle, Frankenhuis and Rickard2013; Stearns & Rodrigues, Reference Stearns and Rodrigues2020; Zietsch & Sidari, Reference Zietsch and Sidari2020) though conceptualizing the development of human life history strategies in terms of contingent developmental adaptations has provided a useful heuristic, and a plethora of research, for understanding human trait covariation and their antecedents (see Del Giudice, Reference Del Giudice2020 for a review). For instance, pubertal maturation and the onset of sexual behaviors is a definitive characteristic of adolescent development. There is considerable variability, however, in the timing of pubertal development and sexual behaviors. Some children will experience puberty earlier than others, initiate romantic and sexual behaviors sooner, and may even have more sexual partners during adolescence and young adulthood. Other children will experience later puberty, likely delay sexual behavior onset, and have fewer sexual partners (Deardorff et al., Reference Deardorff, Gonzales, Christopher, Roosa and Millsap2005; Heywood et al., Reference Heywood, Patrick, Smith and Pitts2015; Ibitoye et al., Reference Ibitoye, Choi, Tia, Lee and Sommer2017). From a within-species perspective, the tendency for these sexual development outcomes to cluster together is the hypothesized result of evolved trade-offs that are conditionally sensitive to early developmental experiences (Belsky et al., Reference Belsky, Steinberg and Draper1991; Ellis, Reference Ellis2004; Ellis & Del Giudice, Reference Ellis and Del Giudice2019; Ellis et al., Reference Ellis, Figueredo, Brumbach and Schlomer2009). Importantly, these trade-offs reflect both heritable individual differences as well as conditional adaptations. Regarding the latter, variability in pubertal development and sexual behaviors are hypothesized to result, at least in part, from early experiences that indicate either a faster (accelerated pubertal development, earlier age at first sex, less restrictive mating behavior) or slower (delayed pubertal development, later age at first sex, more discriminating mating behavior) life history strategy is the likely evolutionarily optimal developmental trajectory.
Advancements in the application of life history theory to human development have identified two fundamental dimensions of environmental risk that are key for life history strategy development (Ellis et al., Reference Ellis, Figueredo, Brumbach and Schlomer2009). Given adequate bioenergetic resources, environmental harshness and unpredictability are hypothesized to uniquely influence the development of human life history characteristics. Empirical studies of harshness and unpredictability have generally supported this hypothesis. Specifically, a small corpus of research designed to study the mutual influence of early experiences with harshness (operationalized as income-to-needs/socioeconomic status) and unpredictability (operationalized as paternal transitions, employment transitions, and/or household moves) has shown that these cues are related to sexual development, including pubertal development in girls (age at menarche), onset of sexual intercourse, and sexual partner number (Belsky et al., Reference Belsky, Schlomer and Ellis2012; Hartman et al., Reference Hartman, Sung, Simpson, Schlomer and Belsky2018; Simpson et al., Reference Simpson, Grisskevicius, Kuo, Sung and Collins2012; Sung et al., Reference Sung, Simpson, Griskevicius, Kuo, Schlomer and Belsky2016; though see Colich et al., Reference Colich, Rosen, Williams and McLaughlin2020). None of these studies have tested the hypothesis that pubertal development and age at first sexual intercourse could be mechanisms by which harshness and unpredictability influence sexual partner number (Belsky et al., Reference Belsky, Schlomer and Ellis2012). Numerous studies of life history strategy development focused on related aspects of the early environment, such as maternal harshness (Belsky et al., Reference Belsky, Steinberg, Houts and Halpern-Felsher2010) and father absence (James et al., Reference James, Ellis, Schlomer and Garber2012; Richardson et al., Reference Richardson, La Guardia and Klay2018) suggests this hypothesis is plausible (see Ellis & Del Giudice, Reference Ellis and Del Giudice2019 for a review). Last, studies of harshness and unpredictability have generally not accounted for possible genetic confounding or evaluated different forms of gene-environment interplay, such as gene-by-environment interactions (GxE). Addressing genetics is an important consideration given life history characteristics are a hypothesized product of conditional adaptations and heritable individual differences (Ellis et al., Reference Ellis, Figueredo, Brumbach and Schlomer2009) and failing to account for genetic correlations in this literature can result in spurious associations (Barbaro et al., Reference Barbaro, Boutwell, Barnes and Shackelford2017). The purpose of the current study was to build on prior life history-informed research by examining a mediational model that included developmental cues to harshness and unpredictability, girls’ age at menarche (AAM), age at first sexual intercourse, and adolescent and young adult sexual partner number. In addition, possible genetic confounding among these associations and potential GxE were tested using a polygenic score (PGS; Day et al., Reference Day, Thompson, Helgason, Chasman, Finucane, Sulem, Ruth, Whalen, Sarkar, Albrecht, Altmaier, Amini, Barbieri, Boutin, Campbell, Demerath, Giri, He, Hottenga, Karlsson and Perry2017).
The influence of harshness and unpredictability on sexual development
In a comprehensive analysis of the evolution and development of life history strategies, Ellis et al. (Reference Ellis, Figueredo, Brumbach and Schlomer2009) identified harshness and unpredictability as distinct environmental dimensions that are key for regulating life history strategy development. Environmental harshness refers to sources of extrinsic morbidity and mortality, or factors that increase risk for disability and death. Environmental unpredictability refers to spatial/temporal variation in extrinsic morbidity and mortality. From an evolutionary perspective, delaying maturation and the onset of reproductive behaviors – characteristics of a slower life history strategy - permit an organism to better utilize environmental resources during development (such as parental investment) resulting in increased reproductive competitiveness. However, when environmental harshness is high, the more optimal evolutionary strategy is to accelerate pubertal development and begin reproduction sooner – characteristic of a faster life history strategy – or risk losing the opportunity to reproduce entirely via death or disability. Similarly, early environmental unpredictability can undermine effective prediction of the adult reproductive environment. Developing a faster life history strategy under these conditions represents a more optimal strategy given the risk that investing in a slower life history strategy could ultimately be wasted if adult mortality turns out to be high. Taken together, harshness and unpredictability are theorized to each uniquely influence the development of faster life history strategies, characterized by earlier pubertal maturation, earlier onset of sexual behaviors, and possible increased sexual partner numbers (Ellis et al., Reference Ellis, Figueredo, Brumbach and Schlomer2009).
In a modern societal context, what are relevant cues to harshness and unpredictability? Drawing from contemporary research, perhaps the most salient indicator of extrinsic morbidity/mortality in a Western society is socioeconomic status. Across diverse societal structures, lower socioeconomic status is robustly linked with increased risk of mortality (Debiasi & Dribe, Reference Debiasi and Dribe2020; Elo, Reference Elo2009; Marmot, Reference Marmot2014). Socioeconomic status (and related constructs such as poverty and income inequality) is a major social determinate of health (Marmot, Reference Marmot2014) and is related to other more direct measures of environmental harshness (Colich et al., Reference Colich, Rosen, Williams and McLaughlin2020), such as exposure to violence (Foster et al., Reference Foster, Brooks-Gunn, Martin, Flannery, Vazsonyi and Waldman2007) and harsh child-rearing practices (Belsky, Reference Belsky1993). With regard to unpredictability, research has focused on indicators of instability and change in the early developmental environment, such as changes in parental figures, changes in parental employment, and residential moves given the centrality of the family and parenting on early child development (Belsky et al., Reference Belsky, Schlomer and Ellis2012; Simpson et al., Reference Simpson, Grisskevicius, Kuo, Sung and Collins2012). Empirical studies have generally found support for the propositions that harshness and unpredictability influence the development of faster life history strategies, though the evidence appears stronger for unpredictability. In a study of the number of oral and sexual intercourse partners among 15-year-olds, Belsky et al. (Reference Belsky, Schlomer and Ellis2012) found that exposure to harshness and unpredictability during early development undermined maternal parenting sensitivity via increased maternal depression, which lead to increased sexual partner number among adolescents. Environmental unpredictability, which was comprised of paternal transitions, employment changes, and residential changes, showed a direct association with sexual partner number while harshness (income-to-needs) did not. Simpson et al. (Reference Simpson, Grisskevicius, Kuo, Sung and Collins2012) examined a number of life history related traits, including age at first sex and number of sexual partners by age 23. Using a measure of unpredictability indicated by maternal life stress associated with employment changes, residential changes, and changes in male co-habitation, the authors found that unpredictability forecast more sexual partners and, at least for boys, earlier sexual onset. Harshness, operationalized as socioeconomic status, was not related to these outcomes. To determine which features of early unpredictability might be most salient for predicting life history traits, Hartman et al. (Reference Hartman, Sung, Simpson, Schlomer and Belsky2018) examined paternal transitions, employment changes, and residential changes separately. Using two different samples, the authors found that paternal transitions most strongly and consistently predicted faster life history characteristics, including age at sexual onset and number of sexual partners. Last, in the only study that has jointly examined whether harshness and unpredictability are related to pubertal development, Sung et al. (Reference Sung, Simpson, Griskevicius, Kuo, Schlomer and Belsky2016) found that harshness, not unpredictability, was associated with earlier AAM. Taken together, the existing evidence suggests that unpredictability, and in particular paternal transitions, may be relatively more important for predicting the onset of sexual behaviors and number of adolescent and young adult sexual partners. The evidence also indicates, however, that harshness rather than unpredictability may be a stronger predictor of AAM. This latter conclusion should be considered in light of recent meta-analytic evidence, however, that showed socioeconomic status was not related to pubertal development (male and female) though other aspects of harsh environments, such as abuse, violence, and assault, were related (Colich et al., Reference Colich, Rosen, Williams and McLaughlin2020). Socioeconomic status as an indicator of harshness may not be strongly related to AAM but may influence life-history related traits via other means.
The hypothesis that AAM and age at first sexual intercourse will mediate associations between harshness and unpredictability and sexual partner number stems from two observations. First is the finding that that earlier AAM is robustly related to sexual behaviors (Belsky et al., Reference Belsky, Steinberg, Houts and Halpern-Felsher2010). A recent meta-analysis, for example, found earlier AAM was associated with earlier sexual onset, more sexual behaviors (e.g., petting, oral sex), and increased risky sexual behaviors (e.g., noncondom/contraception use, substance use during sex; Baams et al., Reference Baams, Dubas, Overbeek and van Aken2015). Second, perhaps the most well-studied association in the human development life history literature is the link between biological father absence, AAM, and sexual behavior outcomes. Dozens of studies over a number of years have established a reliable association, with some exceptions due to energetic deprivation (see Sear et al., Reference Sear, Sheppard and Coall2019; Sohn, Reference Sohn2017), between early (the first 5–7 years) biological father absence and younger AAM and sexual behaviors (see Ellis, Reference Ellis2004; Ryan, Reference Ryan2015; Webster et al., Reference Webster, Graber, Gesselman, Crosier and Orozco Schember2014). At least one study has examined AAM as a possible mediator of the relations between socioeconomic status and father absence on sexual onset and sexual risk taking (i.e., James et al., Reference James, Ellis, Schlomer and Garber2012). Results showed that both socioeconomic status and father absence were related to sexual onset and sexual risk taking via a multiple mediated pathway that included (1) their association with family relationship quality during early adolescence and (2) the link between family relationship quality and AAM. Though father absence is not equivalent to unpredictability, it is worth highlighting that paternal transitions, which is linked to biological father absence, appears to be the primary driver of unpredictability associations with life history characteristics (Hartman et al., Reference Hartman, Sung, Simpson, Schlomer and Belsky2018). As noted above, however, the one study that examined the associations between harshness and unpredictability with AAM did not find that unpredictability was linked to earlier AAM. One goal of this study is to further examine these associations and test a model using AAM as a possible mediator of harshness and unpredictability associations with sexual behavior outcomes.
Addressing gene-environment interplay in research on harshness and unpredictability
Putative environmental associations between cues to harshness and unpredictability with life history phenotypes may actually reflect gene-environment correlations (Jaffee & Price, Reference Jaffee and Price2007). Passive gene-environment correlations, for example, occur when parents provide a rearing environment that is correlated with their children’s inherited characteristics, passed down from parent to child. As a result, childhood experiences with harshness and unpredictability may actually reflect a shared genetic liability between parents and children. Genetic factors that might dispose parents toward, for example, unstable relationships, frequent employment transitions, and low socioeconomic status could be passed down to their children, who may then be more generally liable toward earlier AAM and risky sexual behaviors. This proposition has been extensively studied in research on the association between father absence and AAM (e.g., Ellis et al., Reference Ellis, Schlomer, Tilley and Butler2012; Mendle et al., Reference Mendle, Turkheimer, D’Onofrio, Lynch, Emery, Slutske and Martin2006, Reference Mendle, Harden, Turkheimer, Van Hulle, D’onofrio, Brooks-Gunn, Rodgers, Emery and Lahey2009; Tither & Ellis, Reference Tither and Ellis2008) but less so in research on harshness and unpredictability. One exception is research by Sung et al. (Reference Sung, Simpson, Griskevicius, Kuo, Schlomer and Belsky2016) who used maternal AAM as a partial genetic control for the association between harshness and unpredictability with AAM. Although such studies have provided insight into the possible unique effects of early environmental experiences by controlling or partially controlling for genetic confounding, additional research is needed to examine other possible forms of gene-environment interplay that could be important for life history strategy development, such as GxE.
As genotyping costs have decreased over time, a growing body of research has begun to include DNA derived genotypes in research on links between early experiences, AAM, and sexual behaviors. Much of this work has relied on candidate gene approaches, which reflects an important first step in integrating DNA information into this area of research (e.g., Hartman et al., Reference Hartman, Widaman and Belsky2015; Schlomer & Cho, Reference Schlomer and Cho2017; Schlomer et al., Reference Schlomer, Murray, Yates, Hair and Vandenbergh2019). This approach to genetic research has several known limitations, however, such as low statistical power (see Dick et al., Reference Dick, Agrawal, Keller, Adkins, Aliev, Monroe, Hewitt, Kendler and Sher2015) and the biological reality that complex phenotypes, such as AAM and sexual behaviors, are not single-gene phenotypes and are underlain by many genes of small effect (Fisher, Reference Fisher1918; Lohmueller et al., Reference Lohmueller, Pearce, Pike, Lander and Hirschhorn2003). As a result of these and other limitations, genomic studies in human development have moved toward utilizing genomic aggregates, the most powerful of which are PGSs, which are based on findings from genome-wide association studies (GWAS; Belsky & Israel, Reference Belsky and Israel2014). By adding many variants of small effect together, PGSs can predict more phenotypic variance than any individual variant on its own. Recent primary and replication research using a PGS consisting of AAM predictive variants (Day et al., Reference Day, Thompson, Helgason, Chasman, Finucane, Sulem, Ruth, Whalen, Sarkar, Albrecht, Altmaier, Amini, Barbieri, Boutin, Campbell, Demerath, Giri, He, Hottenga, Karlsson and Perry2017) found that the father absence/AAM association could not be explained by genetic variability captured in the PGS (Gaydosh et al., Reference Gaydosh, Belsky, Domingue, Boardman and Harris2018; Schlomer & Marceau, Reference Schlomer and Marceau2021); nor was the association moderated by the PGS, reflecting a lack of GxE (Schlomer & Marceau, Reference Schlomer and Marceau2021). Potential genetic confounding or GxE of harshness and unpredictability associations have not been tested using genomic information in general or the AAM PGS specifically. Further, it is unknown if links between AAM, age at first sexual intercourse, and sexual partner number are driven by genetic variants specific to AAM.
The current study
The purpose of the current study was to build on prior research on harshness and unpredictability by testing whether cues to these factors are directly related to adolescent and young adult sexual behaviors and/or indirectly related via AAM and age at first sexual intercourse. Further, this study addressed potential genetic confounding and GxE among these associations using a PGS comprised of AAM-related genetic variants (Day et al., Reference Day, Thompson, Helgason, Chasman, Finucane, Sulem, Ruth, Whalen, Sarkar, Albrecht, Altmaier, Amini, Barbieri, Boutin, Campbell, Demerath, Giri, He, Hottenga, Karlsson and Perry2017). To do so, a path model was conducted that included the mutual influence of cues to harshness and unpredictability on AAM, age at first sexual intercourse, and lifetime and past year sexual partner number. AAM and age at first sexual intercourse were included in the model as mediators of harshness and unpredictability outcomes enabling both direct and indirect effects of harshness and unpredictability to be evaluated. Potential genetic confounding of these associations was addressed by subsequently including the AAM PGS as a covariate in the model. Last, GxE was tested by determining if path model associations differed across levels of the PGS.
Methods
Participants
Data for this study were drawn from the Avon Longitudinal Study of Parents and Children (ALSPAC). ALSPAC is a population-based, longitudinal cohort study that initially enrolled 14,541 pregnant women in the United Kingdom whose expected date of delivery was between April 1st, 1991 and December 31st, 1992 (Boyd et al., Reference Boyd, Golding, Macleod, Lawlor, Fraser, Henderson, Molloy, Ness, Ring and Davey Smith2013; Fraser et al., Reference Fraser, Macdonald-Wallis, Tilling, Boyd, Golding, Davey Smith, Henderson, Macleod, Malloy, Ness, Ring, Nelson and Lawlor2013). Of these initial pregnancies, there was a total of 14,676 fetuses, resulting in 14,062 live births and 13,988 children who were alive at 1 year of age. Data collection began when mothers were 8 weeks pregnant and includes 68 data collection time points between birth and child age 18. An additional 913 children were enrolled over the course of the study. Thus, the total sample available after age 7 years included 15,454 pregnancies, resulting in 15,589 fetuses of which 14,901 were alive at 1 year of age. Written informed consent was obtained from all study participants. Ethical approval for the study was obtained from the ALSPAC Ethics and Law Committee and the Local Research Ethics Committees. Consent for biological samples has been collected in accordance with the Human Tissue Act (2004). Approval for this study was obtained from the University at Albany, SUNY Institutional Review Board.
Specific to the current study, the ALSPAC sample included N = 7,262 cases who self-identified as female. Of the female sample, genome-wide data were available on N = 4,355 participants. By design, the ALSPAC genomic sample is entirely European ancestry. Cases were excluded from the current study if they did not have a score on the PGS (e.g., genotyping failures, quality controls, N = 92) or if they were missing data across model variables (N = 618). The final analytic sample for the current study was N = 3,645 female participants.
Measures
Please note that the ALSPAC study website contains details of all the data that are available through a fully searchable data dictionary and variable search tool: http://www.bristol.ac.uk/alspac/researchers/our-data/.
Harshness – Financial difficulties
When children were 8, 33, and 85 months of age mothers responded to a series of questions about difficulty affording food, clothing, rent, heat, and things for the child (0 = very difficult, 3 = not difficult). Item reliabilities at each wave were good (α’s = .889 to .894) and items were averaged within wave. Mean financial difficulties were .594 (SD = .700), .546 (SD = .663), and .326 (SD = .470) and ranged from 0 to 3 across waves, respectively. The measures were moderately correlated (mean inter-item correlation = .534). The three assessments were averaged and coded so that higher values equal more financial difficulties. Of the N = 3,645 analytic sample, N = 3,382 (92.78%) cases had data on financial difficulties. The composite scale ranged from 0 to 3 (M = .512, SD = .546).
Unpredictability – Paternal transitions
Environmental unpredictability was operationalized as paternal transitions, which previous research indicates is more consistently and strongly related to life history characteristics relative to other aspects of unpredictability (i.e., residential and employment transitions; Hartman et al., Reference Hartman, Sung, Simpson, Schlomer and Belsky2018). Information about whether the mother had a male live-in partner was assessed when the mother was pregnant and when the target child was 1 year and 9 months, 2 years and 9 months, 3 years and 11 months, 6 years and 1 month, and 7 years and 1 months (6 assessments total). To create a paternal transitions variable, cases were coded at each assessment as 1 = Yes, male live-in partner or 0 = No male live-in partner (or no partner at all). Starting with the prenatal assessment, a paternal transition was counted if the mother’s partner status changed between adjacent time points (5 total transition points). When data were missing on a prior adjacent time point, the most recent available data was used. For example, when assessing if there was a transition between time 4 (3 years, 11 months) and time 5 (6 years, 1 month), if data were present for time 5 but not time 4, data from time 3 (3 years, 11 months) were used. Using this coding, a transition occurred if the mother went from having a male live-in partner to not having a male live-in partner and if the mother went from not having a live-in male partner to having a live-in male partner. The paternal transitions variable thus codes for the entrances and exits of male partners into and out of the household during the child’s first 7 years of life. To help account for additional transitions between time points, information about the length of the relationship between the mother and the live-in male partner was also used. For example, if a case was coded has having no transition between time 1 (birth) and time 2 (1 year, 9 months), but the length of the relationship between the mother and the male live-in partner was less than 1 year and 9 months, the case was recoded to reflect that a transition had occurred. In addition, if at one assessment the male live-in father figure was identified as the biological father and at the subsequent assessment the male live-in father figure was not the biological father, a transition was counted.
The five transition variables were averaged and the resulting paternal transition variable ranged from 0 to 1. Based on the available data (N = 2,721, 74.7%), 13.6% of the sample children experienced at least one transition in their first 7 years of life (M = .05, SD = .13).
Age at menarche
Between the ages of 8 and 17 years of age, daughter’s AAM was prospectively assessed up to nine times during approximately annual assessments. Surveys addressed if menstruation had begun and if so at what age. AAM responses from the first available wave were used in the current study to minimize recall bias (Culpin et al., Reference Culpin, Heron, Araya, Melotti, Lewis and Joinson2014; Gaydosh et al., Reference Gaydosh, Belsky, Domingue, Boardman and Harris2018; Schlomer & Marceau, Reference Schlomer and Marceau2021). Based on the N = 3,011 available cases (82.6%), mean AAM in months was 152.00 (SD = 14.03) and ranged from 91 to 203.
Lifetime and past year sexual partners
ALSPAC adolescents completed sexual behavior assessments at age 17.5, which included questions about whether they have ever had sexual intercourse and the number of people they have had sexual intercourse with in their lifetime and within the past year. For lifetime sexual partners (N = 1,800, 49.38%), adolescents were coded as 0 if they had never had sex. Within the current analytic sample, mean lifetime sexual partners was 2.27 (SD = 3.19) and ranged from 0 to 20 (20 was the maximum number they could report). To reduce the influence of outliers and skew, this variable was recoded such that 0 = Never had sexual intercourse, 1 = 1 to 2 partners, 2 = 3 to 5 partners, 3 = 6 to 10 partners, and 4 = more than 10 partners (M = 1.12, SD = 1.03; range = 0 to 4). A similar approach was taken for past year sexual partners (N = 1,802, 49.44%), which ranged from 0 to 20 partners (M = 1.26, SD = 1.73). This variable was recoded so that 0 = None or never had sex, 1 = 1 partner, 2 = 2 partners, 3 = 3 to 4 partners, and 4 = 5 or more partners (M = 1.10, SD = 1.08; range = 0 to 4).
Age at first sex
Age at first sex (N = 2,223, 60.99%) was collected retrospectively at age 21 and at age 23 (Northstone et al., Reference Northstone, Lewcock, Groom, Boyd, Macleod, Timpson and Wells2019). To maximize sample size and reduce missing data, responses from both assessments, which were highly correlated (r = .94), were used to create an age at first sex variable. Missing data on the age 23 assessment was recovered by using data from the age 21 assessment. The age 23 assessment served as the “base” measure since it includes more cases with later ages at first sex and is therefore less right-censored. Mean of at first sex was 16.42 years (SD = 1.86; range = 12–24). Similar to prior studies (e.g., Bingham & Crockett, Reference Bingham and Crockett1996; James et al., Reference James, Ellis, Schlomer and Garber2012; Schlomer et al., Reference Schlomer, Murray, Yates, Hair and Vandenbergh2019) age at first sex was binned to help account for right censoring. To preserve as much continuity as possible, age at first sex was coded such that 1 = 12–13 years, 2 = 14–15 years, 3 = 16–17 years, 4 = 18–19 years, 5 = 20–21 years, 6 = 22–24 years (only one case was at 24), and 7 = never had sex (M = 3.03, SD = 1.11, range = 1–7).
Genotyping
DNA were collected on approximately 11,000 children at age 7 using blood, cell line, and mouthwash samples and genotyped using the Illumina HumanHap550 platform (Boyd et al., Reference Boyd, Golding, Macleod, Lawlor, Fraser, Henderson, Molloy, Ness, Ring and Davey Smith2013; Pembrey, Reference Pembrey2004). Raw genotype data were subjected to quality control measures and participants were excluded based on gender mismatch, minimal or excessive heterozygosity, >3% missingness, cryptic relatedness (IBD > .10), and non-European ancestry. Genomic data were imputed using 1000 Genomes phase 1, version 3. Following quality controls, a total of N = 4,355 female participants with genome-wide data were available.
Polygenic score
The PGS used in the current study was derived from findings by Day et al. (Reference Day, Thompson, Helgason, Chasman, Finucane, Sulem, Ruth, Whalen, Sarkar, Albrecht, Altmaier, Amini, Barbieri, Boutin, Campbell, Demerath, Giri, He, Hottenga, Karlsson and Perry2017), who conducted the largest (GWAS) to date on AAM, comprising a sample of over N = 300,000 women from the ReproGen Consortium, 23andME, and the UK Biobank. GWAS results revealed 389 non-redundant genome-wide significant (p < .05 × 10−8) signals that collectively explained approximately 7.0% of the variance in AAM (see Day et al., Reference Day, Thompson, Helgason, Chasman, Finucane, Sulem, Ruth, Whalen, Sarkar, Albrecht, Altmaier, Amini, Barbieri, Boutin, Campbell, Demerath, Giri, He, Hottenga, Karlsson and Perry2017).
Supplemental information provided by Day et al. (Reference Day, Thompson, Helgason, Chasman, Finucane, Sulem, Ruth, Whalen, Sarkar, Albrecht, Altmaier, Amini, Barbieri, Boutin, Campbell, Demerath, Giri, He, Hottenga, Karlsson and Perry2017) were used to create the PGS in the ALPAC data. Of the 389 independent SNPs, N = 372 were available in ALSPAC. Following quality control procedures that account for excess missingness (>1%), low minor allele frequency (<1%), and Hardy−Weinberg violations (p < 5 × 10−7), 339 variants remained. The 339 variants were then weighted by their respective effect sizes (Dudbridge, Reference Dudbridge2013) provided by Day et al., and averaged, resulting in a normally distributed PGS. The raw PGS variable was multiplied by 1000 to aid path model convergence (M = −5.30, SD = 1.24; range = −9.27 to −.92)
Population stratification
A unique issue in molecular genetic research is population stratification, which refers to naturally occurring allele frequency differences among populations with different genomic ancestries. When both phenotype and genotype are correlated with genetic ancestry, population stratification can lead to spurious results (Cardon & Palmer, Reference Cardon and Palmer2003). Prior work on the ALSPAC genetic sample showed little population structure (Jones et al., Reference Jones, Stergiakouli, Tansey, Hubbard, Heron, Cannon, Holmans, Lewis, Linden, Jones, Davey Smith, O’Donovan, Owen, Walters and Zammit2016; Okbay et al., Reference Okbay, Beauchamp, Fontana, Lee, Pers, Rietveld, Turley, Chen, Emilsson, Meddens, Oskarsson, Pickrell, Thom, Timshel, de Vlaming, Abdellaoui, Ahluwalia, Bacelis, Baumbach, Bjornsdottir and Benjamin2016). However, to help adjust the results for potential additional undetected population structure, principal coordinates analysis (PCA) was conducted on the genome-wide data using the -- pca command in plink (Chang et al., Reference Chang, Chow, Tellier, Vattikuti, Purcell and Lee2015). Eigenvalues examined in a scree plot (Cattell, Reference Cattell1966) leveled off following the first PC (PC1), suggesting the first PC explained the most variance and that subsequent PCs may reflect overfitting (i.e., error; D’Agostino & Russell, Reference D’Agostino, Russell, Armitage and Colton2005).
Analysis plan
Path models were conducted in R using the lavaan package (Rosseel, Reference Rosseel2012). Maximum likelihood estimation with robust standard errors (MLR) was used to estimate the models, which helps account for variable nonnormality. Full information maximum likelihood within lavaan was implemented to handle missing dataFootnote 1 . Across models, residual covariances were modeled among dependent variables (i.e., age at first sex, lifetime and past year sexual partners). To test the genetic confounding hypothesis, the PGS was subsequently added to the model as a predictor of all variables. Population stratification was addressed by similarly including PC1 within the PGS confound model. Last, a multiple-group path model was conducted that compared participants above and below the PGS median. Between group differences were tested using equality constraints on individual path coefficients and evaluating the chi-square difference between the unconstrained and constrained models. The Yuan−Bentlar correction was used to adjust chi-squares for difference testing when estimating models via MLR (Yuan et al., Reference Yuan, Bentler and Zhang2005).
Results
Preliminary analyses
As can be seen in Table 1, correlations among model variables were generally small (i.e., r < .20) with the exception of the intercorrelations among age at first sex, lifetime, and past year sexual partners, which ranged from .823 to −.536. The only other relatively large correlation was between AAM and the PGS (r = .291). Earlier AAM was related to earlier age at first sex (r = .134) and more lifetime (r = −.101) and past year (r = −.080) sexual partners.
Correlation matrix of analysis variables
Note: LT = lifetime, PY = past year; PC1 refers to the first principal coordinate.
*p < .05.
Both harshness and unpredictability were significantly correlated, in the predicted directions, with all model variables, expect the PGS and PC1 (see Table 1). The PGS was only related to AAM and PC1 showed no associations with any model variables.
The influence of harness and unpredictability cues on AAM on sexual behaviors
As can be seen in Figure 1, a mediational path model was conducted that included lifetime and past year sexual partners as endogenous variables, predicted by age at first sexual intercourse, AAM, harshness, and unpredictability. To determine if there were indirect effects of harshness, unpredictability, and AAM on sexual partner number, harshness, unpredictability, and AAM were additionally used as predictors of age at first sex. Last, harshness and unpredictability were used as exogenous predictors of AAM. The initial model included all possible paths and covariances, resulting in a fully saturated model. Results showed that AAM did not directly predict lifetime (β = −.007, ns) or past year (β = .003, ns) sexual partners nor did paternal transitions (β = .032, ns and β = .013, ns, respectively). These paths were dropped from the model and the model was estimated a second time. Dropping these paths did not appreciably result in model misfit (χ2(5) = 4.83, ns, CFI = 1.00, TLI = 1.00; RMSEA = 0.00, 90% CI = .00 – .024) and the overall model showed good fit to the data. Significant paths depicted in Figure 1 were taken from the second model.
Direct and indirect effects of harshness, unpredictability, AAM, and age at first sexual intercourse on adolescent and young adult girls’ sexual partner number. The dashed lines represents a nonsignificant paths. Covariances among sexual partner variables modeled but not shown. Coefficients are standardized and significant paths were taken from a model that excluded the nonsignificant paths. *p < .05.
Inspecting paths in Figure 1 reveals that earlier age at first sexual intercourse was strongly associated with more lifetime (β = −.615) and past year sexual partners (β = −.531). After controlling for this association, only early harshness showed direct associations with lifetime (β = .069) and past year (β = .051) sexual partners. Though the effect sizes were quite small, more harshness during the first 5–7 years of life forecast more sexual partners during adolescence. Harshness was also related to earlier age at first sexual intercourse (β = −.102) as was and more unpredictability (β = −.110) and earlier AAM (β = .123). Last, harshness showed a small association with AAM (β = −.054) such that more harshness early in development was related to earlier AAM.
Both harshness and unpredictability showed small indirect effects on lifetime (β = .062 and .068, respectively) and past year (β = .054 and β = .058, respectively) sexual partner number via age at first sexual intercourse. Harshness also showed a multiple mediated indirect association via AAM and age at first sex (β = .004). The total effect of harshness (direct plus indirect effects) on lifetime and past year sexual partner number was β’s = .135 and .109, respectively. Indirect effects of AAM on lifetime and past year sexual partner number via age at first sex were β = .076 and β = .065, respectively.
Taken together, these findings indicate that only harshness had direct and indirect effects on sexual partner number. Harshness showed additional indirect effects via both age at first sexual intercourse and AAM. In addition, harshness, unpredictability, and AAM were each directly and uniquely related to age at first sex and both unpredictability and AAM showed indirect effects (only) on sexual partner number via age at first sexual intercourse. Notably, only harshness significantly predicted earlier AAM (and not unpredictability). Associations among model variables tended to be small, however.
PGS confounding and GxE
To test whether the PGS could account for the associations depicted in Figure 1, an additional model was conducted wherein the PGS was added as an exogenous predictor of all model variables. To help account for possible population stratification, PC1 was likewise included as a predictor of model variables. Results showed the PGS significantly predicted AAM (β = .290, p < .05) and was uncorrelated with all other model variables (β’s ranged from −.023 to .009, ns). As suggested by the raw correlations (see Table 1), PC1 was uncorrelated with all model variables (β’s ranged from −.005 to .038, ns).
With regard to the path coefficients, adding the PGS (and PC1) did little to diminish model associations, all of which remained statistically significant (exact coefficients can be seen in Table 2, PGS Control). The largest effect size change was observed for the association between harshness and AAM, which was reduced by .007 (from −.054 from −.047; a 13.0% reduction). All other associations changed by ≤ |.005| or stayed the same. These findings indicate that the PGS could not substantially account for the significant associations depicted in Figure 1.
Coefficient comparison for PGS confound and multiple group path models
Note: AAM = age at menarche, LT = lifetime, PY = past year. Table values are standardized betas (β). Δχ 2(1) is for the coefficient difference between Below/Above Median.
*p < .05, † p < .10.
Last, to determine whether the PGS moderated any of the significant model associations, a multiple group analysis was conducted comparing participants below (N = 1,822) and above (N = 1,823) the PGS median. Model coefficients for each group as well as their respective chi-square difference tests are provided in Table 2. Overall, no significant differences were detected between participants above or below the PGS median, though it is worth noting that associations between harshness and AAM, lifetime, and past year sexual partners were significant among participants above the median and not significant among those below the median. Although there was not a significant main effect of unpredictability on AAM, potential moderation of this association by the PGS was also tested, though not detected (Δχ 2(1) = .25, ns). A marginal group difference (p < .10) was found for the harshness AAM association such that the regression coefficient was somewhat stronger among participants above the PGS median compared to below.
Discussion
The purpose of the current study was threefold. First, based on prior research and theory, a model was evaluated to determine if AAM and age at first sexual intercourse, key facets of women’s life history strategy development, mediated associations between cues to harshness (financial difficulties) and unpredictability (paternal transitions) used in previous research (e.g., Belsky et al., Reference Belsky, Schlomer and Ellis2012) with number of lifetime and past year sexual partners. Second, because putative environmental associations between harshness and unpredictability with life history characteristics might alternatively reflect gene-environment correlations, an AAM PGS was used to determine if genetic confounding could account for any of the model associations. Third, because other forms gene-environment interplay could be important for these associations, possible GxE was tested using the AAM PGS. Results of the path model (Figure 1) showed that only harshness was directly related to more lifetime and past year sexual partners; both harshness and unpredictability were directly related to earlier age at first sexual intercourse; and that only harshness was related to earlier AAM. These associations were notably small. In addition, earlier age at first sexual intercourse strongly predicted more lifetime and past year sexual partners and earlier AAM modestly predicted earlier age at first sex. Both harshness and unpredictability showed small indirect effects on lifetime and past year sexual partner number, largely mediated by age at first sexual intercourse. Tests of genetic confounding showed that the AAM PGS could not account for any of the model associations. Tests of GxE were suggestive of possible stronger associations among participants above the PGS median compared to below, though no significant differences were detected.
Partial support for AAM and age at first sexual intercourse as mechanisms for harshness and unpredictability
Perhaps the most interesting finding from the path model was the lack of association between unpredictability, operationalized as paternal transitions, and AAM. Because decades of research have linked biological father absence to earlier AAM, it might be expected that paternal transitions would likewise be associated with earlier AAM. However, research by Sung et al. (Reference Sung, Simpson, Griskevicius, Kuo, Schlomer and Belsky2016) who, using a measure of unpredictability that also included residential and employment changes, also did not find an association between unpredictability and AAM.
Indeed, a prior study using the ALSPAC data, the focal sample of the current study, did detect the predicted father absence/AAM association (Schlomer & Marceau, Reference Schlomer and Marceau2021). What could account for this difference? The current findings and those by Sung et al. (Reference Sung, Simpson, Griskevicius, Kuo, Schlomer and Belsky2016) highlight the fact that unpredictability, paternal transitions, and biological father absence, though intercorrelated, capture different dimensions of early development. Consider for example, a child whose biological father was absent from birth and whose mother never co-habitated with another male during early development. Also consider a child whose biological father was present from birth and remained in the home throughout development. If indexing father absence, the former child would be considered father absent and the latter father present. If indexing paternal transitions, these two children would be identical, both experiencing zero paternal transitions during development. Although children who are father absent may also experience paternal transitions, the current findings suggest that absence of the biological father, perhaps in a dose-dependent manner, may be more impactful for pubertal development than transitions of father figures into and out of the household. In line with this possibility, paternal investment theory (see Ellis, Reference Ellis2004) emphasizes the importance of paternal care as a key regulator of pubertal timing. It may be that even intermittent paternal care, as a function of paternal transitions, may buffer the impact of biological father absence on accelerated AAM. Research is needed to test this hypothesis.
While paternal transitions showed no significant association with AAM, a small association was detected between harshness – operationalized as financial difficulties – and AAM. Consistent with a recent meta-analysis, however, that found no association between SES and AAM (Colich et al., Reference Colich, Rosen, Williams and McLaughlin2020) the effect size observed in this study could easily go undetected (not significant) in studies with smaller samples than the one currently used. Combined with the observation that a somewhat larger (though still small in absolute terms) association was found for age at first sexual intercourse that was similarly found for unpredictability, this finding suggests that harshness and unpredictability, as operationalized in this study, has a larger impact on age at first sexual intercourse than on AAM. The primary mechanism by which harshness and unpredictability influenced sexual partner number was via age at first sexual intercourse. This pattern of results is somewhat similar to that reported by James et al. (Reference James, Ellis, Schlomer and Garber2012) wherein SES and father absence effects on sexual risk taking were mediated by timing of sexual debut; direct effects on AAM were not detected when the quality of family relationships were taken into account. With regard to specific cues to harshness and unpredictability, stronger tests of life history strategy development should include more proximate mechanisms of these associations given the culmination of evidence suggests direct associations are likely to be relatively week. Along these lines, the meta-analysis by Colich et al. (Reference Colich, Rosen, Williams and McLaughlin2020) did find that more proximal and specific aspects of harshness, such as violence exposure and child maltreatment, were associated with AAM. Future studies explicitly focused on the mutual influence of harshness and unpredictability should carefully consider how these constructs are operationalized and the implications for interpreting results.
Results from the path model also indicate that when controlling for age at first sexual intercourse harshness continued to a show small direct association with lifetime and past year sexual partner number. These findings suggest that, in addition to age at first sexual intercourse, there are other mechanisms that need to be identified, particularly with regard harshness. These mechanisms are likely many and include moderators as well. Though not intended to be a comprehensive list, some possibilities based on prior literature include parenting behaviors and quality of family relationships (James et al., Reference James, Ellis, Schlomer and Garber2012; Warren & Barnett, Reference Warren and Barnett2020), externalizing behavior problems and substance use (Doom et al., Reference Doom, Vanzomeren-Dohm and Simpson2016), self-regulation and cognitive abilities (Li et al., Reference Li, Liu, Hartman and Belsky2018), and attachment (Sung et al., Reference Sung, Simpson, Griskevicius, Kuo, Schlomer and Belsky2016). Additional, research on the possible set of mechanisms that link harshness with the development of life history characteristics – including specific measurement characterizations – are needed to test these and other possible intervening mechanisms.
Harshness, unpredictability, and AAM associations with sexual onset and sexual partner number are not genetically confounded
Adding the AAM PGS (as well as PC1) did little to effect the model associations. These findings provide initial evidence, from the standpoint of measured DNA, that the evaluated associations are not genetically confounded. This conclusion should be contextualized, however, by highlighting that these associations held while controlling for genetic variation that is specific to AAM. It is possible that additional, unmeasured genetic variation could account for these associations. For example, the AAM PGS showed modest to no correlations with age at first sexual intercourse and sexual partner number. It is possible that a PGS developed specifically for predicting age at sexual onset or sexual partner number might also be correlated with AAM and possibly account for these associations. Along these lines, the AAM PGS explained 8.5% of the variation in AAM in the current sample (.2912). Behavioral genetic estimates of AAM heritability suggest that, on average, approximately 50% of the variability in AAM is explained by genetic variation (Towne et al., Reference Towne, Czerwinski, Dermerath, Blangero, Roche and Siervogel2005). This means that the present PGS only captures about 17% of the total genetic variance in AAM. The tendency for variance explained by PGSs to fall short of heritability estimates is common, and known as the missing heritability problem (Manolio et al., Reference Manolio, Collins, Cox, Goldstein, Hindorff, Hunter, McCarthy, Ramos, Cardon, Chakravarti, Cho, Guttmacher, Kong, Kruglyak, Mardis, Torimi, Slatkimn, Valle, Whittemore, Boehnke and Visscher2009). Several reasons for this discrepancy have been proposed though one possibility is that current genome-wide approaches are not yet powerful enough to adequately identify all of the relevant genetic information for a complex phenotype. It is possible that the PGS used in this study does an adequate job of capturing genetic variation unique to AAM but not components of the variation that are shared between multiple phenotypes, such as between AAM and sexual behaviors. Additional research using even more powerful PGS generative methods will be needed to test this possibility.
No evidence for GxE among harshness, unpredictability, and AAM associations with sexual onset and sexual partner number
Tests of PGS moderation revealed that associations did not differ significantly between participants below compared to above the PGS median. It is notable, however, that associations between harshness and AAM and sexual partner number were significant among participants above but not below the PGS median. Though any interpretation of these possible differences should be met with caution, they may reflect an overall tendency for genetic effects to be suppressed under environmental adversity (Carlson et al., Reference Carlson, Mendle and Harden2014). For example, the association between harshness and AAM was more strongly negative among participants who were above the PGS median. This association was driven by participants above the PGS median who showed later AAM when exposed to lower harshness and earlier AAM when exposed to higher harshness. Because higher values on the PGS are related to later AAM, exposure to low harshness may permit better expression of the higher number of variants linked to later AAM (by virture of being higher on the PGS). Nonetheless, these null findings may reflect a more general difficulty with detecting GxE using a PGS approach. By definition, PGSs are aggregates of genetic variants that each show a strong-enough association with a given phenotype to reach statistical significance at a given sample size. Owing to the very small effects typically found for any individual genetic variant, hundreds of thousands or more participants are needed to achieve the statistical power required for genome-wide methods. With increasing sample size also comes increasing sample heterogeneity, the sources of which are both known and unknown, that can result in diminished effect sizes if associations vary across different types and combinations of heterogeneity. Genetic variants that are detected via genome-wide association may be biased toward variants that are indifferent to heterogeneity (genetic penetrance). As a result, it may be difficult to detect moderation using PGSs that, by design, include variants that have limited conditionality.
Limitations and conclusions
Limitations of the current study include the ancestrally homogenous sample, which was comprised entirely of self-identified Caucasians from the UK. Generalizability of the current findings may be limited to this group. An ancestrally homogenous sample can be viewed as a strength, however, in genetic research given the findings are less likely to be a spurious result of population stratification. Combined with the ancestral control used in the PGS analysis (i.e., PC1), the current genetic findings are unlikely related to population stratification. It should also be pointed out that ALPSAC was part of the initial GWAS discovery sample used by Day et al. (Reference Day, Thompson, Helgason, Chasman, Finucane, Sulem, Ruth, Whalen, Sarkar, Albrecht, Altmaier, Amini, Barbieri, Boutin, Campbell, Demerath, Giri, He, Hottenga, Karlsson and Perry2017) to create the AAM PGS. As a result, PGS associations, particularly with AAM, could be inflated in this current study. This issue was addressed in a similar study and the degree of bias was found to be minimal (Horvath et al., Reference Horvath, Knopik and Marceau2019). In addition, the choice to test PGS moderation by comparing participants below and above the median could be seen as a limitation. However, this approach was chosen given its simplicity and potential trade-offs between conducting a single, multiple group model versus several additional models that each include a novel product term. Beyond potential issues with model convergence in product-term path models, the additional analyses could be viewed as increasing alpha inflation. Last, although this study was couched in life history theory using the harshness and unpredictability dimensions, the operationalization of unpredictability was limited to paternal transitions. It is possible that, for example, the null association between unpredictability and AAM would have been significant if additional components of unpredictability, which presumably would increase the effect size, were included in the measure. At least one study suggests this may not be the case (Sung et al., Reference Sung, Simpson, Griskevicius, Kuo, Schlomer and Belsky2016), though additional research is warranted.
Acknowledgments
We are extremely grateful to all the families who took part in this study, the midwives for their help in recruiting them, and the whole ALSPAC team, which includes interviewers, computer and laboratory technicians, clerical workers, research scientists, volunteers, managers, receptionists, and nurses.
Funding statement
The UK Medical Research Council and Wellcome (Grant ref: 217065/Z/19/Z) and the University of Bristol provide core support for ALSPAC. This publication is the work of the authors and GLS will serve as guarantor for the contents of this paper. A comprehensive list of grants funding is available on the ALSPAC website (http://www.bristol.ac.uk/alspac/external/documents/grant-acknowledgements.pdf). GWAS data was generated by Sample Logistics and Genotyping Facilities at Wellcome Sanger Institute and LabCorp (Laboratory Corporation of America) using support from 23 and Me.
Conflicts of interest
None.
