The association of histone modification changes with autism spectrum disorder (ASD) has not been systematically examined. We conducted a histone acetylome-wide association study (HAWAS) by performing H3K27ac chromatin immunoprecipitation sequencing (ChIP-seq) on 257 postmortem samples from ASD and matched control brains. Despite etiological heterogeneity, ≥68% of syndromic and idiopathic ASD cases shared a common acetylome signature at >5,000 cis-regulatory elements in prefrontal and temporal cortex. Similarly, multiple genes associated with rare genetic mutations in ASD showed common “epimutations.” Acetylome aberrations in ASD were not attributable to genetic differentiation at cis-SNPs and highlighted genes involved in synaptic transmission, ion transport, epilepsy, behavioral abnormality, chemokinesis, histone deacetylation, and immunity. By correlating histone acetylation with genotype, we discovered >2,000 histone acetylation quantitative trait loci (haQTLs) in human brain regions, including four candidate causal variants for psychiatric diseases. Due to the relative stability of histone modifications postmortem, we anticipate that the HAWAS approach will be applicable to multiple diseases.

To address this lack of knowledge, we globally interrogated the histone acetylomes of enhancers in a large cohort of ASD and control samples by analyzing tissue from three brain regions postmortem: prefrontal cortex (PFC), temporal cortex (TC), and cerebellum (CB). These brain regions were chosen due to the role of frontal and temporal lobe in social cognition and the cerebellar dysfunction observed in some animal models of ASD (). H3K27ac was selected as the representative acetylation mark because it highlights active enhancers and promoters () and is also correlated with gene expression and transcription factor binding (). We used the data to define aberrantly acetylated enhancer and promoters in ASD brain and thereby characterize commonly altered pathways, upstream regulatory factors, and developmental dynamics of affected loci. In addition, we used chromatin immunoprecipitation sequencing (ChIP-seq) reads to call SNPs within enhancers and promoters. We then used the genotype-independent signal correlation and imbalance (G-SCI) test () to detect haQTLs in regulatory regions and assessed their relationship to known psychiatric disease-associated variants. This dataset from post-mortem human brains will provide a rich resource for future molecular analyses of ASD and serve as proof of principle for the HAWAS approach, which can be applied to a wide variety of human diseases.

Large-scale gene expression studies on ASD postmortem brain regions, as well as the prenatal and postnatal developing brain, suggest that alterations in common molecular pathways such as transcriptional regulation, synaptic function, and immunity may occur during brain development and contribute to ASD pathophysiology (). How the genetic and environmental heterogeneity is translated into shared molecular pathways is not well understood. In addition, most of the efforts have so far focused on gene expression and genetic changes in coding regions. Many of these coding variants are extremely rare and account only for a small proportion of ASD cases (). Therefore, it has been proposed that epigenetic changes caused by non-coding genetic variation or by environmental insults might contribute to ASD (). An attempt to characterize epigenomic changes in patients is thus likely to provide novel insights into the etiology of ASD (). Thus far, epigenome-wide association studies (EWAS) of psychiatric and other diseases have mostly focused on DNA methylation (). In contrast, little is known about histone modification changes in psychiatric disease () or the genetics of population variation in histone modification ().

Autism spectrum disorder (ASD) is a collection of neuro-developmental disorders characterized by deficits in social interaction and social communication, along with restricted and repetitive behavior patterns. DNA sequence variation affecting the function of several hundred genes has been implicated in the etiology of ASD at various levels of significance (). These genetic changes include copy number variation, chromosomal rearrangements, and also rare single-base pair mutations in coding genes (). In addition, environmental risk factors such as chemical toxins and maternal infection during gestation are thought to play a role in the occurrence of ASD (). Thus, the etiology of ASD is complex and multifactorial, including both genetic and environmental components.

