BDNF rs6265 methylation and genotype interact on risk for schizophrenia

Abstract

Epigenetic mechanisms can mediate gene-environment interactions relevant for complex disorders. The BDNF gene is crucial for development and brain plasticity, is sensitive to environmental stressors, such as hypoxia, and harbors the functional SNP rs6265 (Val⁶⁶Met), which creates or abolishes a CpG dinucleotide for DNA methylation. We found that methylation at the BDNF rs6265 Val allele in peripheral blood of healthy subjects is associated with hypoxia-related early life events (hOCs) and intermediate phenotypes for schizophrenia in a distinctive manner, depending on rs6265 genotype: in ValVal individuals increased methylation is associated with exposure to hOCs and impaired working memory (WM) accuracy, while the opposite is true for ValMet subjects. Also, rs6265 methylation and hOCs interact in modulating WM-related prefrontal activity, another intermediate phenotype for schizophrenia, with an analogous opposite direction in the 2 genotypes. Consistently, rs6265 methylation has a different association with schizophrenia risk in ValVals and ValMets. The relationships of methylation with BDNF levels and of genotype with BHLHB2 binding likely contribute to these opposite effects of methylation. We conclude that BDNF rs6265 methylation interacts with genotype to bridge early environmental exposures to adult phenotypes, relevant for schizophrenia. The study of epigenetic changes in regions containing genetic variation relevant for human diseases may have beneficial implications for the understanding of how genes are actually translated into phenotypes.

Keywords:

• BDNF
• DNA methylation
• epigenetics
• hypoxia
• obstetric complications
• prefrontal cortex
• rs6265
• schizophrenia
• working memory

Introduction

Differentially methylated regions often contain single nucleotide polymorphisms (SNPs) critical for human diseases.¹ CpG methylation state is recognized as a major determinant of natural genetic variation,² and a strong genetic component underlies inter-individual variation in DNA-methylation profiles.³ While it has been argued that many SNPs contribute to gene-expression changes and phenotypes relevant for diseases via epigenetic mechanisms,^3,4 the possibility that DNA methylation changes may compensate and/or modulate the effect of genetic variation has been less studied. Such event may be important in reconciling variable penetrance of genetic variants associated with human diseases. Our research has previously approached this issue by analyzing variable methylation of CpGs associated with functional SNPs, which likely originate during evolution from spontaneous point mutation of 5-methyl-cytosine into thymine. Specifically, we previously reported that methylation of a CpG created by a widely studied functional SNP in the COMT gene was sensitive to environmental experience and associated with brain phenotypes only in individuals homozygous for the ancestral allele, but no association was found in heterozygotes.⁵

In the present work we focus our analysis on another popular gene in behavioral genetics, highly important in neurodevelopment and in system-level neural phenotypes and sensitive to environmental factors, the gene coding for the Brain Derived Neurotrophic Factor (BDNF). BDNF is a highly regulated protein, a crucial factor for the development of the feto-placental unit⁶ and of the brain,⁷ as well as for neural plasticity, energy metabolism, learning, episodic and working memory (WM) ^7–12 and cancer.¹³ The effects of BDNF on synaptic plasticity and neuron survival strongly suggest a role for this factor in schizophrenia, a neurodevelopmental disorder whose risk is heritable and characterized by physiological prefrontal cortex (PFC) dysfunction during WM,^14–16 as well as reduced prefrontal levels of BDNF.¹⁷ However, BDNF has not been associated with schizophrenia in recent large scale genome wide association studies.¹⁸ BDNF expression is sensitive to early-life environment^19–21 and specifically hypoxic-stressors,^22,23 which, in turn, are critical factors involved in the pathophysiology of schizophrenia.^24–27 The failure to find association between BDNF and schizophrenia in large-scale case control studies may conceivably reflect a complex interaction of genotype and environmental experience that alters methylation status.

The SNP rs6265 (G > A, Val⁶⁶Met), in the human Pro-BDNF sequence, is an example of a well-characterized functional SNP in this gene, which influences intracellular trafficking, activity-dependent secretion, as well as memory related behavior and brain activity in transgenic mice and in humans.^8–11 However, while cell and animal models are consistent in showing a functional effect, studies in humans reveal discrepancies in associating the Met or the Val allele to phenotypes relevant for schizophrenia and different psychiatric disorders.^{8–10,28–33} Interestingly, hypoxia-related events during pre-, peri- and early post-natal life (hOCs) have been reported to interact with BDNF in affecting risk for schizophrenia, independently from rs6265. Specifically, birth hypoxia has been associated with reduced BDNF levels in cord samples, taken at delivery, of individuals who developed schizophrenia later in life,³⁴ and other genetic variation of the BDNF signaling interact with early life hypoxia in affecting risk for the disorder.³⁵ Surprisingly, these studies failed to show a role of rs6265 in interacting with obstetric complications in affecting risk for schizophrenia. A possible explanation for these results is that epigenetic mechanisms related to these environmental factors modulate the effect of the most likely relevant functional SNP of BDNF. Indeed, the G>A (Val>Met) substitution of the rs6265 SNP creates/abolishes a CpG site, so that the Val allele has a CpG at this position while the Met allele does not. Since this CpG is differentially methylated in humans,³⁶ we hypothesize that the epigenetic changes able to modulate the effect of rs6265 may correspond to a different methylation status of the Val allele in ValVals and ValMets. In other words, while in experimental models the rs6265 genotype is clearly associated with certain outcomes, in living people who experience actual life and its many stressors DNA methylation specifically of one allele of this SNP may change in response to those factors, such as hOCs, and may modulate specific phenotypes, thus weakening the effect of the rs6265 genotype.

