Multi-omics integration identifies key upstream regulators of pathomechanisms in hypertrophic cardiomyopathy due to truncating MYBPC3 mutations
Clinical Epigenetics volume 13, Article number: 61 (2021)
Hypertrophic cardiomyopathy (HCM) is the most common genetic disease of the cardiac muscle, frequently caused by mutations in MYBPC3. However, little is known about the upstream pathways and key regulators causing the disease. Therefore, we employed a multi-omics approach to study the pathomechanisms underlying HCM comparing patient hearts harboring MYBPC3 mutations to control hearts.
Using H3K27ac ChIP-seq and RNA-seq we obtained 9310 differentially acetylated regions and 2033 differentially expressed genes, respectively, between 13 HCM and 10 control hearts. We obtained 441 differentially expressed proteins between 11 HCM and 8 control hearts using proteomics. By integrating multi-omics datasets, we identified a set of DNA regions and genes that differentiate HCM from control hearts and 53 protein-coding genes as the major contributors. This comprehensive analysis consistently points toward altered extracellular matrix formation, muscle contraction, and metabolism. Therefore, we studied enriched transcription factor (TF) binding motifs and identified 9 motif-encoded TFs, including KLF15, ETV4, AR, CLOCK, ETS2, GATA5, MEIS1, RXRA, and ZFX. Selected candidates were examined in stem cell-derived cardiomyocytes with and without mutated MYBPC3. Furthermore, we observed an abundance of acetylation signals and transcripts derived from cardiomyocytes compared to non-myocyte populations.
By integrating histone acetylome, transcriptome, and proteome profiles, we identified major effector genes and protein networks that drive the pathological changes in HCM with mutated MYBPC3. Our work identifies 38 highly affected protein-coding genes as potential plasma HCM biomarkers and 9 TFs as potential upstream regulators of these pathomechanisms that may serve as possible therapeutic targets.
Hypertrophic cardiomyopathy (HCM), characterized by thickening of the myocardium that is not explained by abnormal loading conditions, is the most common inherited cardiac disease . More than 1,500 associated mutations, primarily in genes encoding parts of the sarcomere, such as MYBPC3, MYH7, and TNNT2, have been identified. The remodeled myocardium is characterized by cardiomyocyte hypertrophy and disarray, extensive fibrosis, and reduced capillary density . HCM is a heterogeneous disease as the onset, the disease phenotype, and the severity of the clinical presentations differ greatly among mutation carriers . Additionally, different mutated genes exhibit distinct biological impact on cellular contractility, energy metabolism, and sarcomeric protein expression in HCM hearts . However, the driving pathological mechanisms underlying the heterogenous HCM remain largely unknown.
Next-generation sequencing technologies (NGS) are instrumental in understanding disease etiology, delivering a clinical diagnosis, and discovering new treatment options. Multiple studies employed NGS to reveal the epigenetic modifications in heart failure, including DNA methylation and histone (de-)acetylation, which provided insights and identification of the driving mechanism underlying the disease [5, 6]. We previously showed the influence of histone acetylation changes on QRS complex-related GWAS loci in HCM [7, 8]. Additionally, studies from our group and others have shown that H3K27ac corresponds to the gene expression in human cardiac tissues as well as human cardiomyocytes [9, 10], highlighting it as a promising epigenetic mark for predicting gene expression. Studies also employed RNA sequencing to identify the affected transcription factor-mediated upstream regulatory events and the distinct gene expressions that define heart failure [9, 10]. Furthermore, data-piling studies are now connecting proteomics to NGS to get comprehensive information on the disease biology for precision medicine .
In this study, we aim to understand the pathomechanisms driving HCM by employing a multi-omics approach, including chromatin immunoprecipitation sequencing (ChIP-seq), RNA sequencing (RNA-seq), and proteomics, using myocardial tissue obtained from clinically well-phenotyped HCM patients with truncating MYBPC3 mutations and compared these with non-failing donor hearts. We revealed altered histone acetylome, transcriptome, and proteome profiles in HCM versus control hearts and studied affected biological functions. Besides, we identified key factors that may play a critical role in regulating the pathomechanisms underlying HCM. We also evaluated the contribution of histone acetylation and transcription signals in 11 cell types in the heart. Combined, this multi-omics study gives insight into the underlying disease pathways driving HCM and identifies promising candidates for therapeutic strategies.
Pairwise comparison between HCM and control hearts reveals distinct histone acetylome profiles
H3K27ac ChIP-seq was used to capture acetylated DNA regions in each sample and to compare the acetylation levels between 13 HCM (n = 13) and 10 control hearts (Fig. 1a and Additional file 9: Table S1) using DESeq2. In total, we identified 4,226 regions presented higher acetylation levels (hyperacetylated regions) and 5,084 regions presented lower acetylation levels (hypoacetylated regions) in HCM versus control hearts (Padj < 0.05, Fig. 1b, Additional file 10: Table S2A). Examples of hyper- and hypoacetylated regions are shown in Additional file 1: Figure S1. Region-to-gene annotation using either a 5 kb or 50 kb window from the transcription start site (TSS) revealed genes in the hyperacetylated regions are mostly involved in muscle contraction and extracellular matrix (ECM) development-related processes, whereas genes in the hypoacetylated regions are mostly involved in metabolic processes (Additional file 10: Table S2B-S2E).
Pairwise comparison between HCM and control hearts reveals distinct transcriptome profiles
Following RNA-seq of the same biopsies, the transcriptional profiles between 13 HCM (n = 13) and 10 control hearts were compared using DESeq2. In total, we identified 936 up-regulated genes and 1,097 down-regulated genes in HCM hearts compared to controls (Padj < 0.05, Fig. 1c, Additional file 11: Table S3A). The top biological processes enriched by up-regulated genes are involved in the muscle system process and energy production. The top biological processes enriched by down-regulated genes are involved in lipid metabolism and cell adhesion (Additional file 11: Table S3B and S3C).
Pairwise comparison between HCM and control hearts reveals distinct proteome profiles
We also performed proteomics in 11 HCM samples and compared their protein expression levels to another control group (n = 8) using DESeq2. In total, we identified 216 up-regulated proteins and 225 down-regulated proteins in HCM hearts compared to controls (P < 0.05, Fig. 1d, Additional file 12: Table S4A). The top enriched biological processes by up-regulated proteins are involved in muscular and ECM development. The top enriched biological processes by down-regulated proteins are involved in metabolism (Additional file 2: Figure S2, Additional file 12: Table S4B, and S4C).
Integrating histone acetylome, transcriptome, and proteome changes in HCM
We also integrated H3K27ac ChIP-seq and RNA-seq data in an unsupervised manner using O2PLS. Notably, without predefining the patient and the control groups, the first joint component of both ChIP-seq and RNA-seq plots discriminated HCM from control hearts (Fig. 1e, f). DNA regions and genes that contributed to the separation between HCM and control hearts are shown in Additional file 13: Table S5. Next, we aimed to identify key genes that underlie the separation between two groups and show altered expression levels in HCM hearts versus controls. Therefore, we overlapped 2,000 genes that discriminate HCM hearts and controls from the integrative analysis (O2PLS) with the differentially expressed genes from the pairwise comparison with the top and obtained 618 up-regulated genes and 308 down-regulated genes (Fig. 1g and Additional file 14: Table S6A). These overlapping genes are enriched for biological processes in the circulatory system and muscle contraction and pathways involved in the ECM formation and complement system (Fig. 1h, Additional file 14: Table S6B).
Since only 11 out of 13 HCM samples were included in the proteomics experiment and a different set of control samples was used in comparison with the HCM group, we could not apply O2PLS to integrate the proteomic data with either H3K27ac ChIP-seq or RNA-seq data. We, therefore, overlapped differentially expressed genes supported by DESeq2 and O2PLS analyses with proteomic data and identified 36 up-regulated protein-coding genes and 17 down-regulated protein-coding genes in HCM versus control hearts (Table 1). Notably, the protein levels of 38 protein-coding genes are detectable in the plasma and 5 out of them were consistently changed in the same direction in HCM versus control hearts at DNA (5 kb from TSS), RNA, and protein levels, including the up-regulation of ASPN, FMOD, MCAM, and NPPA and the down-regulation of AASS, highlighting them as promising candidates for biomarker discovery in HCM (Table 1).
Combined omics analyses identify ECM, muscle contraction, and metabolism as the main mechanisms altered in HCM
Genes annotated to differentially acetylated regions in the histone acetylome data, differentially expressed genes in the transcriptome data, and differentially expressed proteins in the proteome data consistently pointed to biological functions involved in enhancement of ECM, enhancement of muscle contraction, and suppression of metabolism in HCM versus control hearts. To further evaluate these biological functions, we collected gene sets that are known to regulate ECM remodeling, muscle contraction, and metabolism and performed GSEA to study whether they are primarily found in the up-regulated or down-regulated data sets in our study or they are randomly distributed. We observed that genes involved in ECM and cardiac muscle contraction were positively correlated with genes annotated to the hyperacetylated regions in the ChIP-seq data, the up-regulated genes in the RNA-seq data, and genes encoding the up-regulated proteins in the proteomics data (Fig. 2a, b). Genes that are related to fatty acid metabolism were significantly correlated with genes annotated to the hypoacetylated regions in the ChIP-seq data, the down-regulated genes in the RNA-seq data, and genes encoding the down-regulated proteins in the proteomics data (Fig. 2c, Additional file 3: Figure S3). Combined, we confirmed that pathways involved in ECM, muscle contraction, and metabolism were affected in HCM.
Consistently, protein-coding genes (36 up-regulated and 17 down-regulated ones) were enriched for comparable biological processes and pathways as shown in previous analyses, such as ECM and muscle contraction (Additional file 14: Table S6C). We also studied protein networks among them and observed that ECM and muscle contraction were the most enriched pathways by the up-regulated protein-coding genes and metabolism was the most enriched pathway by the down-regulated protein-coding genes (Fig. 3a, b).
Enriched transcription factor binding motifs in differentially acetylated regions
To define the actual factors that regulate altered RNA/protein expression in HCM hearts, we studied the putative transcription factors (TFs) encoded by enriched transcription factor binding motifs (TFBMs) in differentially acetylated regions by scanning through the HOCOMOCO motif database of 769 human primary and alternative binding models. We obtained 125 TFBMs in the hyperacetylated regions and 115 TFs were predicted to bind to them. In the hypoacetylated regions, 120 TFBMs were enriched and they encoded for 111 TFs. Out of those, 68 annotated TFs were enriched in both hyper- and hypo-acetylated regions (Additional file 4: Figure S4).
Notably, 2 TFs (KLF15 and ETV4) that were encoded by enriched TFBMs in the hyperacetylated regions showed higher mRNA levels in HCM hearts than controls, and 7 TFs (AR, CLOCK, ETS2, GATA5, MEIS1, RXRA, and ZFX) that were encoded by enriched TFBMs in the hypoacetylated regions showed lower mRNA levels in HCM hearts than controls (Fig. 3c). These TFs are enriched for biological functions involved in muscle hypertrophy and lipid metabolism (Additional file 15: Table S7), suggesting their potential roles in regulating the upstream signaling in HCM. Next, we examined their expressions in cardiomyocytes using the public single-nuclear RNA-seq data from adult human hearts . We observed that they were all expressed in ventricular cardiomyocytes (Additional file 7: Figure S7). Interestingly, KLF15, AR, CLOCK, MEIS1, and ZFX were more abundantly expressed than ETV4, ETS2, GATA5, and RXRA. We also examined their mRNA expressions in genome-engineered cardiomyocytes with mutated MYBPC3 compared to the control cardiomyocytes from a published study using RNA-seq (Additional file 18: Table S10), which provided insights into the early biological changes due to mutated MYBPC3 in cardiomyocytes . Interestingly, in line with our findings obtained from the bulk level at the severe stage, the mRNA level of ETV4 was significantly increased (P value = 9.56E−9) in mutant cardiomyocytes compared to the controls, whereas the mRNA level of RXRA was significantly decreased (P value = 0.001) in mutant versus control cardiomyocytes. However, they showed comparable expression levels between HCM and control hearts in the proteome data. We further examined their protein expressions using western blot and confirmed that they were not significantly changed in HCM versus control hearts (Fig. 3d). We performed additional immunofluorescence staining and examined the expression of KLF15 in human-induced pluripotent stem cell-derived cardiomyocytes (hiPSC-cardiomyocytes) with and without mutated MYBPC3. We observed the expression of KLF15 in hiPSC-cardiomyocytes (Fig. 6a); however, more cell lines derived from different patients and healthy individuals are needed to study whether the KLF15 signal differs between diseased and control cardiomyocytes.
A schematic overview of all analysis steps is shown in Fig. 3e.
Genes discriminating HCM from controls show consistent changes on various levels
We further selected an example of one down-regulated and one up-regulated gene and examined them in more detail to demonstrate the strength of our integrative omic analysis. ATP2A2, one protein-coding gene from the overlapping candidates, encodes the sarcoplasmic reticulum Ca2+ -ATPase pump SERCA2a and plays a critical role in the regulation of calcium handling . We identified ATP2A2 as one of the major candidates in discriminating HCM hearts from controls in our integrated H3K27ac ChIP-seq and RNA-seq analysis (Fig. 4a). Besides, its mRNA and protein levels were significantly lower in HCM versus control hearts in the pairwise comparison, and the suppressed protein level was also validated using western blot (Fig. 4b).
HSPA2, another protein-coding gene from the overlapping candidates, is involved in protein quality control . We demonstrated that HSPA2 acts as another key player in discriminating HCM hearts from controls (Fig. 4a). Its mRNA and protein levels were significantly higher in HCM versus control hearts in the pairwise comparison, and the enhanced protein level was also validated using western blot (Fig. 4c). The upstream regulatory region of HSPA2 also showed a higher acetylation level in HCM versus control hearts (Table 1). Moreover, we showed a profound HSPA2 intensity in HCM compared to the control heart using the immunohistochemistry staining at low magnification (Fig. 4d). Next, we employed immunofluorescence staining at high magnification and observed occasional HSPA2 aggregates in cardiomyocytes from the HCM heart, whereas no HSPA2 aggregates were shown in the control heart (Fig. 4e).
Allelic imbalance of MYBPC3 in HCM hearts is observed at both DNA and RNA levels
To further explore the potential of the produced data, we investigated the contribution of both MYBPC3 alleles in the sequencing datasets. In the patient cohort, three heterozygous truncating mutations were present in MYBPC3, namely c.2373dupG (n = 5), c.2827C > T (n = 6), and c.927-2A > G (n = 2). We observed that the average acetylation ratio of MYBPC3 with c.2373dupG, c.2827C > T, and c.927-2A > G mutation to wildtype allele was 50%, 25%, and 66.6%, respectively (Additional file 5: Figure S5A). The average mRNA expression ratio of MYBPC3 with c.2373dupG, c.2827C > T, and c.927-2A > G mutation to wildtype allele was 6.7%, 19.7%, and 43.4%, respectively (Additional file 5: Figure S5B). The acetylation and mRNA levels of three mutations were not observed in control hearts (Additional file 5: Figure S5C and S5D). It has to be noted that it is not possible to effectively distinguish between wildtype and mutant alleles in the proteomics data.
Cellular fraction sub-analysis of bulk cardiac tissue transcriptome indicates cardiomyocyte enrichment in a cellular-specific response in HCM
Since there are multiple cell types present in the heart samples, we collected cell-type-specific markers for cardiomyocytes and 11 non-myocyte cell types as revealed by recent single-cell studies in the heart [16, 17], ranging from 8 to 19 markers per cell type (Additional file 17: Table S9), and we examined their expression levels in our bulk sequencing data. Cardiomyocyte-specific markers showed higher histone acetylation and mRNA levels than markers of 11 non-myocyte cell types (Fig. 5a, b), regardless of health and disease. The mRNA expression of several cell-type-specific markers was significantly different between HCM and control hearts (Fig. 5c). The upstream regulatory region of three cardiomyocyte-specific markers and one T-cell-specific marker also showed significantly higher acetylation activities in HCM hearts versus controls (Fig. 5d).
Identified candidates supported by multi-omics data in human-induced pluripotent stem cells-derived cardiomyocytes (hiPSC-cardiomyocytes)
We further evaluated the presence of obtained candidates, which were supported by the acetylome, transcriptome, and/or proteome data (Table 1), using hiPSC-cardiomyocytes with and without mutated MYBPC3. We first employed the transcriptome profile comparing the genome-engineered cardiomyocytes harboring mutated MYBPC3 with the controls (Additional file 18: Table S10) from a published study . Out of 53 candidates, 18 protein-coding genes (34%) from our bulk data also showed significantly affected mRNA levels between diseased and control cardiomyocytes (Table 1). Next, we examined the expression levels of several candidates between hiPSC-cardiomyocytes with and without mutated MYBPC3 using immunofluorescence staining, including ACTN2 and ATP2A2. Consistent with the suppressed mRNA and protein levels of ACTN2 in HCM versus control hearts, we observed that unlike the well-aligned cytoskeleton structure in control cardiomyocytes, HCM cardiomyocytes seemed to have a disrupted cytoskeleton as revealed by ACTN2 staining (Fig. 6b). Besides, we again observed the inhibited expression of ATP2A2/SERCA2a in HCM cardiomyocytes when compared with the controls (Fig. 6c). Combined, we showed that candidates screened by multi-omics approaches from the cardiac tissues also expressed in cardiomyocytes and showed similar changing direction in HCM versus control cardiomyocytes.
In the present study, we studied changes in histone acetylome, transcriptome, and proteome between HCM and control hearts. Integrating these multi-omics data, we present for the first time a set of DNA regions and genes that differentiate HCM from control hearts and identified 53 protein-coding genes as the major contributors. Since these comprehensive datasets consistently indicated altered biological functions involved in ECM, muscle contraction, and metabolism, we further studied and identified the upstream TFs that could play a critical regulatory role. We also observed an abundance of acetylation signals and transcripts signals from cardiomyocyte-specific markers compared to markers of 11 non-myocyte populations in the heart.
By integrating ChIP-seq, RNA-seq, and proteomics data we identified a list of deciding elements in discriminating HCM hearts from controls. Among these, ATP2A2/SERCA2a, the cardiomyocyte-specific marker, transports cytosolic calcium into the sarcoplasmic reticulum and subsequently regulates the contraction and relaxation of cardiomyocytes . Suppressed ATP2A2/SERCA2a expression is associated with impaired relaxation in cardiomyocytes and contributes to diastolic dysfunction in HCM [14, 19]. Consistent with previous findings, we also showed lower mRNA and protein levels in HCM hearts compared to controls. Additionally, we employed hiPSC-derived cardiomyocytes and confirmed the suppressed ATP2A2/SERCA2a in HCM versus control cardiomyocytes.
Another element HSPA2 (heat shock-related 70 kDa protein 2), which protects cardiac integrity by correcting misfolded proteins upon stress , showed higher upstream acetylation level, mRNA level, and protein level in HCM hearts versus controls. Up-regulated HSPA2 protein level has been observed at both the initial stage of diseased cardiomyocytes and at the end stage of diseased heart with diastolic dysfunction [21, 22]. We previously also showed elevated HSPA2 protein level in HCM with mutations in different sarcomeric proteins when compared to controls, and it was negatively correlated with MYBPC3 peptide counts . Here we further demonstrated a cardiomyocyte-specific elevation of HSPA2 accompanied by the occasional HSPA2 aggregates in the HCM heart, which may contribute to the disease progression by hampering the protein quality control system.
Interestingly, both ATP2A2/Serca2a- and HSPA2-mediated activities are ATP-dependent, highlighting the strength of our sequencing data that align with previous observations showing energy deficiency as a key pathomechanism in HCM development. A reduced myocardial external efficiency, the ratio between cardiac work and oxygen consumption, was observed in asymptomatic MYBPC3 mutation carriers compared to the control hearts . We also provided a review summarizing the metabolic changes in HCM hearts . Isolated cardiomyocytes from HCM hearts with MYBPC3 mutations showed a reduced super-relax state, which may result in an increased sarcomere energy utilization . Likewise, hiPSC-derived cardiomyocytes with mutated MYH7, another secomeric gene that is often mutated in HCM patients, also showed a reduced super-relax state, an increased contractility, and a greater basal oxygen consumption rate when compared to the controls, suggesting a higher ATP demand . Taken together, previous studies show mutation-mediated contractile dysfunction with associated disturbed energetics. Here, besides the enriched metabolism using the ChIP-seq, RNA-seq, and proteomic data, we also identified several TFs as potential key upstream factors in regulating metabolism in HCM hearts, such as KLF15, AR, and RXRA [27,28,29]. Furthermore, we confirmed the expression of KLF15, a key regulator in cardiac metabolism , in hiPSC-cardiomyocytes using immunofluorescence staining, highlighting the potential involvement of KLF15 in diseased cardiomyocytes. Drugs targeting energy metabolism, such as trimetazidine, have been investigated in the clinical trial in HCM patients . Combined, identified TFs could be potential therapeutic targets in restoring metabolic homeostasis in HCM patients and future studies in hiPSC-cardiomyocyte models are warranted to define the mutation-mediated sequence in mitochondrial function and cellular metabolism.
Although those identified TFs showed the same changing direction at the acetylation and mRNA levels in HCM versus control hearts, their protein levels were comparable between HCM and control hearts. Similarly, some of those identified key protein-coding genes, which changed in the same direction at the mRNA and protein levels, showed comparable acetylation levels between HCM and control hearts. These findings indicate that the central dogma of cellular information transition from DNA, RNA, to protein is oversimplified . Numerous factors, such as the post-translational modifications on gene expression and TF activity, mRNA degradation, and protein degradation [33,34,35,36], could lead to a poor correlation between DNA, mRNA, and protein expressions. Additionally, both H3K27ac ChIP-seq and RNA-seq data were generated from the same 13 HCM and 10 control samples, whereas only 11 out of 13 HCM samples were included in the proteomics analysis in comparison with 8 different control samples. Thus, different donor samples and group sizes might obscure the integrative analysis of the multi-omics data. Therefore, it is important to find an optimal way to integrate multi-omics data rather than a simple correlative. The discovery of new types of machinery beyond the central dogma is also needed in future studies.
Previous studies from our group and others showed around a 30% ratio of mutant to wildtype MYH7- and MYBPC3-alleles in HCM patients at the mRNA and the protein levels, respectively [37, 38]. In this study, we present the imbalanced mRNA expression of three truncating mutations in MYBPC3. Notably, we showed for the first time the presence of such an imbalance at the DNA level by evaluating allele-specific histone acetylation signals. The mutant/wildtype expression ratio differs between truncating mutations and between DNA and mRNA levels, suggesting mutation-specific imbalance and additional machinery between DNA and mRNA. It is important to point out that H3K27ac ChIP-seq is sensitive in capturing upstream regulatory regions of a gene and is restored in covering the core gene body , and the allelic expression based on the histone acetylation level could be limited by the nature of the technique. Nevertheless, this finding highlights another utilization of omics data in detecting allelic imbalance for genetic diseases.
As the cell composition changes in diseased hearts, we investigated cell-type-specific signals from cardiomyocytes and 11 non-myocyte cell types and observed more ChIP-seq signals and transcripts in cardiomyocyte-specific markers than all non-myocyte markers, suggesting a considerable portion of cardiomyocytes-derived signals in the bulk sequencing data. Gene set enrichment analysis showed enriched pathways in ECM, which is in line with previous reports indicating fibrosis as a hallmark of HCM . Apart from fibroblasts, cardiomyocytes also express ECM-related genes, including two fibroblast-specific markers (DKK3 and ABI3BP) . Because cell-type-specific markers might also be expressed by other cell types and the limitation of capturing pure cell populations from snap-frozen tissues without damaging the DNA and RNA integrity, we can only speculate that the majority of the signals were likely obtained from cardiomyocytes and it remains a challenge to separate them in bulk sequencing. Additionally, we used cell-type-specific markers obtained from the murine models [16, 17], which may not translate one-to-one across species. Therefore, future studies are required to validate the cellular-specific response in HCM by isolating single-cell populations from human cardiac tissues.
In conclusion, our study presents detailed information in HCM hearts with truncating MYBPC3 mutations. These data showed altered ECM, muscle contraction, and metabolism in HCM. Integrative analyses further identified a subset of protein-coding genes and upstream TFs that could drive these pathophysiological mechanisms and serve as promising diagnostic and/or therapeutic targets. We also showed a considerable amount of cardiomyocyte-derived signals compared to non-myocyte cell types, providing cardiomyocyte-specific insights to better understand HCM in future studies.
Materials and methods
Human cardiac samples
The study protocol was approved by the local medical ethics review committees, including the Biobank Research Ethics Committee of University Medical Center Utrecht (protocol number 12/387), the Local Ethics Committee of the Erasmus MC (2010-409), the Washington University School of Medicine Ethics Committee (Institutional Review Board), and the Sydney Heart Bank (HREC Univ. Sydney 2012/030). Since septal thickening is a typical feature of HCM, septal myectomy is commonly performed to relieve left ventricular outflow tract obstruction in HCM patients [2, 42, 43]. Therefore, we collected cardiac tissue of the interventricular septum of 13 HCM patients. Genetic analyses of all patients revealed three truncating pathogenic heterozygous mutations in MYBPC3, namely c.2373dupG in 5 HCM patients, c.2827C > T in 6 patients, and c927-2A > G in 2 HCM patients. Control tissues of 18 non-failing donors were obtained from the Biobank of University Medical Center Utrecht, the Washington University School of Medicine, and the Sydney Heart Bank. Informed consent was obtained from each patient prior to surgery or was waived by the ethics committee when acquiring informed consent was not possible due to the death of the donor. Samples were collected and snap-frozen in liquid nitrogen and stored at − 80 ˚C up until analysis. Detailed clinical characteristics are shown in Additional file 9: Table S1, and an overview of samples included in the following experiments is shown in Fig. 1a.
H3K27ac chromatin immunoprecipitation and sequencing (ChIP-seq)
We performed chromatin immunoprecipitation and sequencing (ChIP-seq) using the H3K27ac mark to study the differences of the histone acetylome between patient and control samples as previously described . To study the differences of the histone acetylome between patient and control samples, we performed chromatin immunoprecipitation and sequencing (ChIP-seq) using the H3K27ac mark. Briefly, all cardiac samples were sectioned at a thickness of 10 µm, and chromatin was isolated using the MAGnify™ Chromatin Immunoprecipitation System kit (Life Technologies) according to the manufacturer’s instructions. The anti-histone H3K27ac antibody (ab4729, Abcam) was used for immunoprecipitation. Captured DNA was purified using the ChIP DNA Clean & Concentrator kit (Zymo Research). Libraries were prepared using the NEXTflex™ Rapid DNA Sequencing Kit (Bioo Scientific). Samples were PCR amplified, checked for the proper size range and absence of adaptor dimers on a 2% agarose gel, and barcoded libraries were sequenced 75-bp single end on an Illumina NextSeq500 sequencer. Sequencing reads were mapped against the reference genome (hg19 assembly, NCBI37) using the BWA package (mem –t 7 –c 100 –M –R)1. Multiple reads mapping to the same location and strand were collapsed to a single read, and only uniquely placed reads were used for peak calling. Peaks/regions were called using Cisgenome 2.02 (–e 150 -maxgap 200 –minlen 200). Region coordinates from all samples were stretched to at least 2000 base pairs and collapsed into a single common list. Overlapping regions were merged based on their coordinates. Only regions supported by at least 2 independent datasets were further analyzed. Autosomal sequencing reads from each ChIP-seq library were overlapped back with the common region list to set the H3K27ac occupancy for every region-sample pair. Differentially acetylated regions between HCM and control hearts were identified using DESeq2 under the default setting in the Galaxy environment .
Annotating genes in the vicinity of differentially acetylated regions
Region-to-gene annotation was performed to study potentially affected genes in the vicinity of DNA regions with altered acetylation levels in HCM hearts when compared with controls. Differentially acetylated regions located within either ± 5 kb or ± 50 kb window from the transcription start site (TSS) of all genes were obtained, and the nearest genes of these regions were collected.
Predicting transcription factor binding motifs in differentially acetylated regions
To study the putative upstream signaling, we studied the enriched transcription factor binding motifs (TFBMs) by differentially acetylated regions and motifs-encoded transcription factors (TFs). DNAse I hypersensitivity regions in human cardiac samples, which play a key role in transcription factor footprinting , were collected from the ENCODE project and overlapped with differentially acetylated regions in this study . Overlapping DNA sequences between differentially acetylated regions and DNAse I hypersensitivity regions were used to studying the enriched transcription factor binding motifs and motifs-encoded TFs using MEME Suite AME tool (HOCOMOCO Human v11 Full, average odds scoring method, and Fisher’s exact test) .
We also performed RNA sequencing (RNA-seq) and obtained the transcriptome landscapes in all samples. Briefly, RNA was isolated using the RNeasy Micro Kit (Qiagen) or ISOLATE II RNA Mini Kit (Bioline) according to the manufacturer's instructions. Sample quality was assessed using the 2100 Bioanalyzer with an RNA 6000 Pico Kit (Agilent), and sample quantity was measured using Qubit Fluorometer with an HS RNA Assay (Thermo Fisher). Afterward, libraries were prepared using the NEXTflex™ Rapid RNA-seq Kit (Bioo Scientific) and sequenced by the Nextseq500 platform (Illumina). Sequenced reads were aligned to the human reference genome GRCh37 using STAR v2.4.2a . Reads per kilobase million reads sequenced (RPKMs) were calculated with edgeR’s RPKM function . To identify a list of differentially expressed genes between HCM and control hearts at Padj < 0.05, we employed d DESeq2 to process all the raw counts per sample per group in the Galaxy environment .
Allele-specific expression of MYBPC3 in ChIP-seq and RNA-seq data
RNA-seq and ChIP-seq reads were processed with TrimGalore (version 0.6.5) to detect and clip off adaptor sequences, reads with a remaining length of 20 or more nucleotides, and an average read quality q > 20 were selected. The reads were aligned to the human genome (GRCh37) with transcript annotation from Ensembl (version 74) using the STAR aligner (version 2.7.1a). To facilitate the equal alignment of reads from wildtype and mutant alleles, alignments were made without clipping off end sequences (STAR option alignEndsType set to EndToEnd), and only best scoring alignments were selected (STAR option outSAMprimaryFlag set to AllBestScore). Alignments were selected to remove duplicated reads by an in-house HTSeq based python (version 2.7.10) script . Reads that aligned to the same genomic interval and that aligned with identical bases to the two non-indel MYBPC3 SNPs (rs397516082 and rs387907267) of interest, were considered duplicates. ChIP-seq reads matching to one of the three SNPs of interest were assigned to the wildtype or mutant allele based on the following rules: for the single-base polymorphic SNPs (rs397516082 and rs387907267), based on having the wildtype or mutant nucleotide aligned to the SNP position; for the indel SNP (rs397515963) based on having the exact wildtype or mutant sequence of the SNP with 10 surrounding bases in a read aligned to the SNP. In addition, ChIP-seq reads with wildtype or mutant assignment were required to have > 90% of the bases aligned to the genome. For RNA-seq reads, the assignment to be wildtype or mutant allele derived was complicated by changed splice patterns for two of the three SNPs. For rs387907267 (c.2827C > T) that has no changed splice patterns, reads were classified based on having the wildtype or mutant nucleotide aligned to the SNP position. For rs397515963 (c.2373dupG) that disrupts the correct splicing of intron 23, reads covering the splice junction with the correct intron spliced out were counted as wildtype reads, whereas reads using either the correct intron donor or acceptor site, but not both, were counted as mutant reads. For rs397516082 (c.927-2A > G) that disrupts the correct splicing of intron 11, reads covering the splice junction with the correct intron spliced out were counted as wildtype reads, whereas reads that were running through the intron 11–exon 12 splice site and into the SNP position were counted as wildtype or mutant based on having the wildtype or mutant nucleotide aligned to the SNP position.
Integrating ChIP-seq and RNA-seq data with Two-way Orthogonal Partial Least Squares
To find common parts between RNA-seq and ChIP-seq data simultaneously across all genes and regions, a data integration approach is considered using two-way orthogonal partial least squares (O2PLS) . O2PLS decomposes both RNA-seq and ChIP-seq datasets into joint, omic-specific, and residual parts. The joint subspaces contain variations that are correlated to one another. The joint principal components (JPCs) that span the joint subspaces are obtained by finding linear combinations of genes and regions that maximize the covariation. Omic-specific subspaces capture variation unrelated to another omics dataset, enabling JPCs to better estimate the underlying system. Here, the rows of the ChIP-seq and RNA-seq data should represent the same samples. Note that O2PLS uses all genes and regions in the datasets and does not rely on prior information about the position or function of these features. Furthermore, O2PLS is unsupervised, and its algorithm is implemented in the OmicsPLS R package and freely available from CRAN . Prior to the analysis, genes with expression lower than 10 counts in at least 22 samples were removed. Samples in both datasets were matched, and 23 overlapping samples were retained. The expression data are normalized. Both datasets were log-transformed and quantile normalized across samples. The dimensionality of the preprocessed datasets is 23 by 15,882 (RNA-seq) and 23 by 33,642 (ChIP-seq).
We performed proteomics using cardiac samples from the same patients (n = 11) and 8 non-failing donor samples from the Sydney Heart Bank (Additional file 9: Table S1). Briefly, proteins were loaded to a 4–12% NuPAGE Novex Bis–Tris 1.5 mm mini gel (Invitrogen) for separation, followed by fixation (50% ethanol and 3% phosphoric acid) and staining using 0.1% Coomassie brilliant blue G-250 solution. In-gel digestion was performed, and the samples were concentrated in a vacuum centrifuge as described previously . Nano-LC–MS/MS was performed as described previously . Briefly, separated peptides were separated using an Ultimate 3000 Nano LC–MS/MS system (Dionex LC-Packings, Amsterdam, The Netherlands) and trapped. Eluting peptides were ionized into a Q Exactive mass spectrometer (Thermo Fisher, Bremen, Germany), and intact masses were measured in the orbitrap. Among all, the top 10 peptide signals were selected and analyzed using the MS/MS. MS/MS spectra were searched against a Uniprot human reference proteome FASTA file (Swissprot_2017_03_human_canonical_and_isoform.fasta, 42,161 entries) using MaxQuant version 18.104.22.168. Differentially expressed proteins between HCM and control hearts at P < 0.05 were identified using DESeq2 .
Functional enrichment analysis
Gene ontology (GO) enrichment analysis: To study the enriched biological functions, GO enrichment analysis was performed using the ToppFun ToppGene Suite under the default settings (FDR correction, P value cutoff of 0.05, and gene limits between 1 and 2,000) .
Gene set enrichment analysis (GSEA): Established gene sets involved in the most enriched biological functions in HCM versus control hearts were collected from Molecular Signature Database v7.1 (Additional file 16: Table S8) [56, 57]. Each gene set per biological function was studied for its positive or negative correlation with genes annotated from the altered acetylated levels in the ChIP-seq data, genes with altered mRNA levels in the RNA-seq data, and protein-coding genes with altered protein levels in the proteomic data in this study under the following setting: Number of permutations: 1,000; Phenotype labels: up_versus_down; Collapses dataset to gene symbols: false; Permutation type: gene_set; Enrichment statistics: weighted; Metric for ranking genes: Signal2Noise; Gene list sorting mode: real; Gene list ordering mode: descending; Max size: 500; Min size: 1.
Protein–protein interaction (PPI) networks: Protein networks were performed using the STRING Version 11.0 under the following settings: the meaning of network edges: confident, minimum required interaction score: medium confidence (0.400) .
Human-induced pluripotent stem cells-derived cardiomyocytes (hiPSC-cardiomyocytes)
Clones from one control hiPSC line and one hiPSC line derived from an HCM patient with mutated MYBPC3 were differentiated to cardiomyocytes as previously described . The contractile function of differentiated control and HCM cardiomyocytes was examined as previously described . Cells were seeded at the density of 20,000 cells per well to the 96-well plate. After culturing for 27 days, cardiomyocytes were fixed and stained for interested proteins.
Western blot was performed as described previously . Primary antibodies, including mouse anti-HSPA2 (heat shock-related 70 kDa protein 2, 66291-1, Proteintech Group), mouse anti-GAPDH (10R-G109a, Fitzgerald Industries International), anti-ATP2A2 (sarcoplasmic reticulum Ca2 + -ATPase 2, also known as SERCA2), anti-KLF15 (AV32587, Sigma-Aldrich), anti-RXRA (ab125001, Abcam), and anti-β-actin, were used.
Immunohistochemistry and immunofluorescence (IF) staining
Tissue sections were thawed and left at RT for 20 min inside a closed box. Then the sections were treated with peroxidase blocking solution (3% H2O2 in MeOH) and blocked with 1% BSA. HSPA2 (Anti-HSPA2 rabbit-Polyclonal antibody Prestige Antibodies HPA000798) primary antibodies were incubated for 1 h at RT (1:100 for IF and 1:200 for IHC). Sections were washed and incubated with the secondary antibodies for 30 min. Vector Vectastain Universal Elite ABC Kit (PK-6200) was used to enhance the brightfield staining (IHC) and WGA (for cell membranes) and DAPI (nuclei) were added to counterstain fluorescence staining. The brightfield slides were washed, treated with DAB, counterstained with Mayers Hematoxylin, dehydrated and mounted with DPX. Fluorescence staining was mounted with Mowiol. Images were acquired with the Vectra Polaris Scanner (brightfield/IHC) or Confocal Nikon A1 (IF). The analysis was performed with QuPath, Fiji, and NIS Nikon software.
HiPSC-cardiomyocytes were fixed using 4% paraformaldehyde solution for 10 min at room temperature and incubated overnight with the primary antibodies, including ACTN2 (A7811, Sigma-Aldrich) and ATP2A2 (MA3-910, Thermo Fisher Scientific). Afterward, cells were washed three times with PBS and incubated for 1 h with Hoechst 33,342 and the secondary antibodies (Invitrogen), including Alexa Fluor 488 and Alexa Fluor 568. Images were taken by Leica confocal microscope at 20× magnification.
Availability of data and materials
All relevant data are available within the article and the supplementary files. Because of the sensitive nature of the data collected for this study (13 patient samples and 10 control samples), requests to access the raw sequencing dataset from qualified researchers trained in human subject confidentiality protocols may be sent to the corresponding authors. It is important to note that processed RNA-seq and H3K27ac ChIP-seq data of 11 patient samples and 7 control samples are published in a GWAS study  and available in Additional file 11: Table S3 and Additional file 13: Table S5, respectively. Raw proteomics data within the article can be found at the ProteomeXchange Consortium via the PRIDE partner repository with the dataset identifier PXD012467 .
chromatin immunoprecipitation and sequencing
differentially acetylated regions
transcription start site
transcription factor binding motifs
reads per kilobase million reads sequenced
two-way orthogonal partial least squares
joint principal components
gene set enrichment analysis
Elliott PM, Anastasakis A, Borger MA, Borggrefe M, Cecchi F, Charron P, et al. ESC guidelines on diagnosis and management of hypertrophic cardiomyopathy. The task force for the diagnosis and management of hypertrophic cardiomyopathy of the European Society of Cardiology (ESC). Eur Heart J Oxf Acad. 2014;35:2733–79.
Llam N, Bollen IAE, Niessen HWM, Dos Remedios CG, Michels M, Poggesi C, et al. Sex-specific cardiac remodeling in early and advanced stages of hypertrophic cardiomyopathy. PLoS ONE. 2020;15:e0232427–e0232427.
Sabater-Molina M, Pérez-Sánchez I, Hernández Del Rincón JP, Gimeno JR. Genetics of hypertrophic cardiomyopathy: a review of current state. Clin Genet. 2018;93:3–14.
Wijnker PJM, van der Velden J. Mutation-specific pathology and treatment of hypertrophic cardiomyopathy in patients, mouse models and human engineered heart tissue. Biochim Biophys Acta Mol Basis Dis. 2020;165774.
Liu C-F, Tang WHW. Epigenetics in Cardiac Hypertrophy and Heart Failure. JACC Basic Transl Sci. 2019;4:976–93.
Zhang W, Song M, Qu J, Liu G-H. Epigenetic modifications in cardiovascular aging and diseases. Circ Res. 2018;123:773–86.
Hemerich D, Pei J, Harakalova M, van Setten J, Boymans S, Boukens BJ, et al. Integrative functional annotation of 52 genetic loci influencing myocardial mass identifies candidate regulatory variants and target genes. Circ Genom Precis Med. 2019;12:e002328.
Manduchi E, Hemerich D, Setten J, Tragante V, Harakalova M, Pei J, et al. A comparison of two workflows for regulome and transcriptome-based prioritization of genetic variants associated with myocardial mass [Internet]. Genet Epidemiol. 2019. https://doi.org/10.1002/gepi.22215.
Gilsbach R, Schwaderer M, Preissl S, Grüning BA, Kranzhöfer D, Schneider P, et al. Distinct epigenetic programs regulate cardiac myocyte development and disease in the human heart in vivo. Nat Commun. 2018;9:391.
Pei J, Harakalova M, Treibel TA, Thomas Lumbers R, Boukens BJ, Efimov IR, et al. H3K27ac acetylome signatures reveal the epigenomic reorganization in remodeled non-failing human hearts. Clin Epigenetics BioMed Central. 2020;12:1–18.
Low TY, Mohtar MA, Ang MY, Jamal R. Connecting proteomics to next-generation sequencing: proteogenomics and its current applications in biology. Proteomics. 2019;19:e1800235.
Litviňuková M, Talavera-López C, Maatz H, Reichart D, Worth CL, Lindberg EL, et al. Cells of the adult human heart. Berlin: Nature Publishing Group; 2020. p. 1–7.
Helms AS, Tang VT, O’Leary TS, Friedline S, Wauchope M, Arora A, et al. Effects of MYBPC3 loss-of-function mutations preceding hypertrophic cardiomyopathy. JCI Insight [Internet]. JCI Insight; 2020 [cited 2020 Dec 28];5. Available from: https://pubmed.ncbi.nlm.nih.gov/31877118/
Kresin N, Stücker S, Krämer E, Flenner F, Mearini G, Münch J, et al. Analysis of contractile function of permeabilized human hypertrophic cardiomyopathy multicellular heart tissue. Front Physiol. 2019;10:239.
Dorsch LM, Schuldt M, dos Remedios CG, Schinkel AFL, de Jong PL, Michels M, et al. Protein quality control activation and microtubule remodeling in hypertrophic cardiomyopathy. Cells. 2019. https://doi.org/10.3390/cells8070741.
Skelly DA, Squiers GT, McLellan MA, Bolisetty MT, Robson P, Rosenthal NA, et al. Single-Cell Transcriptional Profiling Reveals Cellular Diversity and Intercommunication in the Mouse Heart. Cell Rep. 2018;22:600–10.
Gladka MM, Molenaar B, de Ruiter H, van der Elst S, Tsui H, Versteeg D, et al. Single-cell sequencing of the healthy and diseased heart reveals cytoskeleton-associated protein 4 as a new modulator of fibroblasts activation. Circulation. 2018;138:166–80.
New perspectives on the role of SERCA2’s Ca2+ affinity in cardiac function. Biochimica et Biophysica Acta (BBA) - Molecular Cell Research. Elsevier; 2006;1763:1216–28.
Somura F, Izawa H, Iwase M, Takeichi Y, Ishiki R, Nishizawa T, et al. Reduced myocardial sarcoplasmic reticulum Ca(2+)-ATPase mRNA expression and biphasic force-frequency relations in patients with hypertrophic cardiomyopathy. Circulation; 2001 [cited 2020 Jul 10];104. Available from: https://pubmed.ncbi.nlm.nih.gov/11489771/
Ranek MJ, Stachowski MJ, Kirk JA, Willis MS. The role of heat shock proteins and co-chaperones in heart failure. Philos Trans R Soc Lond B Biol Sci [Internet]. The Royal Society; 2018 [cited 2020 Apr 9];373. Available from: https://www.ncbi.nlm.nih.gov/pmc/articles/PMC5717530/
Birket MJ, Raibaud S, Lettieri M, Adamson AD, Letang V, Cervello P, et al. A human stem cell model of fabry disease implicates LIMP-2 accumulation in cardiomyocyte pathology. Stem Cell Rep. 2019. https://doi.org/10.1016/j.stemcr.2019.07.004.
Li W, Rong R, Zhao S, Zhu X, Zhang K, Xiong X, et al. Proteomic analysis of metabolic, cytoskeletal and stress response proteins in human heart failure. J Cell Mol Med; 2012 [cited 2020 Jul 10];16. Available from: https://pubmed.ncbi.nlm.nih.gov/21545686/
Güçlü A, Knaapen P, Harms HJ, Parbhudayal RY, Michels M, Lammertsma AA, et al. Disease stage-dependent changes in cardiac contractile performance and oxygen utilization underlie reduced myocardial efficiency in human inherited hypertrophic cardiomyopathy. Circ Cardiovasc Imaging. 2017;10. Available from: http://dx.doi.org/https://doi.org/10.1161/CIRCIMAGING.116.005604
van der Velden J, Tocchetti CG, Varricchi G, Bianco A, Sequeira V, Hilfiker-Kleiner D, et al. Metabolic changes in hypertrophic cardiomyopathies: scientific update from the Working Group of Myocardial Function of the European Society of Cardiology. Cardiovasc Res. 2018;114:1273.
McNamara JW, Li A, Lal S, Martijn Bos J, Harris SP, van der Velden J, et al. MYBPC3 mutations are associated with a reduced super-relaxed state in patients with hypertrophic cardiomyopathy. PLoS ONE. 2017;12:e0180064.
Toepfer CN, Garfinkel AC, Venturini G, Wakimoto H, Repetti G, Alamo L, et al. Myosin sequestration regulates sarcomere function, cardiomyocyte energetics, and metabolism, informing the pathogenesis of hypertrophic cardiomyopathy. Circulation. 2020;141:828–42.
Gonthier K, Poluri RTK, Audet-Walsh É. Functional genomic studies reveal the androgen receptor as a master regulator of cellular energy metabolism in prostate cancer. J Steroid Biochem Mol Biol; 2019 [cited 2020 Jul 5];191. Available from: https://pubmed.ncbi.nlm.nih.gov/31051242/
Subbarayan V, Mark M, Messadeq N, Rustin P, Chambon P, Kastner P. RXRα overexpression in cardiomyocytes causes dilated cardiomyopathy but fails to rescue myocardial hypoplasia in RXRα-null fetuses. J Clin Investig. 2000. https://doi.org/10.1172/jci8150.
Fisch S, Gray S, Heymans S, Haldar SM, Wang B, Pfister O, et al. Kruppel-like factor 15 is a regulator of cardiomyocyte hypertrophy. Proc Natl Acad Sci USA. 2007;104:7074–9.
Prosdocimo DA, Anand P, Liao X, Zhu H, Shelkay S, Artero-Calderon P, et al. Kruppel-like factor 15 is a critical regulator of cardiac lipid metabolism. J Biol Chem. 2014;289:5914–24.
van Driel BO, van Rossum AC, Michels M, Huurman R, van der Velden J. Extra energy for hearts with a genetic defect: ENERGY trial. Neth Heart J. 2019;27:200.
Franklin S, Vondriska TM. Genomes, proteomes, and the central dogma. Circ Cardiovasc Genet. 2011;4:576.
Valencia-Sanchez MA, Liu J, Hannon GJ, Parker R. Control of translation and mRNA degradation by miRNAs and siRNAs. Genes Dev. 2006;20:515–24.
Sha Z, Zhao J, Goldberg AL. Measuring the overall rate of protein breakdown in cells and the contributions of the ubiquitin-proteasome and autophagy-lysosomal pathways. Methods Mol Biol. 2018;1844:261–76.
Teixeira D, Sheth U, Valencia-Sanchez MA, Brengues M, Parker R. Processing bodies require RNA for assembly and contain nontranslating mRNAs. RNA. 2005;11:371–82.
Zhang Y, Li S, Xiang Y, Qiu L, Zhao H, Hayes JD. The selective post-translational processing of transcription factor Nrf1 yields distinct isoforms that dictate its ability to differentially regulate gene expression. Sci Rep. 2015;5:1–30.
Montag J, Kowalski K, Makul M, Ernstberger P, Radocaj A, Beck J, et al. Burst-like transcription of mutant and wildtype -alleles as possible origin of cell-to-cell contractile imbalance in hypertrophic cardiomyopathy. Front Physiol. 2018;9:359.
Parbhudayal RY, Garra AR, Götte MJW, Michels M, Pei J, Harakalova M, et al. Variable cardiac myosin binding protein-C expression in the myofilaments due to MYBPC3 mutations in hypertrophic cardiomyopathy. J Mol Cell Cardiol. 2018;123:59–63.
Sanchez GJ, Richmond PA, Bunker EN, Karman SS, Azofeifa J, Garnett AT, et al. Genome-wide dose-dependent inhibition of histone deacetylases studies reveal their roles in enhancer remodeling and suppression of oncogenic super-enhancers. Nucleic Acids Res. 2018;46:1756.
Freitas P, Ferreira AM, Arteaga-Fernández E, de Oliveira AM, Mesquita J, Abecasis J, et al. The amount of late gadolinium enhancement outperforms current guideline-recommended criteria in the identification of patients with hypertrophic cardiomyopathy at risk of sudden cardiac death. J Cardiovasc Magn Reson. 2019;21:1–10.
Heras-Bautista CO, Mikhael N, Lam J, Shinde V, Katsen-Globa A, Dieluweit S, et al. Cardiomyocytes facing fibrotic conditions re-express extracellular matrix transcripts. Acta Biomater. 2019;89:180–92. https://doi.org/10.1016/j.actbio.2019.03.017.
Reant P, Captur G, Mirabel M, Nasis A, Sado MD, Maestrini V, et al. Abnormal septal convexity into the left ventricle occurs in subclinical hypertrophic cardiomyopathy. J Cardiovasc Magn Reson. 2015;17:64.
Dearani JA, Danielson GK. Septal myectomy for obstructive hypertrophic cardiomyopathy. Semin Thorac Cardiovasc Surg Pediatr Card Surg Annu; 2005 [cited 2020 Dec 14]; Available from: https://pubmed.ncbi.nlm.nih.gov/15818363/
Afgan E, Baker D, Batut B, van den Beek M, Bouvier D, Cech M, et al. The Galaxy platform for accessible, reproducible and collaborative biomedical analyses: 2018 update. Nucleic Acids Res. 2018;46:W537–44.
Ouyang N, Boyle AP. TRACE: transcription factor footprinting using DNase 1 hypersensitivity data and DNA sequence. bioRxiv. 2019. https://doi.org/10.1101/801001v1.abstract.
Rosenbloom KR, Sloan CA, Malladi VS, Dreszer TR, Learned K, Kirkup VM, et al. ENCODE data in the UCSC genome browser: year 5 update. Nucleic Acids Res. 2013;41:D56-63.
McLeay RC, Bailey TL. Motif Enrichment Analysis: a unified framework and an evaluation on ChIP data. BMC Bioinform. 2010. https://doi.org/10.1186/1471-2105-11-165.
Dobin A, Davis CA, Schlesinger F, Drenkow J, Zaleski C, Jha S, et al. STAR: ultrafast universal RNA-seq aligner. Bioinformatics. 2013;29:15–21.
Robinson MD, McCarthy DJ, Smyth GK. edgeR: a Bioconductor package for differential expression analysis of digital gene expression data. Bioinformatics. 2010;26:139–40.
Anders S, Pyl PT, Huber W. HTSeq–a Python framework to work with high-throughput sequencing data. Bioinformatics. 2015;31:166–9.
el Bouhaddani S, el Bouhaddani S, Houwing-Duistermaat J, Salo P, Perola M, Jongbloed G, et al. Evaluation of O2PLS in Omics data integration. BMC Bioinform. 2016. https://doi.org/10.1186/s12859-015-0854-z.
el Bouhaddani S, el Bouhaddani S, Uh H-W, Jongbloed G, Hayward C, Klarić L, et al. Integrating omics datasets with the OmicsPLS package. BMC Bioinform. 2018. https://doi.org/10.1186/s12859-018-2371-3.
Warmoes M, Jaspers JE, Pham TV, Piersma SR, Oudgenoeg G, Massink MP, et al. Proteomics of mouse BRCA1-deficient mammary tumors identifies DNA repair proteins with potential diagnostic and prognostic value in human breast cancer. Mol Cell Proteomics; 2012 [cited 2020 Jul 10];11. Available from: https://pubmed.ncbi.nlm.nih.gov/22366898/
Piersma SR, Broxterman HJ, Kapci M, de Haas RR, Hoekman K, Verheul HM, et al. Proteomics of the TRAP-induced Platelet Releasate. J Proteomics [Internet]. J Proteomics; 2009 [cited 2020 Jul 10];72. Available from: https://pubmed.ncbi.nlm.nih.gov/19049909/
Chen J, Bardes EE, Aronow BJ, Jegga AG. ToppGene Suite for gene list enrichment analysis and candidate gene prioritization. Nucleic Acids Res. 2009;37:W305–11.
Subramanian A, Tamayo P, Mootha VK, Mukherjee S, Ebert BL, Gillette MA, et al. Gene set enrichment analysis: A knowledge-based approach for interpreting genome-wide expression profiles. Proc Natl Acad Sci USA. 2005;102:15545–50.
Liberzon A, Birger C, Thorvaldsdóttir H, Ghandi M, Mesirov JP, Tamayo P. The Molecular Signatures Database (MSigDB) hallmark gene set collection. Cell systems. 2015;1:417.
Snel B. STRING: a web-server to retrieve and display the repeatedly occurring neighbourhood of a gene. Nucleic Acids Res. 2000. https://doi.org/10.1093/nar/28.18.3442.
van den Brink L, Brandão KO, Yiangou L, Mol MPH, Grandela C, Mummery CL, et al. Cryopreservation of human pluripotent stem cell-derived cardiomyocytes is not detrimental to their molecular and functional properties. Stem Cell Res. 2020;43:101698. https://doi.org/10.1016/j.scr.2019.101698.
Birket MJ, Ribeiro MC, Kosmidis G, Ward D, Leitoguinho AR, van de Pol V, et al. Contractile defect caused by mutation in MYBPC3 revealed under conditions optimized for human PSC-cardiomyocyte function. Cell Rep. 2015;13:733–45.
Perez-Riverol Y, Csordas A, Bai J, Bernal-Llinares M, Hewapathirana S, Kundu DJ, et al. The PRIDE database and related tools and resources in 2019: improving support for quantification data. Nucleic Acids Res. 2019;47:D442–50.
We would like to thank Pedro Espinosa for performing immunohistochemistry and immunofluorescent stainings. Grateful thanks to Cris dos Remedios and the Sydney Heart Bank for providing non-failing donor tissue.
This work was supported by the Netherlands Foundation for Cardiovascular Excellence (to C.C.), the NWO VENI Grant (No. 016.176.136 to M.H.), three NWO VIDI Grants (No. 91714302 to C.C., No. 016096359 to M.C.V and No. 91715303 to R.P.D.), ZonMW-NWO VICI Grant 91818902 (to J.V.), the Erasmus MC fellowship grant (to C.C.), the RM fellowship grant of the UMC Utrecht (to C.C.), Wilhelmina Children’s Hospital research funding (No. OZF/14 to M.H.), the Netherlands Cardiovascular Research Initiative: An initiative with the support of the Dutch Heart Foundation (CVON2014-40 DOSIS to J.V., M.H., F.W.A., CVON2014-11 RECONNECT to C.C., M.C.V., Dekker 2015T041 to A.F.B. and Queen of Heart to C.C. and M.C.V.), UCL Hospitals NIHR Biomedical Research Centre (to F.W.A.), and the Starting Grant (STEMCARDIORISK) from the European Research Council under the European Union’s Horizon 2020 Research and Innovation Program (H2020 European Research Council; Grant Agreement 638030) to R.P.D.
Ethics approval and consent to participate
The study protocol was approved by the local medical ethics review committees, including the Biobank Research Ethics Committee of University Medical Center Utrecht (protocol number 12/387), the Local Ethics Committee of the Erasmus MC (2010-409), the Washington University School of Medicine Ethics Committee (Institutional Review Board), the Sydney Heart Bank (HREC Univ. Sydney 2012/030) and the Leiden University Medical Center (LUMC). Informed consent was obtained from each patient prior to surgery or was waived by the ethics committee when acquiring informed consent was not possible due to the death of the donor.
Consent for publication
The authors declare that they have no competing interests.
Springer Nature remains neutral with regard to jurisdictional claims in published maps and institutional affiliations.
Examples of differentially acetylated regions between HCM and control hearts.
(A) Protein-protein interactions between up-regulated proteins in HCM versus control hearts. (B) Protein-protein interactions between down-regulated proteins in HCM versus control hearts
(A) Gene set enrichment report showing the correlation between extracellular matrix (EC)-related genes and genes annotated to differentially acetylated regions. (B) Gene set enrichment report showing the correlation between ECM-related genes and differentially expressed genes. (C) Gene set enrichment report showing the correlation between ECM-related genes and genes encoding differentially expressed proteins. (D) Gene set enrichment report showing the correlation between cardiac muscle contraction-related genes and genes annotated to differentially acetylated regions. (E) Gene set enrichment report showing the correlation between cardiac muscle contraction-related genes and differentially expressed genes. (F) Gene set enrichment report showing the correlation between cardiac muscle contraction-related genes and genes encoding differentially expressed proteins. (G) Gene set enrichment report showing the correlation between fatty acid metabolism-related genes and genes annotated to differentially acetylated regions. (H) Gene set enrichment report showing the correlation between fatty acid metabolism-related genes and differentially expressed genes. (I) Gene set enrichment report showing the correlation between fatty acid metabolism-related genes and genes encoding differentially expressed proteins
(A) Enriched transcription factor binding motifs in the hyperacetylated regions in HCM versus control hearts. (B) Enriched transcription factor binding motifs in the hypoacetylated regions in HCM versus control hearts
The acetylation and mRNA levels of the wildtype and mutant MYBPC3 alleles in all samples
Principal component analysis (PCA) plot showing the clustering between control samples and HCM with different truncating mutations in the MYBPC3 gene
The expression levels of obtained transcription factors in ventricular cardiomyocytes using published single-cell sequencing data
The expression levels of cardiomyocyte-specific markers in ventricular cardiomyocytes using published single-cell sequencing data
An overview of the clinical characteristics in all cardiac samples
(A) Differentially acetylated regions between HCM and control hearts from the H3K27ac ChIP-seq data. (B) Genes annotated to differentially acetylated regions using the 5 kb window of the transcription start site. (C) Genes annotated to differentially acetylated regions using the 50 kb window of the transcription start site. (D) Enriched GO terms and pathways by genes annotated to the hyperacetylated regions (50 kb window). (E) Enriched GO terms and pathways by genes annotated to the hypoacetylated regions (50 kb window)
(A) Differentially expressed genes between HCM and control hearts from the RNA-seq data. (B) Enriched GO terms and pathways by up-regulated genes in HCM versus control hearts. (C) Enriched GO terms and pathways by down-regulated genes in HCM versus control hearts
(A) Differentially expressed proteins between HCM and control hearts from the proteomics data. (B) Enriched GO terms and pathways by up-regulated proteins in HCM versus control hearts. (C) Enriched GO terms and pathways by down-regulated proteins in HCM versus control hearts
(A) Unsupervised analysis (O2PLS) ranking DNA regions that contributed to the separation between HCM and control hearts. (B) Unsupervised analysis (O2PLS) ranking genes that contributed to the separation between HCM and control hearts
(A) The overlapping genes between differentially expressed genes using DESeq2 and the top 2,000 genes that discriminate HCM from control hearts using O2PLS. (B) Enriched GO terms and pathways by the overlapping genes. (C) Enriched GO terms and pathways by a subset of the overlapping genes, which showed the same changing direction at their protein levels in HCM hearts when compared with controls
Enrichment results of identified transcription factors that were supported by both ChIP-seq and RNA-seq data in HCM versus control hearts
Gene sets collected from the Molecular Signature Databases
The acetylation and mRNA expression levels of included cell-type-specific markers
A published transcriptome data comparing genome-engineered cardiomyocytes with and without mutated MYBPC3 using RNA sequencing
About this article
Cite this article
Pei, J., Schuldt, M., Nagyova, E. et al. Multi-omics integration identifies key upstream regulators of pathomechanisms in hypertrophic cardiomyopathy due to truncating MYBPC3 mutations. Clin Epigenet 13, 61 (2021). https://doi.org/10.1186/s13148-021-01043-3