Heritability of cortisol production and metabolism throughout adolescence: a twin study

Context: Inter-individual differences in cortisol production and metabolism emerge with age, and may be explained by genetic factors. Objective: To estimate the relative contributions of genetic and environmental factors to inter-individual differences in cortisol production and metabolism throughout adolescence. Design: Prospective follow-up study of twins. Setting: Nationwide register. Participants: 218 mono-and dizygotic twins (N=109 pairs) born between 1995-1996, recruited from the Netherlands Twin Register. Cortisol metabolites were determined in 213, 169 and 160 urine samples at the ages of 9, 12 and 17y, respectively.

A c c e p t e d M a n u s c r i p t 4 cortisol metabolism, distinct patterns of genetic and environmental influences were observed, with heritability that either increased with age or peaked at age 12y.

Précis
In a longitudinal twin study, we found that the contribution genetic factors to cortisol production decreased over adolescence, while it increased or peaked at age 12y for cortisol metabolism.

Introduction
Cortisol, the main product of the hypothalamus-pituitary-adrenal (HPA) axis, is a crucial steroid hormone in the physiological stress response following homeostasis disturbance (1).Dysregulation of HPA axis activity has been associated with cardiovascular diseases and psychiatric conditions, e.g., major depressive disorder, post-traumatic stress disorder, panic disorder, and chronic anxiety (2)(3)(4).It has been recognized that experiences in early life may induce permanent alterations in the settings of several endocrine systems, including the HPA axis (5,6).
There is a paucity of data on the magnitude of the contribution of genetic factors to variance in HPA axis activity.Several studies in monozygotic (MZ) and dizygotic (DZ) twins addressed the relative contributions of genetic and environmental factors on serum and salivary cortisol levels.However, the magnitude of the heritability estimates varies (Table 1).A meta-analysis of five twin studies published before 2001, of which 4 were conducted in adults and 1 in both children and adults, estimated the heritability of basal cortisol assessed in serum or saliva at 62% (7).The total sample size in the meta-analysis was small (209 MZ and 190 DZ pairs), and power analyses included in the paper indicated that the statistical power to distinguish between genetic and shared environmental influences was low.Later studies reported lower heritability estimates for salivary cortisol levels, and showed that those estimates differed at different time points of the circadian cortisol rhythm (8)(9)(10)(11).
In adults, the heritability of cortisol level was 32-34% in samples obtained directly after awakening or 30 minutes post-awakening (9).In the evening samples obtained from these adults, the heritability was found to equal zero.In children, the heritability estimates were 28%, 60%, and 8%, directly after awakening, 30 minutes post-awakening, and in the evenings, respectively (8).In summary, previous twin research consistently found that variation in morning salivary, serum and hair cortisol levels appears to be at least partially heritable (Table 1).Cortisol secretion varies with age, and interindividual differences in diurnal cortisol levels are presumed to emerge during the second decade of life (12).Earlier studies have focused on cortisol levels in serum or salivary.These cortisol levels represent the net effect of cortisol production and metabolism.The contributions of genetic and environmental factors to individual differences in cortisol production or metabolism, as determined by cortisol metabolite excretion in urine, have yet to be studied.Day-to-day excretion of cortisol in the urine is moderately stable (13).Cortisol is metabolized to cortisone by 11β-hydroxysteroid dehydrogenase (HSD) type 2 in the kidney, while the reverse reaction occurs by 11β-HSD type 1 in liver and adipose tissue.Cortisol is also metabolized irreversibly by the A-ring reductases (5α-reductase and 5βreductase), and cytochrome P450 (CYP) 3A4 in the liver.Cortisol metabolism is stable across the menstrual cycle (14).The aim of the current study was to focus on indices of cortisol production and metabolism across adolescence, and to estimate the relative contributions of genetic and environmental factors to cortisol production and metabolism in a sample of children who were registered at birth in the Netherlands Twin Register (NTR).The twins took part in a longitudinal study and were seen at 9, 12 and 17 years of age (15).