We hypothesize that BDNF rs6265 methylation is sensitive to environmental (hOCs) exposure both in ValVal and ValMet subjects and is also associated with phenotypes relevant for schizophrenia. To this purpose, we analyze (Fig. 1):

Figure 1. Research design.

1. the interaction between hOCs exposure, a known risk factor for schizophrenia, BDNF rs6265 genotype and methylation in blood;

2. the interaction between hOCs exposure, BDNF rs6265 genotype and methylation on intermediate phenotypes of schizophrenia, i.e., working memory performance and related neural activity measured with fMRI;

3. the heritability of rs6265 methylation and the interaction between BDNF rs6265 genotype and methylation on genetic risk for schizophrenia;

4. the interaction between BDNF rs6265 genotype and methylation on serum BDNF levels, which are known to be altered in patients with schizophrenia;

5. the relationship between brain and blood methylation in a sample of transgenic mice;

6. the effect of the rs6265 genetic variation on DNA-protein binding, potentially relevant in explaining the relationship between methylation and BDNF levels.

We found that hOCs exposure is associated with opposite methylation changes in ValVal and ValMet subjects. These methylation changes are associated with intermediate phenotypes of schizophrenia and being partially inherited predict genetic risk for the disorder, differently in ValVals and ValMets. Moreover, the rs6265 genotype affects binding of the transcription factor BHLHB2, potentially explaining the opposite relationship between methylation and serum BDNF levels that we detected in ValVals and ValMets. The correlation between brain and blood methylation in mice further supports the possibility of using peripheral methylation as a proxy of epigenetic changes relevant for the brain.

Results and discussion

Relationship between rs6265 genotype, methylation and hypoxia-related complications (hOCs)

Given that brain methylation cannot be studied in vivo in humans, and that epigenetic changes relevant for developmental programming can be part of a general response of the entire organism,^20,37 we used pyrosequencing to analyze rs6265 methylation in peripheral blood mononuclear cells (PBMCs) of 259 healthy humans. This site is differentially methylated in human PBMCs and, as predictable, methylation was affected by rs6265 genotype, given that only the Val allele has a CpG in this position (ANOVA: N=259; F_2,256=1601.2; P < 0.001***; ValVal > ValMet > MetMet; Fig. S1). ValVal subjects have greater methylation compared with ValMets, while MetMet subjects have no methylation. All analyses were therefore performed only in ValVal and ValMet individuals, and methylation values were normalized in order to test for interactions between genotype and methylation.

First, we analyzed in living healthy humans whether BDNF rs6265 methylation in PBMCs is associated with exposure to hypoxia-related events during prenatal, perinatal, and early postnatal life (hOCs), a known risk factor for schizophrenia. Using the McNeil-Sjöström Scale,³⁵ we assessed exposure to hOCs in 168 healthy humans demonstrating that the relationship between hOCs exposure and BDNF rs6265 methylation is affected by rs6265 genotype. In fact, we found an interaction between rs6265 genotype and hOCs on methylation (Factorial ANCOVA, with age and sex as covariates: N = 169: 110 ValVal, 59 ValMet; F_5,163 = 9.70; P = 0.002**; Fig. 2), so that rs6265 methylation is greater in the presence of hOCs in ValVal homozygotes (P = 0.03*), while it is attenuated in ValMet subjects exposed to hOCs (P < 0.02*; see Supplemental file 1 for detailed results). On the other hand, BDNF rs6265 methylation was not associated with age and sex, in the whole sample, and in the 2 genotype groups (P > 0.28). These results suggest that exposure to hOCs is associated with opposite and long-lasting changes of rs6265 methylation in ValVal and in ValMet subjects.

Figure 2. Relationship between methylation of BDNF rs6265 in PBMCs and hOCs exposure in healthy humans. Boxplot of rs6265 methylation (T scores) as a function of hOCs exposure: Val/Val homozygotes (N = 110) exposed to hOCs have greater methylation compared with ValVal not exposed, while Val/Met subjects (N = 59) exposed have reduced methylation compared with ValMet not exposed. See text and Supplemental file 1 for statistics.

Relationship between rs6265 genotype, methylation, hOCs, and intermediate phenotypes for schizophrenia

We then analyzed the relevance of these methylation changes for schizophrenia. Since we found that rs6265 methylation is related with hOCs which are a risk factor for schizophrenia,^25,35 we assessed the relationship between hOCs exposure, BDNF rs6265 genotype and methylation in PBMCs, and WM performance and related prefrontal activity, 2 well established intermediate phenotypes for this neurodevelopmental disorder.¹⁵ Consistent with the above results, we found that, in a sample of 212 healthy humans, BDNF rs6265 methylation is correlated with 2-back WM performance both in ValVal and ValMet subjects with a qualitatively different direction. More in details, a multiple regression of BDNF rs6265 genotype, methylation and hOCs and their interactions, with 2-back WM accuracy as the dependent variable (age and sex as covariates), revealed an interaction between rs6265 genotype and methylation on WM performance (N = 211: 145 ValVal, 66 ValMet; t = 2.63, P = 0.009**; Fig. 3), so that BDNF rs6265 methylation is associated with 2-back WM performance both in ValVal and ValMet subjects with a qualitatively different direction; consistent with previous studies,^8,38 this analysis also confirmed lower WM performance in ValMet compared with ValVal (t = −2.51, P = 0.01*). Post-hoc analyses of the interaction indicated that greater methylation is associated with lower accuracy in ValVal homozygous (t = −2.46, P = 0.01*), while the opposite is found in ValMet (t = 2.55, P = 0.01*; Fig. 3; see Supplemental file 1 for detailed results). Taken together, these results suggest that the relationship between methylation of the Val allele and WM performance is opposite between the 2 genotype groups. However, the opposite methylation changes associated in ValVals (increased methylation) and ValMets (blunted methylation) with hOCs exposure were similarly associated with decreased WM accuracy.