GWAS analyses have not so far uncovered statistically significant ASD-associated variants that have been replicated in independent analysis at a genome-wide level. We therefore intersected the haQTL set with genome-wide significant (p ≤ 5e-8) variants known to be associated with shared aspects of five psychiatric disorders: schizophrenia, bipolar disorder, major depressive disorder, ASD, and attention-deficit/hyperactivity disorder (). While this GWAS set was too small to test for statistical enrichment near haQTLs, we did uncover two instances where brain haQTLs were strongly linked (R≥ 0.8) to disease-associated variants ( Table S7 ). Most notably, an haQTL (rs4765905) in an intron of the syndromic ASD gene CACNA1C was linked to multiple psychiatric disease-associated SNPs within the locus ( Figure 5 B). Based on Hi-C data from GM12878 cells (), the putative enhancer containing this haQTL SNP was predicted to form a long-range loop to the CACNA1C promoter, suggesting that it could exert its influence on psychiatric disease by modulating the chromatin state of CACNA1C. In addition, we intersected haQTLs with 128 SNPs associated with schizophrenia in a recent large-scale meta-analysis of schizophrenia (). This analysis revealed two additional haQTLs strongly linked to psychiatric disease-associated variants ( Table S7 ). For example, we found that the haQTL SNP rs8054791 was linked to the schizophrenia-associated variant rs9922678 in an intron of GRIN2A, a glutamate receptor gene that has also been associated with ASD ( Figure 5 C).

To identify haQTLs in the three human brain regions, we used the G-SCI pipeline that was previously validated on lymphoblastoid cell lines (). The pipeline uses ChIP-seq reads to call DNA sequence variants in active regulatory regions, followed by filtering to remove low-confidence variants ( STAR Methods ). By analogy to exome sequencing, this stage of the pipeline can be termed “regulome sequencing.” A unique aspect of the G-SCI method is that called variants need not be explicitly genotyped. Rather, counts and base qualities of reference- and alternative-allele ChIP-seq reads are used to infer genotype likelihoods. These likelihoods are then used to compute the haQTL p value of the variant using the G-SCI test, which maximizes statistical power by combining information from peak height variability and allelic imbalance across all individuals within the cohort. In order to separate the cis effect of regulatory SNPs from more general disease effects, we first adjusted ChIP-seq peak heights by regressing out the diagnosis variable (ASD versus control). We then applied the G-SCI test to called SNPs and identified ∼2,000 haQTLs in each of the three brain regions ( Figure 5 A; Table S7 ). Note that these haQTLs are not specific to ASD. Rather, they represent region-specific regulatory variation in the general population.

(C) A SNP within an intron of GRIN2A (chr16: 9,936,580-9,951,221) is an haQTL (Q = 2.2e-2) in TC and CB and is in LD with a genome-wide significant SNP (rs9922678; R 2 = 0.91) for schizophrenia. Histone acetylation tracks (chr12: 2,343,551-2,354,513), read depth analysis (bar graph), and peak height boxplots indicate that the reference “A” allele has higher histone acetylation than the non-reference “G” allele. All acetylation tracks are plotted on the same fold-enrichment scale (y axis: 0–50).

(B) A SNP within an intron of CACNA1C is an haQTL (Q = 1.5e-4) in CB and is in LD with four genome-wide significant SNPs from GWAS of five psychiatric disorders ( Table S7 ). Histone acetylation tracks (chr12: 2,343,551-2,354,513), read depth analysis (bar graph) and peak height boxplots indicate that the reference “G” allele has higher histone acetylation than the non-reference “C” allele. All acetylation tracks are plotted on the same fold-enrichment scale (y axis: 0–120).

Noncoding genetic variants that affect disease susceptibility potentially act via a gene regulatory mechanism (). Because histone acetylation serves as a measure of gene regulatory function, such variants are also likely to influence acetylation levels. It is therefore instructive to identify histone acetylation QTLs (haQTLs), which are defined as genetic variants that correlate with population variation in histone acetylation (). As we and others have previously shown (), haQTLs can be used to prioritize causal variants within disease-associated loci.

It has been shown that genes upregulated during early postnatal development are often differentially expressed in adolescent and adult ASD brain (). We therefore asked whether early postnatal genes might also be enriched for the ASD-related acetylation changes we detected in older subjects (≥10 years old). Using a database of human RNA-seq profiles (), we defined the 2,000 genes most upregulated at each developmental stage (fold change relative to median expression across all stages). We then tested for enrichment of DA peaks near each such gene set. This analysis was performed separately for PFC and TC, using expression profiles from the corresponding regions of the developing human brain. As expected, ASD-Up DA peaks in the adult (more precisely, ≥10 year) brain were significantly overrepresented near adult-upregulated genes ( Figure 4 ). Surprisingly, however, we found even greater enrichment of ASD-Up DA peaks near genes upregulated at 10–12 months after birth, which corresponds to the stage of synapse formation, and neuronal maturation. In contrast, ASD-Down DA peaks did not show stage-specificity. Thus, although chromatin aberrations in ASD affect genes with a broad variety of developmental specificities, genes upregulated at or near 12 months after birth are particularly strongly associated with increased acetylation in ASD cortex.