Participants
We conducted a prospective follow-up study in mono-and di-zygotic twin pairs.Participants in this study were recruited from the Netherlands Twin Register (NTR) (16,17), and invited to take part in the BrainScale study of cognition, hormones and brain development (15,18).Parents of twins born between 1995 and 1996 were invited by letter four to eight weeks before the ninth birthday of the twins.Two weeks later, the parents were approached by phone, to explain about the study and to ask them if they consented to their children taking part.Of the 214 families who were approached, 109 consented to take part (51%).Seventy-eight and seventy-three percent of the participants took part in the follow-up study, and provided samples at the ages of 12 and 17 years, respectively (Figure 1).(19).Twins were categorized based on zygosity and sex: MZM (monozygotic males), DZM (dizygotic males), MZF (monozygotic females), DZF (dizygotic females), DZOS (dizygotic opposite sex male-female or female-male).To establish zygosity, participants were requested to collect buccal swabs from which DNA was isolated.All DNA samples were tested for genome-wide single nucleotide polymorphic (SNP) markers (20).
BrainScale is a collaborative project between the Netherland Twin Register at the Vrije Universiteit Amsterdam and University Medical Center Utrecht.The project was approved by the Central Committee on Research Involving Human Subjects of The Netherlands (CCMO), and studies were performed in accordance with the Declaration of Helsinki.Parents signed informed consent forms for the children and for themselves.Children signed their own informed consent forms at the third measurement.Parents were compensated for travel expenses, and children received a present or gift voucher at the end of the testing day.In addition, a summary of cognition scores and a printed image of their T1 brain MRI scan, when available, were provided afterwards The current study was approved by the medical ethics committee of the Amsterdam UMC, location VUmc.

Study protocol
At the ages of 9, 12 and 17 years, participants visited the study site for different tests.In the week prior to the study visit, urine samples were collected upon awakening in specially provided tubes.
Participants or their parents were requested to store the tubes in their refrigerator, and to bring them to the study visit.Samples were subsequently stored at -20 and -80 degrees Celsius, and thawed only once just before analysis.Urine samples were available from 47 monozygotic (23 male, 24 female) and 62 dizygotic twin pairs (22 male, 21 female, 19 opposite sex).

Laboratory analysis
Analysis of cortisol metabolites was conducted at the Edinburgh Clinical Research Facility Mass Spectrometry Core Laboratory.Glucocorticoid metabolites were measured by gas chromatographytandem mass spectrometry (GC-MS/MS) (21).Samples were analyzed in fifteen batches.Ratios of cortisol metabolites representing the activities of various enzymes involved in cortisol metabolism were calculated, as depicted in Table 2.

Statistical analysis
Outliers, defined as any value greater than 3 standard deviations above the phenotypic mean, were excluded from statistical analysis (on average six per index).Twin pairs with highly discordant outcomes, as assessed by visual inspection of scatterplots, were also removed (on average 0.76 per index), given their impact on the MZ and DZ correlations.We corrected for batch effects by fitting a random effects model to the twin data, which included batch as a random effect (22).As the sampling unit is twin pairs, we included in the model family as a random effect to estimate the MZ and DZ twin (intraclass) correlations.These analyses were done in R 3.4.2.using the nlme library (Nonlinear Mixed-Effects Models) (23,24).The batch-corrected phenotypic data were subject to subsequent analyses, as described below.
The classical twin design exploits the fact that MZ twins share 100% of their alleles identically by descent (IBD; from the biological parents), and DZ twins on average share 50% of their alleles IBD (25).This difference in genetic relatedness allows us to estimate the relative contributions of genetic We used OpenMx to fit the ACE or ADE model to the batch corrected twin data (26).The genetic and environmental variance components were estimated by maximum likelihood.We included as fixed covariates in the model sex and gestational age, since these factors have previously been associated with HPA axis activity (12,27).Below, we present the estimates of total standardized variance components.The standardized variance components were obtained by dividing by the total variance.The twin correlations and variances for each index are presented in Table 3.The DZ correlations (average .257)were lower than the MZ correlations (average .415) in 17 of the 21 phenotypes, pointing towards a relatively simple additive (AE) genetic model.Thirteen of the phenotypes appeared to be consistent with an ACE model.We note that some correlations were zero (MZ: 5βreductase activity at age 12; DZ: Renal 11β-HSD type 2 activity).We attribute this to sampling fluctuation, given the relatively small sample sizes.

A
The relative contributions of genetic (additive and non-additive), shared environmental, and unshared environmental factors, along with 95% confidence intervals, are displayed in Table 4.At age 9y, the contribution of variation in the cortisol production rate was explained for approximately forty percent by genetic factors and sixty percent by unshared environmental factors.The contribution of genetic factors to the inter-individual differences in cortisol production rate decreased with age.At age 17y, genetic factors were no longer contributing to the cortisol production rate.Variation at this age is mainly (for 79%) explained by unshared environmental factors.For indexes representing A-ring reductases activity, the contribution of genetic factors to inter-individual differences was found to increase with age.The contribution of genetic factors to the inter-individual differences in activities of renal 11β-HSD type 2 activity, cytochrome P450 3A4 and 11β-HSD first increased from 9 to 12 years of age, and then decreased from 12 to 17 years.