Figure 3. Relationship between methylation of BDNF rs6265 in PBMCs and working memory (WM) accuracy in healthy subjects. Scatterplot of the correlations between methylation of rs6265 (T scores) and WM accuracy: in Val/Val subjects (N = 145) increased methylation is associated with impaired accuracy, while in Val/Met heterozygotes (N = 66) blunted methylation is associated with impaired accuracy (See text and Supplemental file 1 for statistics).

We also analyzed the potential relationship between rs6265 genotype, methylation, hOCs, and prefrontal activity during WM measured with fMRI, another intermediate phenotype of schizophrenia, in a group of 141 healthy subjects who underwent fMRI during the N-back task (Fig. 4A-E). Multiple regression of the imaging data in SPM8 demonstrated: a positive correlation between BDNF rs6265 genotype and left dorsolateral PFC activity during 2-back, suggesting greater activity in ValMet subjects compared with ValVal (N = 141: 93 ValVal, 48 ValMet; x = −54, y = 24, z = 32, BA46, k = 31, Z = 3.42, p_{FWE-corrected} = 0.027*); a positive correlation between BDNF rs6265 methylation and left prefrontal activity (x = −54, y = 24, z = 32, BA46, k = 31, Z = 3.29, p_{FWE-corrected} = 0.041*); an interaction between BDNF rs6265 genotype, methylation, and hypoxia exposure on left prefrontal activity (Fig. 4A). More specifically, in one prefrontal locale, greater methylation is associated with attenuated prefrontal activity in subjects exposed to hypoxia compared with subjects not exposed in the context of ValVal genotype but not in ValMet subjects (x = −54, y = 24, z = 32, BA46, k = 26, Z = 3.24, p_{FWE-corrected} = 0.047*; Fig. 4B-C; difference test between ValVal subjects with and without hypoxia exposure: Z = 2.38, P = 0.008**). On the contrary, in the other prefrontal locale, lower methylation is associated with attenuated prefrontal activity in subjects with hOCs compared with subjects without hypoxia exposure in the context of ValMet genotype but not in ValVal subjects (x = −54, y = 32, z = 6, BA46, k = 14, Z = 3.40, p_{FWE-corrected} = 0.044*; Fig. 4D-E; Difference test between ValMet subjects with and without hypoxia exposure: Z = −2.32, P = 0.01*). These results suggest that prefrontal activity during WM is predicted by the interaction of hOCs exposure and methylation to describe opposite relationships in ValVal and ValMet subjects. Since attenuated prefrontal activity for equivalent gains in WM accuracy is considered indicative of increased efficiency of prefrontal cortex,^14,39,40 the opposite methylation changes associated with hOCs exposure in the 2 genotypes (increased methylation in ValVal and blunted methylation in ValMets) were also associated with increased activity, i.e., impaired prefrontal efficiency, which represents a phenotype characteristic of patients with schizophrenia and their relatives.^14,39,40

Figure 4. Interaction between BDNF rs6265 genotype, methylation in PBMCs, early life exposure to hypoxia (hOCs), and prefrontal activity during working memory (WM) in healthy subjects. a: 3D rendering of the interaction between rs6265 genotype, methylation and hOCs on BOLD fMRI response in prefrontal cortex of ValVal (N = 93) and ValMet (N = 48) subjects. Color bar represents F-values. b-c: Scatterplots of the interaction in BA 46 (x = −54, y = 24, z = 32) showing that increased methylation is associated with attenuated prefrontal activity in ValVal subjects exposed to hypoxia (B), while no significant relationship emerged in ValMet subjects (C). d-e: Scatterplots of the interaction in BA 46 (x = −54, y = 32, z = 6) showing that methylation is positively associated with prefrontal activity in ValMet subjects exposed to hypoxia (E), while no significant relationship emerged in ValVal subjects (D). See text for statistics.

Relationship between BDNF rs6265 genotype, methylation and genetic risk for schizophrenia