(A) ASD-Up DA peaks in PFC are most significantly enriched near genes upregulated 1 year after birth. Bar height indicates enrichment Q value (FDR). Numbers above bars indicate fold enrichment (Q ≤ 0.05).

To identify transcription factors (TFs) that potentially mediate aberrant histone acetylation in ASD, we used the HOMER tool () to scan for TF-binding motifs within DA peaks. Most notably, we found strong enrichment of RFX motifs in ASD-upregulated peaks, both in TC and in PFC ( Figure 3 ; DA Up). RFX2 has a DA peak at its promoter and RFX3 contains an intronic DA peak ( Table S4 ). These two TFs are therefore the most likely candidates for driving increased acetylation in ASD. Three other TFs or TF families were enriched in DA Up peaks across both cortical regions: PAR bZIP, AP-1 and MEF2. Among the PAR bZIP candidate TFs, E4BP4 and HLF are the most promising, because their promoters are acetylated in cerebral cortex. Among the MEF2 factors, MEF2C is clearly the most prominent candidate, because the corresponding gene locus hosts five DA peaks from the downregulated list in PFC (p = 1.1e-4; Table S6 ). The nuclear receptor motif enriched in ASD-upregulated peaks in TC could relate to the glucocorticoid or mineralocorticoid receptor (GR or MR), because the corresponding gene promoters are marked by H3K27ac peaks. In contrast to the five to six motifs overrepresented in DA Up peaks in each cortical region, SOX was the only motif enriched in DA Down peaks in TC and ETS the only motif in PFC. These binding site enrichment results potentially indicate the presence of master TFs that drive dysregulation of groups of regulatory elements across the ASD genome.

We hypothesized that DA peaks might act as cis-regulatory elements for some of the genes thought to be causal for ASD. We therefore analyzed a curated list of 296 ASD genes for proximity to DA peaks (SFARI database, gene score ≤4) () ( Table S5 ). In order to control for biases in gene length and intergenic size, this analysis was performed using the same statistical procedure as in the GREAT tool () (hypergeometric test). In both PFC (p = 0.017) and TC (p = 0.025), peaks showing increased acetylation in ASD were significantly enriched near the ASD gene set ( Table S5 ). ASD-downregulated peaks in TC were also enriched, though not significantly (p = 0.055) and downregulated peaks in PFC showed no enrichment for known ASD loci. Cerebellar DA peaks were too few in number to analyze in this manner. Cortical DA peaks as a whole (union of the four DA peak sets) were clearly overrepresented near ASD genes (p = 7.6e-4; fold enrichment = 1.1).

Downregulated DA peaks also showed highly significant enrichment for specific functions ( Figure 2 C). Immune-related terms such as abnormal immune serum protein physiology and lymphatic system disease were most prominently enriched in this peak set in PFC, perhaps reflecting unique microglial () or lymphoid cell states () in ASD cortex. Downregulated peaks in TC showed similar immune-related enrichments ( Table S3 ). DA-Down peaks in the two cortical regions were also enriched near histone deacetylase genes, including HDAC2 and HDAC4 (). In particular, the syndromic autism gene HDAC4 () neighbored 16 downregulated DA peaks in TC and 4 in PFC ( Figure 2 D; Table S4 ; note that some of the DA peaks in Table S4 have been annotated with the names of smaller genes within the introns of HDAC4). Chemokine gene loci were also enriched for DA-Down peaks in both cortical regions. In addition, when all genes in the genome were individually scored for enrichment near DA-Down peaks, the chemokine receptor CX3CR1 ranked among the top five in both PFC and TC ( Table 1 ). Multiple gut-related gene groups such as embryonic digestive tract morphogenesis and digestive/alimentary phenotype showed significantly reduced histone acetylation in ASD cortex. These gene sets included multifunctional morphogenetic genes such as FGFR2, chemokine ligand and receptor genes (including CX3CR1), and the HDAC SALL3 ( Table S3 ). Among the individually enriched genes near downregulated peaks in TC, the behavior-related gene GRB10 () was the most statistically significant (p < 1e-9).