Discussion
In this longitudinal study of twin pairs followed between the ages of 9 and 17 years, we demonstrated the relative contributions of genetic and environmental factors to the indices of cortisol production and metabolism throughout adolescence.The most important finding from our study is that the environment plays a key role in the production of cortisol, evidenced by the A previous meta-analysis of twin studies estimated heritability of basal cortisol, representing the net effect of cortisol production and elimination, at 62% (7).Despite this observation, others demonstrated no significant SNP heritability for plasma or salivary morning cortisol (29).The interpretation of such discrepant findings is complicated, when considering that approximately fifty percent of the variance in salivary cortisol was dependent on day-to-day fluctuations (30).
In our study, cortisol production was mainly determined by unshared environmental factors already at age 9y, and the contribution of unshared environmental factors was found to increase with age.This lends support to previous observations suggesting that the settings of the HPA axis are mainly determined by individual circumstances (31).There is overwhelming evidence from animal experiments and epidemiological studies demonstrating that experiences in early life may program future HPA axis activity, with data linking poorer quality of parental care to increased HPA axis activity along with increases in mental illnesses and cardiovascular diseases (32)(33)(34).These observations could be attributed to increased DNA methylations status and a reduced expression of the glucocorticoid receptor (GR) promotor in hippocampal regions (6,35,36).However, more recent evidence suggests that the time window in which epigenetic programming could occur, may extend into adulthood (37)(38)(39).In male middle-aged twins, saliva cortisol levels showed significant cortisol heritability estimates for laboratory measures, but not for measures in the home situation, suggesting that genetic factors influence cortisol responses to specific environmental stressors (40).
In animal experiments, the influences of early life experiences can be observed by randomly allocating offspring into groups with different amounts of exposure to the factor of interest, whereas in humans determining the long-term impact of life experiences on HPA axis activity is more challenging.Previous research in humans linked early-life experiences like battering, neglect,  (33,(41)(42)(43)(44).The findings obtained from these studies are inevitably confounded by factors associated with both exposure and outcome, such as low socio-economic class and low household income.Twin studies offer a powerful tool to study environmental contributions, by controlling for the family background.
In our study, cortisol metabolism was, in contrast to cortisol production, considerably influenced by genetic constitution, and our findings suggest distinct patterns of genetic and environmental contribution to the different metabolic pathways.The activities of the A-ring reductases are, unlike cortisol production, less influenced by unshared environmental factors, especially later in life.
Variation in A-ring reductases activity was for 0-23%, 0-51% and 51-66% explained by (non-additive and additive) genetic factors at 9, 12 and 17 years of age, respectively, indicating a predominant role of genetic constitution in the regulation of A-ring reductases with age.The influences of genetic constitution on the activities of 11β-HSD isozymes and CYP3A4 increased from 9 to 12 years, and then decreased from 12 to 17 years of age.An explanation for the varying influences found for these enzymes may lie in a complex interplay between cortisol metabolizing enzymes and the hormonal regulators of puberty and the pubertal growth spurt, in particular growth hormone (GH) and insulinlike growth factor 1 (IGF-1).A previous twin study has shown that the pubertal growth spurt, peaking between 12 and 14 years, is strictly genetically regulated (45).In addition, there is strong evidence suggesting that these hormones affect the clearance of cortisol (14,46,47).
Our study has several strengths and limitations.The major strength of our study was the long-term follow-up.Furthermore, participants were recruited from a nationwide twin registry and the numbers lost to follow-up were acceptably low.Consequently, selection bias is unlikely to explain our results.Another strength of our study was the use of GC/MS/MS analysis, providing highly reliable measurements.Thus, the findings as presented are unlikely to be explained by measurement error.
However, our study also has its limitations.Participants were requested to collect early-morning A c c e p t e d M a n u s c r i p t 13 urine samples.Preferably, a 24h urine sample would have been analyzed, since cortisol is secreted in a circadian rhythm.Next, sampling started no earlier than at age 9y, so that there remains a lack of knowledge on the genetic and environmental etiology in early childhood.Moreover, differences between boys and girls were observed, necessitating further research.Additionally, due to limited sample size, resulting in large confidence intervals and some remarkable twin correlations (e.g.MZ correlation<DZ correlation for cortisol production rate at age 17y), our results have to be interpreted cautiously.Finally, the collected samples were not randomly distributed across the different analytical batches, which may have led to systematic error.More specifically, samples from twin pairs were allocated in the same batches, which might overestimate the contribution of shared environmental factors.Therefore, random correction for batch effect was carried out.

Conclusion
Our current findings, along with previous observations, emphasize the significant role of individual circumstances on the settings of the HPA axis.Notably, the contribution of unshared environmental factors on cortisol production was considerable and was found to increase with age, implicating a predominant role of individual circumstances with ageing.In contrast to cortisol production, cortisol metabolism was considerably influenced by genetic constitution, and heritability of A-ring reductases was found to increase with age, resulting in a peak of the genetic contribution at the age of 17 years.
For 11β-HSD isozymes and cytochrome P450 3A4, this peak was found at the age of 12 years.
Downloaded from https://academic.oup.com/jcem/advance-article-abstract/doi/10.1210/clinem/dgz016/5586817 by Edinburgh University user on 04 November 2019 Downloaded from https://academic.oup.com/jcem/advance-article-abstract/doi/10.1210/clinem/dgz016/5586817 by Edinburgh University user on 04 November 2019 A c c e p t e d M a n u s c r i p t 6 Downloaded from https://academic.oup.com/jcem/advance-article-abstract/doi/10.1210/clinem/dgz016/5586817 by Edinburgh University user on 04 November 2019 A c c e p t e d M a n u s c r i p t 8 Participants were instructed to bring any packages of recently used medication to the study visits, revealing that there was hardly any recent use of medication in our sample.
and environmental factors to phenotypic individual differences in terms of genetic and environmental variance components.In practice, either an ACE model or an ADE model is fitted to the twin data, where A stands for additive genetic effects, D for non-additive or dominance (genetic) Downloaded from https://academic.oup.com/jcem/advance-article-abstract/doi/10.1210/clinem/dgz016/5586817 by Edinburgh University user on 04 November 2019 A c c e p t e d M a n u s c r i p t 9 effects, C for shared (or common) environmental effects, and E for unshared environmental effects.MZ and DZ pairs differ in their genetic relatedness, but share the pre-and postnatal environment.The extent to which resemblance in both types of twins is not explained by differential genetic resemblance, is the basis for identification of common environmental factors.An ACE model is fitted to twin data if the MZ phenotypic correlation (rMZ) is smaller than two times the DZ correlation (rDZ; rMZ<2*rDZ).With 2 groups of relatives, i.e.MZ and DZ twins, it is not possible to estimate four variance components, i.e. a model including non-additive genetic factors in addition to common environment is not identified.A choice needs to be made to evaluate an ACE or ADE model.If (rMZ>2*rDZ) an ADE model is chosen.Unshared environmental factors incorporates aspects of the environment that are child specific and results in differences between twin pairs, i.e. the differences within MZ pairs are due to these factors.Unshared environmental factors also includes measurement error.Path diagrams of the ACE and ADE models are shown in Figure2.We fitted a series of univariate models.
Downloaded from https://academic.oup.com/jcem/advance-article-abstract/doi/10.1210/clinem/dgz016/5586817 by Edinburgh University user on 04 November 2019 A c c e p t e d M a n u s c r i p t 11 predominant and increasing contribution with age of unshared environmental factors.In addition, we found distinct patterns of genetic and environmental contribution to the different cortisolmetabolizing pathways.
Downloaded from https://academic.oup.com/jcem/advance-article-abstract/doi/10.1210/clinem/dgz016/5586817 by Edinburgh University user on 04 November 2019 A c c e p t e d M a n u s c r i p t 12 emotional maltreatment, perinatal malnutrition, low birthweight and prematurity with future HPA axis activity Downloaded from https://academic.oup.com/jcem/advance-article-abstract/doi/10.1210/clinem/dgz016/5586817 by Edinburgh University user on 04 November 2019 Downloaded from https://academic.oup.com/jcem/advance-article-abstract/doi/10.1210/clinem/dgz016/5586817 by Edinburgh University user on 04 November 2019 A c c e p t e d M a n u s c r i p t 14

Figure 1 Figure 2
Figure 1 This flowchart presents the enrollment of participants for this study Figure 2 Path diagram representing the ACE (top) and ADE (bottom) model.

Figure 2
Figure 2 Path diagram representing the ACE (top) and ADE (bottom) model.