We further investigated whether BDNF rs6265 methylation can be inherited and is associated with genetic risk for schizophrenia. We analyzed BDNF rs6265 methylation in a sample of families with schizophrenia and we calculated heritability from the regression slope of offspring methylation on the average methylation of parents. Interestingly, we estimated in a cohort of 115 families that BDNF rs6265 methylation is not only sensitive to environmental exposures but can also be partially inherited (h² = 0.2; Fig. S2). Next, we analyzed whether rs6265 genotype and methylation interact on risk for schizophrenia, by comparing healthy subjects with patients with schizophrenia, with siblings and with parents of patients in 3 separate analyses. Consistent with the above results, we found an interaction between rs6265 genotype and methylation on schizophrenia risk, when comparing healthy subjects whether with siblings (N = 384, F_5,378 = 10.39, P = 0.001**), or with parents (N = 461, F_5,455 = 9.47, P = 0.002**) or with patients with schizophrenia (N = 406, F_5,400 = 4.01, P = 0.04*; Fig. 5A; Supplemental file 1). Specifically, rs6265 methylation in ValVal subjects is lower in healthy subjects, compared with siblings (post hoc with Tukey HSD: P = 0.01*), parents (P = 0.02*) and patients (P < 0.01**), so that in ValVal the methylation changes associated with hypoxia are also associated with schizophrenia risk. We also found a significant association between BDNF rs6265 methylation and schizophrenia risk in ValMet subjects. In this case, healthy subjects have greater methylation compared with siblings (P = 0.03*) and parents of patients (P = 0.02*), so that also in ValMets the methylation changes associated with hypoxia are associated with schizophrenia risk. However, ValMet patients are not significantly different from healthy subjects (P = 0.79). Since we cannot exclude that these effects may have been confounded by treatment with antipsychotics which alter BDNF levels in patients,²¹ we also evaluated in vitro the potential effect of treatment with haloperidol 1 μM for 2 and 5 hours on rs6265 methylation, in PBMCs of ValVal and ValMet healthy subjects (Fig. 5B). Notably, rs6265 genotype and haloperidol treatment interact on methylation, so that haloperidol significantly increases BDNF rs6265 methylation in ValMet subjects, but it does not in ValVal (N = 17, 11 ValVal and 6 ValMet; Factorial ANOVA: F_2,44 = 5.65, P = 0.006**; post hoc with Tukey HSD: “ValVal basal” vs. “ValVal Haloperidol 2 h”, P = 0.8; “ValVal basal” vs. “ValVal Haloperidol 5 h”, P = 1, “ValVal Haloperidol 2 h” vs. “ValVal Haloperidol 5 h”, P = 0.9; “ValMet basal” < “ValMet Haloperidol 2 h”, P = 0.02*; “ValMet basal” < “ValMet Haloperidol 5 h”, P = 0.002**, “ValMet Haloperidol 2 h” vs. “ValMet Haloperidol 5 h”, P = 0.92; Fig. 5B). These results suggest that greater levels of methylation in ValMet patients may be related to antipsychotic treatment (Fig. 5B), although further experiments are necessary to address the effect of chronic exposure to this and other antipsychotics on DNA methylation. Moreover, since information on hOCs exposure were not available in our sample of families with schizophrenia, we cannot exclude that the levels of methylation in patients, parents and siblings were also related to these or other risk factors. However, the comparison of rs6265 methylation levels of controls with methylation of parents and siblings of patients, in the context of ValVal and ValMet genotype, indicates that the relationship between methylation of the Val allele and genetic risk for schizophrenia is opposite between the 2 genotype groups. Therefore, our data show that methylation changes in the rs6265 region are associated with hOCs exposure, WM accuracy and related prefrontal activity, and genetic risk for schizophrenia, in a distinctive manner, depending on the rs6265 genotype. Specifically, in individuals with the ValVal genotype hOCs exposure is associated with increased methylation, and enhanced methylation—in these individuals—means impaired WM accuracy, increased prefrontal activity during fMRI (i.e., blunted prefrontal efficiency) and increased genetic risk for schizophrenia. On the other hand, in ValMet subjects the same insult is associated with decreased methylation, but blunted methylation turns out to be associated, in this genotype, with, again, impaired WM accuracy and increased genetic risk for schizophrenia.

Figure 5. BDNF rs6265 methylation in PBMCs and schizophrenia risk. a: Bargraph (mean ± s.e.m.) of rs6265 methylation (T scores) in healthy subjects (168 ValVal, 77 ValMet), patients with schizophrenia (122 ValVal, 40 ValMet), siblings (97 ValVal, 43 ValMet), parents of patients (148 ValVal, 69 ValMet): in ValVal subjects, rs6265 methylation is lower in healthy subjects, compared with siblings, parents and patients with schizophrenia. In ValMet subjects, rs6265 methylation is greater in healthy subjects, compared with siblings and parents of patients, while it is not significantly different compared with patients. b: Bargraph (mean + s.e.m.) of the effect of treatment with haloperidol 1 μM for 2 and 5 hours on rs6265 methylation changes in PBMCs of ValVal (N = 11) and ValMets (N = 6). See text and Supplemental file 1 for statistics.

Relationship between rs6265 genotype, methylation and BDNF levels in serum

To be able to attribute any potential relevance to methylation changes it is essential to assess their relationship with gene expression. Thereby, we analyzed the relationship between rs6265 methylation and serum BDNF levels measured with ELISA in a subsample of healthy subjects.⁴¹ A multiple regression, with BDNF levels as dependent variable and BDNF rs6265 genotype and methylation as predictors, revealed a significant interaction between rs6265 genotype and methylation on total BDNF expression in serum (N = 39: 23 ValVal, 16 ValMet; t = 4.37, P < 0.001***; Fig. 6A). Specifically, greater methylation is associated with greater levels of BDNF in ValVal subjects (t = 3.04, P < 0.01**), while this relationship is opposite in ValMet subjects (t = -3; P = 0.01*). Univariate results also indicate an effect of rs6265 genotype on BDNF levels, which are greater in ValMet compared with ValVal subjects (t = 4.3; P < 0.001***; see Supplemental file 1 for detailed results).

Figure 6. BDNF rs6265 methylation and gene regulation. a: Scatterplots of the correlations between rs6265 methylation (T scores) and Total BDNF in serum measured with ELISA (N = 39: 23 ValVal and 16 ValMet): increased methylation in PBMCs is correlated with increased expression in ValVal subjects and with attenuated expression in ValMet heterozygotes. b: Scatterplot of the correlation between methylation of rs6265 region in PBMCs and in prefrontal cortex in a group of transgenic mice carrying the human mutation (N = 17). c: Bargraph (mean + s.e.m.) of the effects of rs6265 genotype on binding of hypoxia-related proteins to the rs6265 region (N = 20: 11 ValVal and 9 ValMet). See text and Supplemental file 1 for statistics.

