Early life telomeres are influenced by environments acting at multiple temporal and spatial scales

An individual's telomere length early in life may reflect or contribute to key life history processes sensitive to environmental variation. Yet, the relative importance of genetic and environmental factors in shaping early life telomere length is not well understood as it requires samples collected from multiple generations with known developmental histories. We used a confirmed pedigree and conducted an animal model analysis of telomere lengths obtained from nestling house sparrows (Passer domesticus) sampled over a span of 22 years. We found significant additive genetic variation for early life telomere length, but it comprised a relatively small proportion (9%) of the total biological variation. Three sources of environmental variation were important: among cohorts, among breeding attempts within years and families, and among nestmates. The magnitude of variation among breeding attempts and among nestmates also differed by cohort, suggesting that interactive effects of environmental factors across time or spatial scales were important, yet we were unable to identify the specific causes of these interactions. The mean amount of precipitation during the breeding season positively predicted telomere length, but neither weather during a given breeding attempt nor date in the breeding season contributed to an offspring's telomere length. At the residual level, level of individual nestlings, offspring sex, size, and mass at 10 days of age also did not predict telomere length. Environmental effects appear especially important in shaping early life telomere length in some species, and an array of complex more focus on how environmental factors that interactions across scales may help to explain some of the variation observed among studies.Funding provided by: National Science FoundationCrossref Funder Registry ID: http://dx.doi.org/10.13039/100000001Award Number: IBN-9816989Funding provided by: National Science FoundationCrossref Funder Registry ID: http://dx.doi.org/10.13039/100000001Award Number: IBN-0542097Funding provided by: National Science FoundationCrossref Funder Registry ID: http://dx.doi.org/10.13039/100000001Award Number: IOS-1257718Funding provided by: National Science FoundationCrossref Funder Registry ID: http://dx.doi.org/10.13039/100000001Award Number: IOS-1656194Funding provided by: National Science FoundationCrossref Funder Registry ID: http://dx.doi.org/10.13039/100000001Award Number: IOS-1656212Study population: This study made use of a library of blood samples collected from 1993 to 2014 from a population of house sparrows located at the University of Kentucky's North Farm Agricultural Research Station. Westneat et al. (2002, 2009) and Heidinger et al. (2021) provide details about this population and the standard methods used in the field to mark individuals, monitor their reproduction and the collection of and storage of blood samples. Briefly, all nestlings were given uniquely numbered metal bands at 10 days of age, and both independent juveniles and adults were given unique combinations of colored-plastic bands during trapping sessions throughout the year. The present data consists of the families of breeding adults who had hatched in focal nest boxes, were banded as nestlings, and had a blood sample taken when they were 10 days old. While we had samples from any nesting attempt of these focal birds that reached 10 days of age in one of our boxes, we often targeted a subset of these. If the focal adult produced only one brood, we analyzed the DNA of all offspring. If they had multiple broods but only within a single season, we analyzed the brood with the most offspring. If they had broods across more than one season, we sampled the largest brood in each season. Some additional broods were analyzed for other reasons and these were also included in the analysis. We identified breeding partners of target individuals and if a blood sample was available, we measured the telomere length from that sample. Most such partners were captured as adults, so their exact age was unknown, but we knew relative age from the date of capture. This sampling regime resulted in 1,591 individuals with telomere lengths of which 1,500 had been sampled at 10 days of age. Telomere measures: Some of the telomere measures used in this analysis came from the dataset used in Heidinger et al. (2021). Subjects in that analysis had all been sampled at 10 days of age and resampled at least once, but many did not breed. The present analysis included breeders only and also added individuals sampled at 10 days of age that were not recaptured and resampled. The birds used in the Heidinger et al. (2021) dataset were randomly assigned to assays that were carried out by A. Kucera whereas the remaining samples were randomly assigned to assays and processed by R. Young. All assays used the same reference sample (described below). Telomeres were measured in whole blood samples, a highly replicative tissue that can be non-destructively sampled and is well suited for telomere analyses, especially in birds which have nucleated red blood cells. We extracted DNA using Qiagen DNeasy Blood and Tissue DNA Extraction Spin Column Kits and following the manufacturer' instructions. DNA concentration was quantified using a Nanodrop 8000 Spectrophotometer (Thermo Scientific, Waltham, MA, USA). DNA quality was assured using electrophoresis on a 2% agarose gel.We used real-time quantitative polymerase chain reaction (qPCR) on an Mx3000P (Stratagene, San Diego, CA, USA) to measure relative telomere length. The methods followed those of Criscuolo et al. (2009) adapted for house sparrows (Young et al. 2022). To measure relative telomere lengths (T/S ratios) for each sample, we ran separate qPCR reactions for telomeres and the single copy control gene, Glyceraldehyde 3-phosphate dehydrogenase (GAPDH), on different plates. The suitability of GAPDH as a control gene in this species has previously been verified using a melt curve analysis (Young et al. 2022). All the samples were run in duplicate and randomly distributed across plates.Each telomere and GAPDH reaction contained a total volume of 25 µl comprised of 20 ng of DNA and either telomere ((forward tel1b (5'-CGGTTTGTTTGGGTTTGGGTTTGGGTTTGGGTTTGGGTT-3') and reverse tel2b (5'- GGCTTGCCTTACCCTTACCCTTACCCTTACCCTTACCCT-3')) or GAPDH primers (GAPDH - forward (5'-AACCAGCCAAGTACGATGACAT-3') and reverse GAPDH (5'-CCATCAGCAGCAGCCTTCA-3')) at 200 nM concentration mixed with 12.5 ul of perfecta SYBRE green supermix Low Rox (Quantabio). The thermal profiles were 10 min at 95 °C, followed by 27 cycles of 15 s at 95 °C, 30 s at 58 °C, and 30 s at 72 °C, finishing with 1 min at 95 °C, 30 s at 58 °C, and 30 s at 95 °C for the telomere reactions and 10 min at 95 °C, followed by 40 cycles of 30 s at 95 °C and 30 s at 60 °C, finishing with 1 min at 95 °C, 30 s at 55 °C, and 30 s at 95 °C.We calculated relative telomere lengths (T/S ratio) according to the following formula: 2ΔΔCt, where ΔΔCt = (Ct telomere - Ct GAPDH) reference sample – (Ct telomere - Ct GAPDH) focal sample (Stratagene 2007), where the Ct value is the number of PCR cycles necessary to accumulate a sufficient fluorescent signal to cross a threshold. Individuals with longer telomeres cross this threshold more quickly than individuals with shorter telomeres, relative to the same house sparrow reference sample used on all plates. This same reference sample was also used to create a 5-point standard curve (40, 20, 10, 5, 2.5 ng) to ensure that all samples fell within the bounds of the standard curve and to calculate reaction efficiencies. In total we ran 53 plates and in all cases the efficiencies were between the recommended 85-115 % (mean ± SEM: telomere: 96.5 ± 0.85% and GAPDH 100 ± 0.99%). This assay produces highly repeatable T/S ratios (ICC: 0.86-0.88) when samples are run in random well locations across plates (Heidinger et al. 2021, Young et al. 2022).Parentage confirmation and pedigree construction: We confirmed maternity and paternity using genotypes obtained from a multiplex PCR of 5 highly polymorphic microsatellite loci adapted from well-established protocols (Stewart et al. 2006, Dawson et al. 2012). Generally, parents and offspring were organized by family, amplified on the same plate, and products analyzed in the same round on an ABI 3730. If adults had multiple broods across several years or switched mates, they may have been run on a different plate than some of their offspring. In some cases, subjects were analyzed 2-3 times on different plates and in all cases, their allele scores were within 2 base pairs of each other. Adults were excluded as parents if 2 or more loci diverged by more than 1 base pair. In some cases, adults were excluded from being the parent of a whole brood. This sometimes led to the identification of alternative adults that did match, and if these were not among our target individuals, that brood was excluded from further analysis. Only offspring that were genetically related to a target adult were included in the dataset. If the focal adult was a male and the offspring was extra-pair, it was excluded from our analysis. If the focal adult was a female, then any of her offspring sired by an extra-pair male were included in the analysis, but we did not identify the true sire in these cases. The parentage information allowed us to create a pedigree of individuals with known telomere lengths and their relatedness to other individuals.The resulting pedigree of 1,629 individuals had 1,591 telomere measures (both nestlings and adults), 1,500 of which were from individuals sampled at 10 days of age. This included 211 sires and 228 mothers, no cases of inbreeding, and 7 tiers (Table S1, Figure S1) with one lineage having 6 generations. The data set covered 428 nesting attempts with at least one nestling sampled. Most mothers had only one nesting attempt (N = 152) but two were represented by offspring from 10 attempts. Mothers were typically part of only one pair (N = 215), meaning that they only had one male partner, but 44 had 2 or more partners with one having 6. Since extra-pair offspring were not assigned to a sire, these numbers reflect social pairings that produced within-pair offspring. We determined the sex of all nestlings using PCR of extracted DNA, following procedures described in Westneat et al. (2002).Statistical Analysis Goal 1: Variance partitioning. Our initial goal was to partition the variance in telomere lengths among sources given the structure of the data. The initial animal model included all adults with a telomere measure, regardless of when they were sampled. We repeated the analysis restricting the set of adult telomere measures to those that were sampled at 10 days of age. These models were fitted with a set of random effects and no fixed effects. We included the relatedness matrix as a key term for assessing covariances among individuals by relatedness, the assay identity as a measure of lab artifacts, year (cohort) to assess effects of year-to-year environmental variation (note: prior analyses showed no effect of latency between sampling and analysis on telomere length, Heidinger et al. 2021), breeding attempt identity to assess all environmental effects common to a breeding attempt (e.g., date in the season, weather conditions during embryo and nestling development, size of the brood), nesting location (labeled as "barn" since nest boxes were clustered on the sides of farm buildings) to capture local environmental effects common to offspring reared in the same location, and pair identity to capture consistent joint influences of the male and female associated with breeding attempts. We did not include maternal identity alone, as males can provide extensive care during incubation and chick rearing, and most adult females in the dataset were members of only one pair and had only one breeding attempt, so the combination of pairD, motherID, and attemptID (or attemptID, motherID and sireID) would have over-specified the model.We reanalyzed the initial model to include only data from individuals sampled at 10 days of age. Variance partitioning models were fitted twice in the R computing environment (R Core Team 2019), first using brms (Burkner 2017) with 2 chains, default priors, 10,000 iterations with a burnin of 1000 and a sampling every 5 iterations. We used this to assess how well the model performed and get initial insight into the magnitude of variance components. The models took a long time to run in brms, so once we confirmed that models were fitting well (no divergent chains and rhats less than 1.1), we reran models in the program MCMCglmm (Hadfield 2010) using weakly informative priors (V = 1, nu = 0.002, Lemoine 2019), 1.1 x 105 iterations, a burnin of 10,000 and sampling every 50 iterations. We then analyzed the data in more detail to assess three possible influences on genetic variance in telomeres. First, we asked if parent age at the time of fertilization influenced telomere length (sex-specific parent age entered as two fixed effects). Second, we asked if the genetic variance in telomere length was a function of parent age via separate analysis of sire-specific (1-9) or mother-specific (1-6) age as a linear random slope in the "animal' random effect (Class et al. 2019. In these models we dropped any random effects that could not be distinguished from 0 and replaced pair identity with either mother or sire identity depending on which sex was in the random slope term. We allowed the model to estimate covariances between the animal intercept and slope. These models were fitted in brms with a linear slope, cauchy priors (0,2), 10000 iterations with a burnin of 1000, and a sample every 5 iterations. Finally, we assessed the role of offspring sex, both as a fixed effect to repeat a test for any sex differences in telomere length (Heidinger et al. 2021, Le Pepke 2021) and to assess sex-specific heritability. The former was tested by adding a fixed effect of offspring sex to the models described below that explored specific environmental factors. For sex-specific heritability, we fitted a bivariate model with the telomere length in male and female nestlings as two responses to investigate differences in heritability (Olsson et al. 2011, Chick et al. 2022) and estimate genetic correlations between the sexes. We also included assay identity, year, and attempt identity as random effects and modelled the covariance between male and female for each of these as well. Models were fitted in brms with each of three prior types for variance/covariances: default, V = diag(2), nu = 1.002, and V = diag(2)*(0.002/1.002), nu = 1.002. Results were affected only slightly. All models had 10,000 iterations, burnin of 1000, and a sampling of every 5. Goal 2: Influence of environmental effects. We examined types of environmental effects in more detail using a sequenced approach. Step one was to assess the magnitude of variance explained by the "environmental" variance components in the initial animal model described above. Our plan was to then explore specific environmental variables in more detail by adding fixed effects that varied most at the requisite level. For example, we expected from Le Pepke et al (2022b) that cohort (year) might explain some variation. If so, then yearly differences in weather might be relevant, so we gathered year-specific mean temperature and precipitation in either the "spring" months (February and March) preceding each breeding season or "summer" months through the period of breeding (April – August). Similarly, any among-attempt variation might be influenced by date in the season in which attempts were started, brood size at hatch, or the specific mean precipitation that occurred during the 25 days after the first egg was laid in that breeding attempt. Because temperature for a specific nesting attempt is highly correlated with date, and previous analyses of both variables revealed that date was a better variable to include for some traits (Westneat et al. 2009), we left attempt-specific temperature out of the model. Because these models focused on environmental sources of variance, we omitted the pedigree and included year and attempt identity as random effects. Each of these investigations was modelled in brms using default priors, 10,000 iterations, a burnin of 1000, and a sample every 5 iterations. We also explored the potential for interactive effects among environmental factors in two ways. To gather general evidence of interactions, we paired down the mixed model to just the informative and relevant random effects (omitting PairID, Barn, and the pedigree). We then split either the breeding attempt random term or the residual via grouping by a higher order random effect. For example, we calculated among-breeding attempt variance and among-nestling-within-attempt (residual) variance in telomere for each year. We also assessed if among- or within-attempt variance differed among parental ages by converting parental age to a random effect (lumping older ages into a 5+ category) and splitting among attempt and residual variance by age, with each parent tested in separate models. These models were fitted in a frequentist framework in SAS 9.4 (SAS 2015) Proc Mixed given the ease with which SAS coding allows this (see supplementary material). The magnitude of improved fit was assessed using a likelihood ratio test against the base model without group-specific variances. We also tested a set of two-way interactions as fixed effects that we reasoned were likely to explain variance at the appropriate level. For instance, if variation among attempts was important, then date by brood size seemed a likely influence given prior work showing that date affects clutch size (Westneat et al. 2009) and nestling growth (Mock et al. 2009) and brood size influences telomeres in a North Dakota population (Young et al. 2022). Similarly, at the residual level, brood size or parent age by the sex of the nestling might be important. We constructed these models after running the initial random effect models. Because these proposed fixed effects differed dramatically in measurement scale, we standardized ((x-mean)/SD) all values of all variables and analyzed them in a mixed model with a reduced set of random effects (omitting Barn, PairID, and the pedigree) in brms with default priors, 10,000 iterations, a burnin of 1000, and sampling every 5 iterations. Goal 3: Phenotypic and genetic covariance between telomeres and body size. Two opposing results already published encouraged us to investigate links between individual offspring telomere length and their condition at the time of measurement. Young et al. (2022) found that relative offspring size positively predicted day 10 telomere length but Le Pepke et al. (2022a) found that telomere length was negatively correlated with offspring size with age controlled, suggesting that higher growth led to shorter telomeres. We used two bivariate analyses to explore these relationships in our data set. One of the bivariate equations had telomere length as the response, included assay identity as a random effect and the important random effects from the analyses described above (year and attempt identity) along with the relatedness matrix. The other equation had either nestling tarsus length, a standard measure of size, or nestling mass as the response. These equations both included nestling age as a fixed effect since there was some variation in the age at which nestlings were measured, and year, attempt identity, and the relatedness matrix as random effects. We set the models to extract covariances for each random effect including the relatedness matrix to assess environmental and genetic correlations and the residual covariance. These models were fitted using brms using default priors, 10,000 iterations, a burnin of 1000 and sampling every 5. As a check, we reran these models with Cholesky priors (lkj(2) as recommended for correlations among random effects in the set prior notes for brms. We found no substantive differences in results