Although ASD is known to be etiologically highly heterogeneous, we hypothesized that its diverse causal genetic and environmental perturbations could potentially converge on a small set of downstream pathways (). We therefore used the GREAT tool () to test for significantly enriched (fold change ≥ 1.5, Q value ≤ 0.01) gene categories in DA loci. GREAT maps regulatory elements to flanking genes based on proximity and uses the hypergeometric test to determine if the fraction of DA peaks near genes from any particular functional category is greater than expected by chance. Overall, DA peaks in PFC and TC showed very similar functional profiles. In both cortical regions, upregulated DA peaks were significantly enriched in neuronal functions known to be perturbed in ASD, including synaptic transmission, metal cation transport, epilepsy, and the glutamate receptor pathway () ( Figure 2 A; Table S3 ). Known ASD genes from these categories include CACNA1C and GRIN2B (), which flank five and seven DA peaks, respectively ( Figure 2 B). In light of the observation that zinc deficiency is common in ASD (), it is intriguing that upregulated DA peaks were also enriched in loci related to zinc ion homeostasis. Genes contributing to this result include the SLC30A5 zinc transporter gene, which flanks multiple DA peaks and has been shown to harbor rare single nucleotide variants in ASD ().

(D) ASD-Down DA peaks in the GRB10 (chr7: 50,657,760-50,861,159) and HDAC4 (chr2: 239,960,131-240,388,294) gene loci. Only the DA peaks closest to GRB10 and HDAC4 are visible in these genomic windows.

(B) ASD-Up DA peaks in the CACNA1C (chr12: 2,161,809-2,900,900) and GRIN2B (chr12: 13,427,172-14,303,010) gene loci. Only the DA peaks closest to CACNA1C and GRIN2B are visible in these genomic windows.

Having confirmed the robustness and consistency of the DA peak set, we constructed an ASD-specific global acetylome signature (AGAS), defined as the first principal component (PC1) of the corresponding peak height matrix. The strength of this signature in each brain sample was thus given by its PC1 score (X coordinate in Figure 1 C). In all three brain regions, ASD samples had significantly lower AGAS scores, and disease status explained 12%–63% of the score variance ( Figures 1 C and S5 ). Conversely, a simple threshold on the AGAS score could be used to predict disease status for 68%–95% of samples in the three brain regions ( Figure 1 C). Thus, ASD was associated with a coherent global shift in the histone acetylome.

PFC and TC gene expression levels have been comprehensively measured using RNA sequencing (RNA)-seq in a parallel study on the same cohort (). To investigate the consistency between chromatin aberrations and gene expression changes in ASD, we focused on promoter regions of differentially expressed (DE) genes (FDR ≤ 0.05; linear mixed model) (). We used the EFilter tool () to convert promoter histone acetylation profiles into expression estimates and then identified the subset of DE genes whose acetylation-based expression estimates were significantly divergent between ASD and control (Q ≤ 0.05; Wilcoxon test; Benjamini-Hochberg correction), after controlling for covariates as before. At these gene loci, promoter acetylation changes were significantly correlated with expression fold change, both in PFC (R = 0.33; p = 0.016) and in TC (R = 0.38; p = 0.0012) ( Figure 1 F). This analysis was not extended to CB, due to the lack of detectable DE genes in that tissue. Thus, while measurement noise and biological differences between chromatin variation and expression variation may have attenuated the similarity between the two, we nevertheless observed evidence of overall consistency between acetylation aberrations and expression dysregulation in ASD.

Of the 45 cases, 7 had a monogenic form of ASD, duplication 15q syndrome (dup15q; Figure 1 A), while the others had no detectable structural variants and were therefore idiopathic (). It is possible that individuals with syndromic ASD could have unique chromatin aberrations. We therefore defined DA peaks separately for syndromic and idiopathic ASD, relative to the same set of controls ( STAR Methods ). Remarkably, acetylation changes were highly concordant genome-wide between the two forms of ASD ( Figure 1 E, R = 0.88 in PFC, R = 0.87 in TC). To maximize statistical power, we therefore retained the original set of DA peaks based on all ASD samples (syndromic plus idiopathic).