On the other hand, we did not detect a significant relationship between BDNF rs6265 methylation and mRNA expression in PBMCs (P < 0.2, not shown). This may raise the concern that BDNF rs6265 methylation represents something indirectly associated with schizophrenia, rather than a risk factor. However, since the expression of BDNF mRNA does not reflect the rate of protein synthesized,⁴² the relationship of BDNF rs6265 methylation with BDNF protein level may also be the result of the contribution of other elements, such as the expression of non-coding transcripts, affecting the level of BDNF protein by acting at a post-transcriptional level.^43,44 Moreover, since human platelets represent a main source of serum BDNF protein but not of BDNF mRNA, serum BDNF has been postulated not to originate from megakaryocyte precursor cells, while potential sources include CNS^41,45,46; indeed, it has been shown that BDNF can readily cross the brain-blood barrier.⁴⁷ The relationship between rs6265 methylation in PBMCs and serum BDNF levels is therefore compatible with a potential link between rs6265 methylation in PBMCs and in brain.

Relationship between brain and blood BDNF methylation in mice

We estimated the relationship between BDNF rs6265 methylation in PBMCs and in brain in transgenic mice carrying the human rs6265 SNP and the analysis revealed a negative correlation between PFC and PBMCs methylation – greater PFC methylation of the rs6265 Val allele is correlated with lower methylation in PBMCs of both ValVal mice (N = 7; Rho = −0.78 ; P = 0.04; Fig. S3A) and ValMet mice (N = 9; Rho = −0.8 ; P = 0.01; Fig. S3B). Similarly, the correlation between PFC and PBMCs methylation of a CpG close to rs6265 (hg19 position: Chr11: 27,679,922-3) is negative in the whole sample (N = 17; Rho = −0.55 ; P = 0.02; Fig. 6B). These results imply that methylation in PBMCs can also be used as a peripheral proxy of PFC methylation of rs6265, although the direction of the correlations demonstrated above may be reversed when considering rs6265 methylation in PFC. However, further studies are necessary to evaluate whether the relationship between rs6265 methylation in PBMCs and phenotypes relevant for schizophrenia is related to brain-blood correlation of methylation levels or is simply due to rs6265 methylation in PBMCs being an “epigenetic fossil” that keeps trace of early events relevant for brain development.

Relationship between rs6265 genotype and DNA-protein binding

Since a differential effect of methylation can be mediated by direct interference with DNA-protein binding,⁴⁸ which can be affected by genotype-specific changes,⁴⁹ we next addressed the potential relationship between BDNF rs6265 genotype and protein/transcription factors (TFs) binding. The rs6265 region, identified as a DNAse I hypersensitive site in brain,⁵⁰ is potentially interesting for the binding of TFs, since the Val allele shows putative binding sites for HIF1α, BHLHB2 (also known as DEC1 or BHLHE40), and CREB. These sites are disrupted by the G(Val)/A(Met) substitution, which also creates a binding site for MITF on the Met allele. Previous studies have shown how the interaction between all these factors is complex and likely involves BDNF, the demethylating protein GADD45B, and hypoxia-related mechanisms. For example, binding of CREB⁵¹ and GADD45B⁵² in the BDNF region upstream of the rs6265 SNP has been detected. Moreover, BHLHB2 and GADD45B are regulated by hypoxia, since they have a high-stringency HIF1 binding site,^53,54 and BHLHB2 represses expression of MITF,⁵³ a TF stimulating HIF1α expression.⁵⁵ Studies in mice have also proved how the basic helix-loop-helix protein BHLHB2 is regulated by neurotrophins and modulate BDNF transcription.⁵⁶ As a consequence, we performed an experiment of chromatin immunoprecipitation (ChIP), with the aim of verifying whether these proteins and TFs related to hypoxia could bind the Val allele in the rs6265 region, potentially interacting in a different way in ValVal and ValMet subjects. Our data confirmed that this region may bind GADD45B, HIF1α, BHLHB2, CREB, and MITF (Fig. 6C, detailed results in Supplemental file1). In addition, ANOVA showed that BDNF rs6265 genotype was associated with differential binding of BHLHB2, which was greater in the context of ValMet genotype, compared with the ValVal (N = 20: 11 ValVal and 9 ValMet; F_1,18 = 5.71; P = 0.02; Fig. 6C). Other TFs showed a similar trend, without reaching statistical significance for the differential binding in ValVals and ValMets. In addition, the greater binding of MITF in ValMets compared with ValVals is likely less relevant, since the putative binding site for this TF is present only on the Met allele.

All these findings indicate that rs6265 methylation of the BDNF Val allele has a differential relationship with BDNF levels and environmental exposures in ValVal and ValMet subjects, likely through a different interaction with transcription factor binding, which may be also particularly important in early development, when TFs are distributed in a concentration gradient.⁵⁷

Limitations of the study

A limitation of our study is that the assessment of hOCs relied solely on maternal recall as the source of information about the exposure. Unfortunately, this is a limitation shared by most of the literature on obstetric complications.²⁶ However, we believe to have controlled potential bias in our protocol by using a standardize questionnaire developed by previous published reports,⁵⁸ by employing the McNeil–Sjöström scale,²⁷ and by excluding from the analysis all the subjects with uncertain information.

