Statistical analyses
To identify and partition sources of variation in autosomal methylation
levels, we used the R package glmmTMB (Brooks et al. 2017)
to run a generalized linear mixed model (GLMM) with the BFGS algorithm
as an optimizer, using the number of methylated and unmethylated Cs
(represented by the sequenced Cs and Ts) for each position as our
dependent variable, assuming a binomial error distribution (Lea et
al. 2017). As fixed effects we added sex (as a two-level class variable
using males as the reference category) and age (as a covariate). With
respect to the latter, each individual’s age was partitioned into an
‘average age’ and ‘delta age’ component following van de Pol & Wright
(2009). Average age was calculated as the average of all ages at which
we assessed a bird’s autosomal methylation level, while delta age was
calculated as the difference between the bird’s actual age and its
average age (i.e. delta age = age – average age). When adding both age
variables as covariates, average age represents the among-individual,
and delta age the within-individual age effect (van de Pol & Wright
2009). If the among- and within-individual age effects were to differ,
this would indicate that the effect of age among individuals cannot be
explained by changes within individuals, revealing age-specific
selective (dis)appearance of individuals with certain levels of
methylation (van de Pol & Wright 2009). We additionally included the
interaction between sex and delta age in our model to test for
sex-differences in the within-individual age trajectory of methylation.
Random effects included were ‘bird identity’, genomic position
(‘position identity’) and an observation-level random effect
(‘observation’). The latter was added to account for overdispersion,
which was detected using the R package performance (Lüdeckeet al. 2021), whereas the first two were included to account for
pseudoreplication caused by repeated sampling of individuals and genomic
positions, respectively. Model evaluation and summary of parameter
estimates and statistics were conducted using the R packageparameters (Lüdecke et al. 2020). The marginal effects of
the interaction were plotted using the R-package sjplot test
(Lüdecke 2018).
Results