Recruitment and randomization

Families in psychosocially at-risk situations living in the German-speaking part of Switzerland (more specifically, in Canton Zurich) were identified through an interdisciplinary network of family counseling and health professionals (e.g., parent-counseling offices, pediatricians, midwives, social counseling, psychological and psychiatric services), and recruited before or shortly after the birth of the child. Recruitment took place at three project sites in the Canton of Zurich. Families were included in the study if they presented at least two of the following risk factors: (a) parental risks (low level of education, early parenthood, alcohol or drug abuse, sickness and disabilities), (b) familial risks (single parenthood, partnership conflicts), (c) social risks (lack of social integration, dissocial environment), (d) material risks (unemployment, financial problems, confined living space), and (e) child-related risks (high-risk pregnancy, health issues, prematurity). Exclusion criteria were: (a) no permanent residency permit, (b) severe illness or disability of the child, (c) severe illness or disability of the parent requiring inpatient and long-term psychiatric treatment, and (d) other intensive treatments or child protection procedures. The screening was performed using a self-developed questionnaire (Neuhauser et al., Reference Neuhauser, Ramseier, Schaub, Burkhardt, Templer and Lanfranchi2015).

After prospective participants signed informed written consent, approximately one half of the recruited families were assigned to the intervention group (IG), and the other half to the control group (CG). Families signed the informed written consent before knowing to which group they would be assigned. The assignment was performed using stratified block randomization. Strata were project site, cumulative psychosocial risk factors assessed using the Heidelberger Belastungsskala (HBS) (Sidor, Eickhorst, Stasch, & Cierpka, Reference Sidor, Eickhorst, Stasch and Cierpka2012), marital status (single parent: yes/no), and German-language proficiency (interpreter: yes/no) (Schaub et al., Reference Schaub, Ramseier, Neuhauser, Burkhardt and Lanfranchi2019). An additional consent form was signed by parents who agreed to provide saliva of their child in addition to behavioral data. Approximately 50% of the families agreed to provide saliva samples. Mothers in participating families were easier to reach (e.g., by phone calls, emails) (p < .01), and slightly older (M = 30.4 years, standard deviation [SD] = 6) than mothers in families that did not agree to participate in biological analyses (M = 28.9 years, SD = 5.2). However, the two groups did not differ on other relevant study variables, such as SES or maternal depressive symptoms. The resulting IG and CG subgroups did not differ regarding the strata mentioned above. The flowchart of the participants’ progress through the phases of the study is presented in Figure 1. Families assigned to the CG had access to community services, but were not supported by PAT. Incentives such as birthday greeting cards, small birthday presents, and a monetary benefit of 70 CHF for data collection were offered to both experimental groups.

Figure 1. CONSORT flow diagram modified from Schaub et al. (Reference Schaub, Ramseier, Neuhauser, Burkhardt and Lanfranchi2019). ¹Criteria for participation were no longer met after randomization and increasing age of the child. ²Insufficient sequencing coverage.

Participants

The final sample for the present study included N = 132 children (CG; n = 60, IG; n = 72), including seven pairs of twins (CG; n = 3, IG; n = 4). Children were aged 3 years (M = 35.5 months, SD = 1.09), and approximately half were female (n = 72, 54,5%). Families were mostly of European origin (71.2%,), and the remaining were originally from Africa, Asia, and Latin America. Notably, the primary language of 79.5% of the children was non-German. The sociodemographic characteristics of the overall sample and of the individual experimental groups are summarized in Table 1.

Table 1. Demographic and biopsychosocial characteristics of the participating children

Note: N = 132, IG, n = 72; CG, n = 60. Methylation is presented in % of methylated cytosines in the assessed CpGs (NGFI-A binding regions of the NR3C1 exon 1F). IG was coded 1, CG was coded 0.

Abbreviations: ADHD = attention deficit/hyperactivity disorder; CG = control group; CpG = cytosine–guanine dinucleotide; IG = intervention group; M = mean; SES = socioeconomic status; SD = standard deviation

* p < .05.

Parents as Teachers

The aim of PAT, which is to enhance children's development, was pursued using four elements. Firstly, home visits during the first 3 years of life conducted by a qualified parent educator, with a minimum of ten visits per year. In this study, members of the IG received, on average, a total of n = 49.9 (SD = 10.1) home visits before data were collected. During the visit, the parent educator worked with the family on three areas: development-oriented parenting, parent–child interactions, and the well-being of the family. The second element comprised nonmandatory group meetings once a month to allow parents to network with other parents and receive information related to educational practices, Parent×Child interactions, and community services for families. The third element was monitoring of the child's general health and development and the final element was support of parents in their integration in the community and the referral of specific needs to other public institutions.