Another limitation of our study is represented by the analysis of only a single region of the epigenome. Such approach was chosen to detect a potential interaction between genome, epigenome and environment, which would have gone undetected in a common whole-epigenome approach. By providing evidence of an ‘epistatic’ interaction between genetic and epigenetic variation, our data raise the possibility that analyses on whole genome and epigenome may give partial information, not taking into account that epigenetic changes can actually have opposite meaning depending on genotype, and vice versa.

Conclusions

Our data display a dynamic interplay between genome, epigenome, and environment on prefrontal function and risk for schizophrenia. Our main finding is the opposite relationship of BDNF rs6265 methylation with phenotypes and environmental factors relevant for schizophrenia in ValVal and ValMet subjects. More specifically, DNA methylation of rs6265 is differentially associated in ValVal and ValMet subjects with hOCs exposure and with prefrontal behavior and activity. Furthermore, the results in healthy subjects, siblings and parents of patients indicate that rs6265 methylation is associated with genetic risk for schizophrenia differentially in the 2 genotypes and that it can also be partially inherited. Although ambiguous at a first glance, these results are consistent in showing that the opposite methylation changes associated in the 2 genotypes with a risk factor for schizophrenia, i.e., hOCs exposure, are also associated with intermediate phenotypes for schizophrenia and genetic risk for the disorder in a consistent way: specifically, in ValVal subjects enhanced methylation is associated with hOCs exposure, impaired prefrontal cognition and inefficiency of prefrontal cortex, and methylation is also greater in siblings and parents of patients compared with controls, while in ValMet subjects blunted methylation is associated with hOCs exposure, impaired prefrontal cognition and inefficiency of prefrontal cortex, and methylation is lower in siblings and parents of patients compared with controls.

These genotype-dependent in vivo findings suggest that other molecular factors play a role. Indeed, rs6265 genotype affects binding of transcription factors and the relationship between DNA methylation and serum BDNF levels, so that the methylation changes associated with hypoxia are also differentially linked in both genotypes to BDNF levels. In this way, environmentally-sensitive DNA methylation modulates the effect of genetic variation, leading to risk phenotypes for complex disorders. More in general, our results indicate that opposite epigenetic, genotype-dependent changes may allow developmental plasticity to “adapt” the organism to environmental conditions, contributing to modulate complex phenotypes above and beyond genetic variation.

Material and methods

Subjects, methylation analysis and genotyping

Two hundred and forty-four healthy subjects, 162 patients with schizophrenia, 140 siblings and 217 parents of patients entered the study, based on inclusion criteria and protocols specified elsewhere.⁵⁹ Briefly, all subjects were white Caucasians from the region of Puglia and provided written informed consent. The Structured Clinical Interview for DSM-IV was used to confirm diagnosis of schizophrenia for patients and to exclude any Axis I psychiatric disorder for siblings, parents of patients and healthy subjects. Exclusion criteria were presence of any neurological or medical condition, presence of head trauma with loss of consciousness and drug abuse within the past 6 months. All patients were on stable pharmacological treatment with antipsychotics. The Institutional Review Board of University of Bari “Aldo Moro,” Bari (Italy), approved protocols and procedures. DNA was extracted from PBMCs using QIAamp DNA Blood Midi Kit (Qiagen, Valencia, CA, US) and bisulfite treated as described previously.⁶⁰ Methylation analysis was performed with pyrosequencing, using primer sequences previously reported³⁶ and focused on CpG methylation sites in the region of the BDNF gene containing rs6265. A consensus LINE-1 sequence was also analyzed to estimate global DNA methylation,⁶⁰ in order to exclude any global effect (see Supplementary Results for LINE-1 methylation results). Methylation was expressed as percentage of methylated cytosines divided by the sum of methylated and unmethylated cytosines (%5mC).⁶¹ SNP rs6265 was genotyped using the pyrosequencing assay designed to interrogate percentage of methylation in this region. In addition, genotypes were double-checked with direct DNA sequencing, as previously described.⁶² Subjects with MetMet genotype were excluded from further analyses, after verified that, as expected, rs6265 methylation was around 0%. Lack of methylation of the Met allele is consistent with absence of a CpG site in this position, confirming the quality of bisulfite treatment conversion of unmethylated Cytosine in Uracil/Thymine. Supplemental file 1 contains information about the different samples analyzed.

Hypoxia-related Obstetric complications assessment

Obstetric Complications (OCs) refer to conditions occurring not only during labor-delivery but also during pregnancy and neonatal period. Specifically, OCs are here referred as “somatic complications and conditions occurring during pregnancy, labor-delivery and the neonatal period” experienced as an offspring with special focus on the CNS.²⁷ We assessed OCs exposure on the basis of interviews administered to mothers using a standard questionnaire developed from other published reports.⁵⁸ OCs data were rated using the McNeil–Sjöström scale for obstetric complications,²⁷ which assigns each OC a severity score on a scale of 1–6. We determined OCs exposure based on the presence of at least one serious OC. As in a previous report,³⁵ we adopted a strict definition of serious OC, i.e., McNeil-Sjostrom Scale score ≥5 . This score allowed identification of individuals exposed to hypoxia-related (hOCs) potentially harmful obstetric complications. The hOCs reported included bleeding during pregnancy, maternal diabetes, maternal infections, Rh incompatibility, adverse fetal position, cord around neck, delivery problems, extended labor duration, use of high forceps, emergency cesarean section, early gestational age at birth and preterm birth, very low birth weight, respiratory distress at birth, and neonatal anomalies. Since a limitation of our study is that the assessment of hOCs relied solely on maternal recall as the source of information about the exposure, we excluded from the analysis all the subjects with uncertain information.