To further characterize the overall consistency of the DA peak sets, we examined their overlaps. Over 45% of ASD-upregulated regions in PFC overlapped ASD-Up peaks in TC (p ≈ 0; Figure 1 D). The same was true of ASD-downregulated peaks. Moreover, the ASD-versus-control acetylation fold change was highly correlated between PFC and TC (R = 0.86; p ≈ 0). Thus, the chromatin dysregulation signature of ASD was highly consistent between the two cortical regions. Cerebellar DA peaks, on the other hand, showed only ∼5% overlap with same-direction cortical DA peaks.

To evaluate the likelihood of false-positive DA peaks, we repeated the entire procedure (initial DA peaks, discarding atypical samples, final DA peaks) after randomly permuting ASD and control labels. At the same false discovery rate (FDR) threshold (Q ≤ 0.05), permuted datasets generated fewer than 100 DA peaks, on average. Moreover, after 1,000 tries, none of the permutated datasets generated as many DA peaks as the true dataset ( Figure 1 B). Thus, the chromatin changes we detected in ASD samples were far in excess of what would be expected by chance.

In order to define the core set of chromatin aberrations in typical ASD cases, we employed a systematic mathematical criterion to exclude atypical samples ( STAR Methods ). Because the number of excluded samples was relatively small (5%–20% of total cases and 14%–20% of dup15q), acetylation fold changes were not substantially altered (PFC: R = 0.94, TC: R = 0.90, CB: R = 0.98; Figure S4 ). We then used the remaining (typical) samples to define the final set of DA peaks for each brain region. Strikingly, we detected 5,153 DA peaks in PFC and 7,009 in TC, indicating widespread, systematic histone acetylation changes in ASD cerebral cortex ( Figures 1 B and 1C). In contrast, only 247 DA peaks were detected in CB. The limited molecular pathology of ASD cerebellum is consistent with results from transcriptomic studies ().

Acetylation Fold Change between ASD and Control, Calculated Using All Samples Displayed on the Y Axis or Using Only Typical Samples Displayed on the X Axis, Related to Figure 1

The heights (aggregate read counts) of consensus peaks represent acetylation levels of cortical and cerebellar regulatory regions in each sample. We normalized these peak heights for GC content ( Figure S1 ) and distributional skews and then controlled for confounders by regressing out multiple biological covariates such as age, sex, and proportion of neurons and also multiple technical covariates ( STAR Methods Figure S2 ). Corrected peak heights were used to define an initial set of differentially acetylated (DA) loci between ASD and control in each brain region (Wilcoxon rank sum test; Q ≤ 0.05, fold change ≥1.3). Based on acetylation levels at these DA loci, we measured inter-sample divergence and found that a small number of atypical ASD samples showed greater similarity to control and vice versa ( Figure S3 ). This was not surprising, given the tremendous etiological heterogeneity of ASD and previous findings from transcriptomic analysis (). Nevertheless, in the majority of cases, ASD acetylomes resembled each other more than they resembled control and vice versa ( Figure S3 ).

Scatterplot of median divergence between acetylomes in PFC (A), TC (B) and CB (C). In this analysis, the acetylome is defined as the vector of peak heights at DA peaks. x axis: median Euclidean distance to other control acetylomes; y axis: median Euclidean distance to other ASD acetylomes ( STAR Methods ). Red dots: ASD samples; blue diamonds: control samples. Solid line: Y = X; Dotted lines: Y = 1.05X and X = 1.05Y. ASD samples above the Y = 1.05X line and control samples below the X = 1.05Y line were defined as atypical samples.

Pearson correlation coefficient is shown at each grid point. After regressing out correlated confounding factors, the top 5 PCs correlated with none of the covariates except diagnosis. InsertSize: fragment median insert size; Dup: percentage of duplicated reads; Reads: sequencing depth; Peaks: number of peaks; Neuron: neuronal cell fraction.