Data collection

The data collection was conducted at the family home, around the child's third birthday. It consisted of self-report questionnaires, mostly filled out by the mothers, and saliva sampling carried out at the end of the collection session. For mothers that were not proficient in German, questionnaires were available in the mother's first language. Descriptive statistics of the study variables, in the overall sample and based on experimental groups, are reported in Table 1.

Maternal depressive symptom

The Brief Symptom Inventory (BSI-18; Spitzer et al., Reference Spitzer, Hammer, Löwe, Grabe, Barnow, Rose and Franke2011) was used to assess self-reported maternal depressive symptoms over the past week. Although the timeframe of the assessment was short, previous studies have demonstrated that differential maternal mood and behavior may induce rapid changes in NR3C1 methylation during early life in the rat and human offspring, even over a few weeks (Murgatroyd et al., Reference Murgatroyd, Quinn, Sharp, Pickles and Hill2015; Weaver et al., Reference Weaver, Cervoni, Champagne, D'Alessio, Sharma, Seckl and Meaney2004). The BSI-18 was developed as a screening tool for assessing psychological distress (somatization, depression, and anxiety) in clinical and nonclinical populations. The depression subscale includes six items rated on a 5-point Likert-type scale (0 = not at all; 1 = a little bit; 2 = moderately; 3 = quite a bit; 4 = extremely) and assesses thoughts of suicide, feelings of loneliness, melancholy, anhedonia, feeling of hopelessness, and feelings of worthlessness. The internal consistency in the present study was excellent (Cronbach's α = 0.92). The mean score of the six items was calculated to obtain an approximately continuous variable. Scores of the BSI-18 assessed a year earlier were also employed in the present study (Cronbach's α = 0.91).

Parental disagreement

The following two items of the Parental Stress Questionnaire (PSQ; Domsch & Lohaus, Reference Domsch and Lohaus2010) related to parenting attitudes in partnership were used as a proxy of parental stress (termed as “parental disagreement” in this study): “My partner and I totally agree on questions about upbringing” and “My partner and I discuss and decide together about parenting tactics”. The rest of the PSQ items were not assessed in the context of ZEPPELIN. For the two assessed items, the internal consistency in the present study was very good (Cronbach's α = 0.81). Items are rated on a 4-point Likert scale (1 = strongly agree, 2 = agree, 3 = disagree and 4 = strongly disagree) and were answered by both mothers who lived with their partner and mothers living alone with their child (single mothers), mostly of whom were still in contact with the father of the child. The mean score of the two items was calculated to create an approximately continuous variable. The scale for the variable scores was reversed to facilitate the interpretation of the results; high scores are indicative of more risk.

SES

To assess the families’ SES we used the well-established International Socio-Economic Index of Occupational Status (ISEI) (Ganzeboom & Treiman, Reference Ganzeboom and Treiman1996). The ISEI was constructed using data from almost 200.000 men and women in 42 countries. This index is a by-product of the mean educational and income levels for each occupation that is converted into ISEI. The ISEI scale ranges from 0 to 100, with greater scores indicating higher SES. Again, we reversed the variable scores in order to facilitate the interpretation of the results; higher scores correspond to higher risk.

Child behavior problems

Child behavior problems were assessed using the well-established Child Behavior Checklist (CBCL 1½–5; Achenbach & Rescorla, Reference Achenbach and Rescorla2000). Parents were asked to rate the child symptoms “now or within the past 2 months”. We used the five Diagnostic and Statistical Manual of Mental Disorders, fourth edition (DSM-IV)-oriented subscales which measure both children's internalizing and externalizing symptoms. Internalizing symptoms include affective problems (e.g., underactive, lack of energy), anxiety problems (e.g., not wanting to leave home), and pervasive developmental problems (e.g., disturbed by any change in routine). Externalizing symptoms include attention-deficit/hyperactivity problems (e.g., quickly shifting from one activity to another), and oppositional defiant problems (e.g., being aggressive). The five subscales include, collectively, 45 items that are rated on a 3-point scale (0 = not true, 1 = somewhat or sometimes true, 2 = very true or often true). Row scores were standardized in T-scores (M = 50 ± 10).

Saliva sampling

In this study, we used saliva as a source of DNA. Methylation studies that mainly target the human brain rely on postmortem tissues; this presents clear unpractical implications (Rhein et al., Reference Rhein, Hagemeier, Klintschar, Muschler, Bleich and Frieling2015). Therefore, easily accessible peripheral tissue such as blood, buccal cells, and saliva are widely used as a proxy of the brain for assessing DNA methylation in humans. In particular, the use of saliva carries many advantages, including undemanding and minimally invasive collection, as well as easy handling of samples (Yoshizawa et al., Reference Yoshizawa, Schafer, Schafer, Farrell, Paster and Wong2013). Because methylation is known to be cell- and tissue-specific, the choice of tissues can affect results and their interpretation (Fiori & Turecki, Reference Fiori and Turecki2016; Zhang et al., Reference Zhang, Zhou, Lin, Lowdon, Hong, Nagarajan and Wang2013). However, studies assessing methylation–disease associations have demonstrated that methylation in peripheral tissues can reflect results obtained using brain tissue (Turecki & Meaney, Reference Turecki and Meaney2016). For instance, studies assessing NR3C1 promoter methylation in saliva extensively demonstrate a link between ELS, increased methylation, and child or adult diseases (Cicchetti & Handley, Reference Cicchetti and Handley2017; Efstathopoulos et al., Reference Efstathopoulos, Andersson, Melas, Yang, Villaescusa, Rȕegg and Lavebratt2018; Murgatroyd et al., Reference Murgatroyd, Quinn, Sharp, Pickles and Hill2015; Oberlander et al., Reference Oberlander, Weinberg, Papsdorf, Grunau, Misri and Devlin2008; Parade et al., Reference Parade, Ridout, Seifer, Armstrong, Marsit, McWilliams and Tyrka2016; Tyrka et al., Reference Tyrka, Parade, Eslinger, Marsit, Lesseur, Armstrong and Seifer2015), reflecting results obtained by studies using brain tissues (McGowan et al., Reference McGowan, Sasaki, D'Alessio, Dymov, Labonté, Szyf and Meaney2009; Labonte et al., Reference Labonte, Yerko, Gross, Mechawar, Meaney, Szyf and Turecki2012).

Saliva from children in the present study was collected in Oragene tubes (OG-500, DNA Genotek, Ontario, Canada) using the passive-drool method according to the manufacturer's instructions, and stored at 4°C until subsequent biochemical analyses.

DNA methylation

DNA was isolated using the PrepIT-L2P protocol (OG-500, DNA Genotek, Ontario, Canada). Sodium bisulfite conversion was performed with 500 ng of DNA using the EZ DNA methylation kit (Zymo Research, Irvine, CA, USA). We performed amplification of the NR3C1 target sequence (Figure 2) using the following specific primers: frw 5′-TTG AAG TTT TTT TAG AGG G-3′ and rws 5′-AAT TTC TCC AAT TTC TTT TCT C-3′, which included universal primer sequences CS1/CS2 on the 5′ ends (Fluidigm, San Francisco, CA, USA). Cycling conditions were as follows: 95°C for 3 min, then 40× (98°C for 20 s, 60°C for 15 s, 72°C for 15 s), and final elongation step at 72°C for 45 s. Amplicons were then purified using e-gel size selection (Thermo Fisher Scientific, Waltham, MA, USA). Subsequently, indexing with unique single barcodes (Fluidigm, San Francisco, CA, USA) was performed through a second polymerase chain reaction (PCR) (95°C, 3 min, then 10× (98°C, 20 s; 60°C, 15 s; 72°C, 15 s), and final elongation at 72°C for 45 s). Indexed amplicons were pooled and submitted to a final purification to remove dimers and amplification artifacts. The pooled samples were then diluted to a concentration of 2 nM and sequenced on an Illumina Miseq using the v3 kit (Illumina, San Diego, CA, USA).

Figure 2. Upper part: Schematic diagram of noncoding alternative first exon in the NR3C1 promoter region. Lower part: Sequence of the NR3C1 exon 1F, chr5:142,783,586–142,783,903 located in the 5′ untranslated region of the NR3C1. Dashed boxes delineate NGFI-A binding regions described by McGowan et al. (Reference McGowan, Sasaki, D'Alessio, Dymov, Labonté, Szyf and Meaney2009). These CpG sites are the most widely reported in the literature, in both human and animal studies assessing the effect of ELS on NR3C1 methylation. CpG 1–5 correspond to the CpG sites 30–32 and 37, 38 reported by Palma-Gudiel et al. (Reference Palma-Gudiel, Córdova-Palomera, Leza and Fañanás2015). Underlined sequences correspond to the primer's position.

After the sequencing, the software Trimmomatic v0.35 (licensed under GPL V3 and available at http://www.usadellab.org/cms/index.php?page=trimmomatic) was employed to identify and remove low-quality products (Bolger, Lohse, & Usadel, Reference Bolger, Lohse and Usadel2014). The Bismark program (v0.19.0) was used to extract the counts of methylated (cytosines) and unmethylated (thymine) bases. After summing up methylated and unmethylated counts, we only kept samples showing a coverage (=total counts) of at least 100× (Chen et al., Reference Chen, Gross, Lutz, Vaillancourt, Maussion, Bramoulle and Ernst2017). Three samples did not reach this threshold and were therefore excluded. For the remaining samples, coverage ranged from 2,539× to 125,583×, with a mean of 57,956× ± 29,040. In this study, we analyzed five CpG sites located within the NGFI-A binding regions of exon 1F (Figure 2). The NGFI-A is a powerful regulator of the expression of the GR, and there is consistency for increased methylation of the CpG sites located in the NGFI-A binding regions among association studies focusing on NR3C1 methylation and ELS (Turecki & Meaney, Reference Turecki and Meaney2016). Table 1 presents the average methylation and SD for each of the five CpG sites. One participant showed extremely high percentages of methylation at CpG 5 (86.9%) (M = 0.45%, SD = 2.1). This value was not considered as an outlier, as sufficient coverage was reached at this CpG site (9,531×). Moreover, the other CpG sites for this sample, including the proximal CpG 4, had low methylation (0.13%), which suggests successful completion of the bisulfite reaction. Finally, other studies found high methylation at this CpG, such as 52,17% (Efstathopoulos et al., Reference Efstathopoulos, Andersson, Melas, Yang, Villaescusa, Rȕegg and Lavebratt2018). However, we conducted all statistical analyses with and without the participant showing high percentages of methylation at CpG 5. Results were largely unchanged after the exclusion of this participant and conclusions were not affected.

Statistical analyses

We first performed pairwise Pearson's correlation between all variables: ELS, behavior problems, and NR3C1 methylation (Table 2). Linear regression was used for assessing the predictive effect of ELS and PAT on NR3C1 methylation. Because NR3C1 methylation was strongly positively skewed, we performed log transformation in order to approach normality. Partial correlation was used for assessing the association between NR3C1 methylation and child behavior problems. As NR3C1 methylation was not normally distributed, correlation analyses were performed using a nonparametric (Spearman) correlation. In order to assess whether NR3C1 methylation contributed as a mediator of ELS on child behavior problems, path analysis was conducted using PROCESS macro model 4 (Hayes, Reference Hayes2017). Predictor and outcome variables were included in the path analysis if they were significantly associated with the mediator. Current maternal depressive symptoms were not used in the mediation analysis, as the outcome variables (child behavior problems) were essentially measured before the predictor variable, posing a problematic temporal ordering. Child age at the time of saliva sampling (=child age), child age at the beginning of the intervention (=child age at randomization), sex, and geographical origin were included as covariates in all statistical analyses. The “PAT experimental group” variable was also included as a covariate in all statistical analyses, to control for possible effects of the intervention on the results. All methylation analyses were performed using the mean methylation of the five CpG sites within the NGFI-A binding regions (mean NGFI-A methylation) and methylation levels at individual CpG sites (CpG 1 to CpG 5). The mean composite of CpG sites was used in the statistical analyses, despite the lack of intercorrelation among CpG sites (Table 2). Indeed, the five CpG sites are located in the binding regions of the transcription factor NGFI-A. Therefore, physical NGFI-A binding impedance and consequent alteration in gene expression may result from the average level of methylation at the five CpG sites, rather than from the level of methylation at one single CpG site of the NGFI-A binding region. The percentage of missing values across the 14 variables varied between 0 and 16.7%. In total 165 out of 1,683 data points (8.93%) were incomplete. Questions related to the CBCL questionnaire were not answered by 21 families, while the ELS variables showed nonsystemic patterns of missing data. We used multiple imputation to create and analyze 50 multiply imputed datasets to ensure that a reliable estimate of the fraction of missing information was obtained (Madley-Dowd, Hughes, Tilling, & Heron, Reference Madley-Dowd, Hughes, Tilling and Heron2019; Van Buuren, Reference Van Buuren2018). Incomplete variables were imputed under the missing at random (MAR) assumption, using the default settings of SPSS (Enders, Reference Enders2010). We compared the results using the pooled imputed datasets with the results of complete cases and found no difference regarding the findings. Analyses were conducted using SPSS (IBM statistic, version 24.0, Armonk, NY: IBM Corp.). Mediation was considered significant when the 95% asymmetric confidence interval (CI) generated by bootstrapping did not include the value 0. In regression and correlation analyses, a p value below 0.05 was considered statistically significant.

Table 2. Bivariate correlation between ELS variables, child behavior problems and NR3C1 methylation

Note: Pearson's correlation.

Abbreviations: ADHD = attention deficit/hyperactivity disorder; CpG = cytosine–guanine dinucleotide; Pervasive dev. problems = pervasive developmental problems; SES = socioeconomic status

* p < .05.

** p < .001.