Working memory (WM) task

During fMRI, all subjects completed a blocked paradigm of the N-back task.⁵ Briefly, ‘N-back’ refers to how far back in the sequence of stimuli the subject had to recall. The stimuli consisted of numbers (1–4) shown in random sequence and displayed at the points of a diamond-shaped box. There was a visually paced motor task, which also served as a non-memory guided control condition (0-back) that simply required subjects to identify the stimulus currently seen. In the WM conditions, the task required recollection of a stimulus seen 2 stimuli (2-back) previously while continuing to encode additionally incoming stimuli. All subjects were trained on the task before the fMRI session. The stimuli of the task were organized in a simple block design in which each block consisted of 8 alternating 0-back and 2-back WM condition lasting 4 m and 8 s. Stimuli were presented via a back-projection system and behavioral responses were recorded through a fiber optic response box, which allowed measurement of accuracy and reaction time for each trial.

fMRI data acquisition

Blood oxygen level-dependent (BOLD) fMRI was performed on a GE Signa 3T scanner (General Electric, Milwaukee, WI), equipped with a standard quadrature head coil. A gradient-echo planar imaging sequence (repetition time, 2000 ms; echo time, 28 ms; 20 interleaved axial slices; thickness, 4 mm; gap, 1 mm; voxel size, 3.75 × 3.75 × 5; flip angle, 90°; field of view, 24 cm; matrix, 64 × 64) was used to acquire 120 volumes while subjects performed the WM task. The first 4 scans were discarded to allow for T1 equilibration effect.

fMRI data analysis