In total, we performed 257 H3K27ac ChIP-seq assays on tissue samples from PFC, TC, and CB, in 94 individuals aged 10 years and above (45 ASD, 49 control; Figure 1 A; Table S1 ). Forty-eight ChIP-seq profiles were discarded based on data quality, resulting in a final acetylome set comprising 209 profiles ( STAR Methods Table S2 ): 81 from PFC (41 ASD, 40 control), 66 from TC (30 ASD, 36 control), and 62 from CB (31 ASD, 31 control). We used DFilter () to call peaks in each of the 209 ChIP-seq profiles and then defined two consensus peak sets: 56,503 cortical peaks (union of PFC and TC) and 38,069 CB peaks. Each consensus ChIP-seq peak defined a region of focal histone acetylation and thus represented a putative promoter or enhancer region.

(F) Correlation between fold change of differential acetylation in the promoter region (DA) and differential gene expression (DGE) in PFC and TC. Only significantly differential loci are shown (PFC:58 genes, TC:79 genes; STAR Methods ). Statistical significance of concordance in fold change direction was calculated using the hypergeometric test.

(E) Correlation between log2 fold change of dup15q and idiopathic samples in PFC and TC. The correlation coefficient and its p value were calculated as in (D). ChIP-seq peaks within the 15q duplication region are highlighted in red.

(D) Venn diagrams showing Overlap between DA peak sets from the three brain regions. The hypergeometric test was used to calculate p values, with the set of all peaks as background. The density plot shows the log2 fold change in TC versus PFC in the union of DA peaks. The corresponding p value was calculated assuming a t-distributed Pearson correlation coefficient.

(C) PCA of ChIP-seq peak heights (DA peaks only) in the three brain regions. Red dots, ASD samples; blue diamonds, control samples. Unfilled dots and diamonds indicate atypical samples. Variance explained by PC1 and PC2 are shown in parentheses. The vertical line, the threshold on ASD-specific global acetylome signature (AGAS) score.

(B) Number of DA peaks detected in the three brain regions. For comparison, the average number of DA peaks called in 1,000 randomized datasets (permutation of case-control labels) is shown. The p value was computed as the fraction of randomized datasets yielding at least as many DA peaks as the true dataset.

(A) Overview of post mortem tissues in three brain regions used in this study. H3K27ac ChIP-seq was performed on prefrontal cortex (PFC), temporal cortex (TC), and cerebellum (CB) samples from 45 ASD (A) and 49 control (C) individuals. ASD subjects with 15q duplication syndrome are highlighted in green.

Discussion

Voineagu et al., 2011 Voineagu I.

Wang X.

Johnston P.

Lowe J.K.

Tian Y.

Horvath S.

Mill J.

Cantor R.M.

Blencowe B.J.

Geschwind D.H. Transcriptomic analysis of autistic brain reveals convergent molecular pathology. Voineagu et al., 2011 Voineagu I.

Wang X.

Johnston P.

Lowe J.K.

Tian Y.

Horvath S.

Mill J.

Cantor R.M.

Blencowe B.J.

Geschwind D.H. Transcriptomic analysis of autistic brain reveals convergent molecular pathology. Despite etiological heterogeneity, our results indicate that shared aberrations in histone acetylation are widespread in ASD cerebral cortex: over 5,000 enhancer or promoter loci were systematically shifted up or down ( Figure 1 B). The fact that histone acetylation changes were broadly similar between PFC and TC indicates similarity in ASD mechanisms across cortical regions and also suggests that our results on differential acetylation are unlikely to represent methodological artifacts. Note that, as expected for a complex disorder with highly heterogeneous etiology, this global signature of chromatin alteration is not shared by all ASD samples ( Figure 1 C). An earlier transcriptomic study revealed a similar pattern of changes shared by many, but not all, ASD cases (). Nevertheless, the fact that the majority of patients conform to a single global epigenomic pattern indicates that the diverse causal mechanisms of ASD have shared downstream effects on the acetylome. In contrast to cerebral cortex, only 247 loci were found to be perturbed in cerebellum, indicating that the former is affected to a much greater degree. This disparity between ASD cerebrum and cerebellum has also been observed at the transcriptomic level (). Syndromic dup15q cases showed acetylome alterations that were highly correlated with those observed in idiopathic ASD (R ≥ 0.87 in cerebral cortex), suggesting that most chromatin aberrations are shared between idiopathic ASD and this syndromic form.

Krumm et al., 2015 Krumm N.