Analysis of the fMRI data was completed using Statistical Parametric Mapping (SPM8; http://www.fil.ion.ucl.ac.uk/spm). Images, for each subject, were realigned to the first volume in the time series and movement parameters were extracted to exclude subjects with excessive head motion (> 2 mm of translation, > 2° rotation). Images were then re-sampled to a 2 mm isotropic voxel size, spatially normalized into a standard stereotactic space (Montreal Institute on Neurology, MNI, template) and smoothed using a 10 mm full-width half-maximum isotropic Gaussian kernel to minimize noise and to account for residual inter-subject differences. A box car model convolved with the hemodynamic response function (HRF) at each voxel was modeled. In the first-level analysis, linear contrasts were computed producing t statistical maps at each voxel for the 2-back condition, assuming the 0-back condition as a baseline. All individual contrast images were entered in a second level random effects analysis. A multiple regression was performed entering BDNF rs6265 genotype, BDNF rs6265 methylation and hypoxia exposure scores as predictors. Because of our strong a priori hypothesis about the dorsolateral prefrontal cortex (DLPFC), we used a statistical threshold of P < 0.05, family wise error small volume corrected using as volume of interest the WFU_PickAtlas Brodmann's areas in which significant clusters were located (BA46).⁶³ Because we did not have a priori hypotheses regarding brain activity outside of DLPFC, we used a statistical threshold of P < 0.05, FWE - corrected for whole-brain comparisons.

PBMCs stimulation with Haloperidol

To evaluate the potential role of antipsychotic treatment on BDNF rs6265 methylation, we assessed DNA methylation in PBMCs of healthy subjects (11 ValVal, 6 ValMet) following in vitro challenge with haloperidol. Briefly, blood (20 ml) was collected and PBMCs were isolated by Ficoll density gradient (ICN, Biomedical, Inc.), as previously described.⁵ After counting, 3.5 × 10⁶ cells tube were resuspended in fresh RPMI 1640 medium (Gibco) (pH 7.5) with 15% FCS and HEPES (Sigma) 10 mM, and incubated at 37°C with haloperidol 1 µM (Janssen Pharmaceutical) for 0 (baseline), 2, and 5 h. Subsequently, DNA was extracted and BDNF rs6265 methylation was analyzed with pyrosequencing.

BDNF expression measurement

In order to assess a potential relationship between methylation and expression, we measured BDNF protein levels in serum, which are likely derived by CNS sources and are altered in patients with schizophrenia.⁴¹ Blood (10 ml) was collected between 8:00 and 9:30 AM in anticoagulant-free tubes and maintained at RT for 1 h, followed by 1 h at 4°C. After centrifugation at 2000 g for 10 min at 4°C, serum samples were stored up to 1 month at −20°C and then analyzed in triplicate at the same time. Sera were diluted 1:50 in sample buffer and total BDNF was quantified using an ELISA kit (BDNF Emax immunoassay system, Promega) in a microplate reader (Anthos Labtec Instrument) set at 450 nm. We also analyzed the relationship between BDNF rs6265 methylation and mRNA expression, in the same sample. BDNF mRNA levels were assessed using the comparative CT method with β-actin as reference (control) gene, using TaqMan® Gene Expression Assays (Applied Biosystems, Cat. #4331182) specific for the following transcripts: NM_170733.3 (Assay ID: Hs00380947_m1); NM_170731.4 (hs00538277-m1); NM_170732.4 (hs00538278-m1); NM_001709.4 (hs00156058-m1).

Correlation of BDNF methylation in brain and PBMCs

We analyzed the correlation between BDNF rs6265 methylation in PBMCs and in PFC in 7 homozygous (ValVal) and 9 heterozygous (ValMet) mice described by Chen,⁶⁴ which reproduce the phenotypic hallmarks described in humans with the variant allele. Mice were maintained on an inbred C57BL/6 background. The animals were housed under standard conditions (12-h light/dark cycle with food and water available ad libitum) and all studies were performed in adult mice. All animal handling and experimental procedures were performed in accordance with the EC (EEC Council Directive 86/609 1987), the Italian legislation on animal experimentation (Decreto Legislativo 116/92), and the National Institutes of Health Guide for the Care and Use of Laboratory Animals. All efforts were made to minimize animal suffering and to reduce the number of animals used. Samples of trunk blood from each mouse were collected in tubes containing sodium citrate and then rapidly processed for DNA analysis (see below), whereas frontal lobes were dissected, frozen on dry ice and stored for further analyses at −80°C. Genomic DNA from whole plasma was extracted using ReliaPrep™ Blood gDNA Miniprep System kit (Promega, Italy), while DNA from brain tissue was isolated using phenol/chloroform extraction method. Analysis of mice BDNF rs6265 methylation was performed with pyrosequencing, as described in humans.

Quantitative chromatin immunoprecipitation assay (ChIP)

Given the differential genotype-dependent relationship between rs6265 methylation and BDNF expression, we performed a ChIP experiment in order to investigate whether the BDNF region containing rs6265 can bind transcription factors and proteins in a genotype-dependent manner. We selected DNA-binding proteins based on bioinformatics predictions and previous research, as described in the main text. Bioinformatics tools (www.genomatix.de) suggest how rs6265 genetic variation can affect binding of transcription factors so that the G(Val)/A(Met) substitution abolishes a potential binding site for HIF1α, BHLHB2 (also known as DEC1 or BHLHE40), and CREB while it creates a binding site for MITF. PBMCs were isolated from 20 healthy individuals (11 ValVal, 9 ValMet), as previously described.⁵ Protein bound to DNA was cross-linked by PBMCs with 1% formaldehyde at room temperature, stopping the reaction 10 min later with the addition of 2.5 M glycine to a final concentration of 125 mM, followed by 5 min incubation at room temperature. ChIP assays were performed using the EpiQuik^TM chromatin immunoprecipitation kit from Epigentek Group Inc. (Brooklyn, NY) starting from ˜0.5 × 10⁶ PBMC cells. Antibodies used for Protein-DNA immunoprecipitation were: anti-HIF1α (antibody provided from Dr. C. W. Pugh, Center for Cellular & Molecular Physiology, University of Oxford, Oxford, United Kingdom⁶⁵), anti-MITF and anti-GADD45B (Aviva System Biology, San Diego, CA), anti-BHLHB2/DEC1 (Bethyl Laboratories, Inc. Montgomery, TX USA), anti-CREB (Merck Millipore Headquarters, Billerica, MA), and normal mouse IgG as a negative control antibody. DNA from these samples was subjected to quantitative PCR analyses, using Power SYBR® Green PCR Master Mix (Life Technologies Corporation, Carlsbad, California) in a Chromo4 Real Time thermocycler (BIORAD). Amplification of the BDNF promoters fragment was performed using the primers: pBDNFf (forward) 5′-CCAAGGCAGGTTCAAGAGG-3′ and pBDNFr (reverse) 5′-CGAACTTTCTGGTCCTCATCC-3′ amplifying a 90 bp fragment including rs6265 SNP. The quantitative PCR conditions were: 95°C for 10 min followed by 40 cycles of 95°C for 15 s, 62°C for 1 min. All PCR signals from immunoprecipitated DNA were normalized to PCR signals from non-immunoprecipitated input DNA. The signals obtained by precipitation with the control IgG were subtracted from the signals obtained with the specific antibodies. Results are expressed as percentage of the input.⁶⁶ Calculations were performed using the average values of at least 3 independent experiments.

Statistical analysis

All statistical analyses—except for fMRI—were performed in the R environment. We used ANOVAs in order to analyze the effect of rs6265 genotype on methylation and on DNA-binding of proteins and transcription factors; we used ANCOVAs (with age and sex as covariates) to analyze the interaction between rs6265 genotype and hOCs on methylation, and the relationship between rs6265 genotype, methylation and schizophrenia risk and diagnosis. To assess in vitro effects of the challenge with haloperidol on rs6265 methylation, we performed an ANOVA with rs6265 genotype and incubation time with haloperidol as independent categorical variables (0 h/basal, 2 h, 5 h) and rs6265 methylation changes as dependent variable. Heritability of rs6265 methylation was estimated from the regression slope of offspring methylation on the average methylation of the parents.⁶⁷ We used multiple regressions in order to analyze the interaction between rs6265 genotype, methylation and hOCs on WM accuracy; the interaction between rs6265 genotype and methylation on BDNF levels. Finally, we used Spearman correlation in order to analyze the relationship between rs6265 methylation in PBMCs and in PFC in mice. Statistical models for the analyses are also reported in Supplemental file 1.

Disclosure of potential conflicts of interest

Dr. Bertolino is a consultant of Hoffman-La Roche Ltd. All the authors declare no conflict of interest.

Acknowledgments

We thank Dr. Chris W Pug and Dr. Zuzana Bencokova for providing anti-HIFα antibody; Riccarda Lomuscio, BA, Maria Teresa Attrotto, MD, Lucia Colagiorgio, MD, Giuseppe Rizzo, MD, for helping with data acquisition, and Dr. Daniel R. Weinberger for helpful discussions. This work was supported by Fondazione Con Il Sud “Capitale Umano ad Alta Qualificazione” grant (awarded to AB); Brain & Behavior Research Foundation Independent Investigator grant (AB); Regione Campania l.5 and EPIGEN Flagship Project CNR grant.