Turner T.N.

Baker C.

Vives L.

Mohajeri K.

Witherspoon K.

Raja A.

Coe B.P.

Stessman H.A.

He Z.-X.X.

et al. Excess of rare, inherited truncating mutations in autism. To examine the genetic basis of the epigenomic aberrations detected in ASD, we tested all high-coverage SNPs within DA peaks for genetic differentiation between patients and controls (chi-square test). The distribution of genetic differentiation p values was close to uniform (data not shown), suggesting that genetic variation in cis SNPs is not a major contributor to case-control acetylation differences at DA peaks. It is thus likely that ASD-specific differential acetylation is driven mostly by other factors such as environmental influences, SNPs in trans (at a different locus), indels, and larger chromosomal variants ().

Kumar et al., 2013 Kumar V.

Muratani M.

Rayan N.A.

Kraus P.

Lufkin T.

Ng H.H.

Prabhakar S. Uniform, optimal signal processing of mapped deep-sequencing data. Yen and Kellis, 2015 Yen A.

Kellis M. Systematic chromatin state comparison of epigenomes associated with diverse properties including sex and tissue type. Overall, acetylation changes in ASD cerebral cortex were significantly correlated with differential gene expression, consistent with the known functional consequences of these alterations in chromatin structure ( Figure 1 F). However, the majority of DA peaks did not lie next to DE genes. This is consistent with previous studies; we and others have shown that differences in chromatin state between two sample types are only moderately correlated with differential expression (). Differences in the sensitivity of ChIP-seq and RNA-seq at various loci could provide one explanation for this phenomenon. For example, post-mortem RNA degradation or low steady-state mRNA levels could reduce the detectability of DE genes in some cases, while low read mappability or occlusion of the acetylated epitope (for example) could limit the sensitivity of DA peak analysis at other loci. Moreover, noise levels could vary between the mRNA and chromatin readouts at individual loci, resulting in differential statistical power. Finally, although histone acetylation and gene expression are correlated in general, post-transcriptional regulation, other histone modifications, DNA methylation status, and the influence of additional regulatory elements within the same locus could all contribute to genuine biological differences between mRNA fold change and acetylation shifts. Thus, case-control chromatin profiling could serve as a valuable complement to the more common strategy of transcriptomic profiling by highlighting novel disease mechanisms.

Voineagu et al., 2011 Voineagu I.

Wang X.

Johnston P.

Lowe J.K.

Tian Y.

Horvath S.

Mill J.

Cantor R.M.

Blencowe B.J.

Geschwind D.H. Transcriptomic analysis of autistic brain reveals convergent molecular pathology. Bourgeron, 2015 Bourgeron T. From the genetic architecture to synaptic plasticity in autism spectrum disorder. Parikshak et al., 2013 Parikshak N.N.

Luo R.

Zhang A.

Won H.

Lowe J.K.

Chandran V.

Horvath S.

Geschwind D.H. Integrative functional genomic analyses implicate specific molecular pathways and circuits in autism. We found evidence for shared pathways and functional themes among DA loci in ASD cerebral cortex ( Figure 2 ). Among loci with increased H3K27ac, there was strong enrichment for genes related to ion channels, synaptic function, and epilepsy/neuronal excitability, all of which have previously been shown to be dysregulated in this disorder (). Moreover, these adult DA loci were strongly enriched for genes developmentally upregulated at or around 12 months of life ( Figure 4 ), which coincides with the peak of early experience-dependent synaptogenesis. A similar temporal enrichment has also been observed for cerebral DE genes in ASD (). Loci with decreased acetylation in ASD also converged on shared functional categories, such as digestive tract morphogenesis, chemokine signaling, HDAC activity, and immune processes related to microglia. Note that it is possible for functional categories to appear systematically enriched in DA peaks merely because of the contribution of a single highly enriched “jackpot” gene. However, our results are likely to be robust to such artifacts, because we discarded functional terms that had fewer than five genes near DA peaks and then manually inspected the remaining top hits (shown in Figure 2 ) for jackpot effects. While the primary causes of ASD are highly heterogeneous, it appears that they nevertheless converge on shared downstream epigenomic changes associated with specific functions. It is possible that these shared chromatin alterations could in turn drive some of the shared symptoms of ASD.