Novel myocardial markers GADD45G and NDUFS5 identified by RNA-sequencing predicts left ventricular reverse remodeling in advanced non-ischemic heart failure: a retrospective cohort study

Background Left ventricular reverse remodeling (LVRR) has been detected in non-ischemic dilated cardiomyopathy (NIDCM) patients following optimal treatment. However, its prediction with only conventional modalities is often difficult. This study sought to examine whether RNA sequencing (RNA-seq) of myocardium tissue samples could predict LVRR in NIDCM. Methods A total of 17 advanced NIDCM patients with left ventricular ejection fraction (LVEF) below 30% who underwent cardiac biopsy from Left ventricle (LV) were prospectively recruited. They received optimal treatment and followed with echocardiogram every 6 months. Based on LVRR status after 12 months of treatment, patients were divided into the reverse remodeling (RR) or non-RR group. Tissue samples were analyzed by RNA-seq, and a functional analysis of differentially expressed genes was carried out. Results There were eight and nine patients in the RR and non-RR groups, respectively. No difference was found in age, sex, disease duration, LV end-diastolic diameter, and LVEF between the two groups. There were 155 genes that were differentially expressed between the two groups. Nicotinamide adenine dinucleotide ubiquinone oxidoreductase subunit (NDUF)S5 and Growth arrest and DNA-damage-inducible protein (GADD)45G, along with several genes related to the mitochondrial respiratory chain and ribosome, were significantly downregulated in the RR as compared to the non-RR group. Conclusion GADD45G and NDUFS5 are potential biomarkers for LVRR in patients with advanced NIDCM.


Background
Heart failure has a complex etiology and is associated with various conditions, and is a major cause of adult mortality. In particular, non-ischemic dilated cardiomyopathy (NIDCM)-which is characterized by left ventricular (LV) dilation and a severely reduced LV ejection fraction (LVEF) without coronary artery disease-is often refractory to established drugs such as β-blockers, angiotensin-converting enzyme inhibitors, and angiotensin II receptor blockers and to cardiac resynchronization therapy (CRT). Therefore, NIDCM often requires heart transplantation or ventricular assist device (VAD) implantation. However, the former is available only to a limited number of patients, whereas the latter is associated with complications like stroke, bleeding, and infection that can be fatal during the long waiting period for transplantation [1].
On the other hand, many NIDCM patients show recovery of cardiac function after optimal pharmacological or device therapy [2,3]. This phenomenon, known as LV reverse remodeling (LVRR), is associated with improved clinical outcomes. Accordingly, several indices have been identified as predictors of LVRR and heart failure, including blood biomarkers like suppression of tumorigenesis-2 and galectin-3 [4][5][6][7][8], late gadolinium enhancement by cardiac magnetic resonance imaging [9], cardiac fibrosis and cardiomyocyte hypertrophy by histopathological examination [10,11], and washout rate in meta-iodobenzylguanidine myocardial scintigraphy [12]. However, conventional modalities cannot accurately predict the occurrence of LVRR.
Analyzing myocardial gene expression can be more useful for this purpose. Several studies have demonstrated a clear correlation between myocardial expression of certain mRNAs or micro (mi)RNAs and LVRR in NIDCM patients by quantitative reverse transcription PCR [13,14]. Moreover, RNA sequencing (RNA-seq) and next generation sequencing (NGS)-which enables higher throughput and greater accuracy in the evaluation of gene expression than microarray analysis [15]-is increasingly used in clinical settings [16]. Previous NGS studies have revealed differences in gene expression profiles between ischemic and non-ischemic cardiomyopathy [17,18] or pre-and post VAD support [19][20][21]. However, there have been no studies using comprehensive transcriptome analysis to examine specific genes related to LVRR. The present study therefore aimed to identify factors that can predict LVRR occurrence by RNA-seq.

Patient selection
Patients with advanced NIDCM whose LVEF was < 30% by echocardiogram and who underwent cardiac biopsy of the LV at Chiba University Hospital from September 2014 to May 2016 were included in the study. Exclusion criteria were as follows: (i) heart failure with ischemic etiology (defined as past history of myocardial infarction or significant stenosis in the main branches of coronary arteries); (ii) inflammatory or infiltrative heart disease including myocarditis, sarcoidosis, or amyloidosis; (iii) other secondary or metabolic cardiomyopathies such as neuromuscular disease and alcoholic cardiomyopathy; (iv) dependence on hemodialysis; and (v) age younger than 20 years or older than 65 years.

Patient follow-up
Prior to myocardial biopsy, all patients underwent blood sample collection, echocardiogram, and right heart catheterization. After myocardial biopsy, patients continued to receive optimal treatment including guidelinedirected medical therapy with β-blockers and reninangiotensin blockers and CRT. Patients were followed for 6 to 12 months, with echocardiogram performed at each follow up; they were then divided into the reverse remodeling (RR) and non-RR groups according to positive or negative LVRR status, respectively. LVRR was defined as (i) LVEF improvement of > 10% in absolute value, and (ii) LVEF > 30%.
RNA extraction and cDNA library preparation LV myocardium samples were collected at the time of cardiac catheter examination or open heart surgery, including for VAD implantation and valve replacement. Samples were immediately treated with RNAlater RNA stabilization reagent (Qiagen, Valencia, CA, USA) and stored at − 20°C. Total RNA was extracted using mir-Vana miRNA Isolation kit (Life Technologies, Carlsbad, CA, USA) according to the manufacturer's instructions, and cDNA libraries were generated using a NEBNext Ultra RNA Library Prep kit (New England BioLabs, Beverly, MA, USA). Sequencing was performed using a HiSeq1500 system (Illumina, San Diego, CA, USA) with a single-read sequencing length of 60 bp.

RNA-seq data analysis
TopHat v.1.3.2 with default parameters was used to map sequences against the UCSC/hg19 reference genome. Gene expression level was quantified using Cufflinks v.2.0.2 [22,23] as fragments per kilobase of exon per million fragments (FPKM) mapped reads. In the present study, miRNAs and genes expressed at low levels (FPKM < 1 in all samples) were eliminated. To identify novel biomarkers that predict the occurrence of LVRR, we carried out sensitivity, specificity, and receiver operating characteristic (ROC) curve analyses for each differentially expressed gene (DEG). Functional annotation was performed using the Database for Annotation, Visualization, and Integrated Discovery [24][25][26]. Gene expression profiles of two groups were compared by principal component analysis (PCA). RNA-sequencing data was deposited in the DNA Data Bank of Japan (accession number PRJDB5368).

Statistical analysis
For patients' clinical data, continuous variables are expressed as the mean ± standard deviation and were compared with the unpaired t test, whereas categorical data are expressed as a percentage and were compared with the χ 2 test. P < 0.05 was considered statistically significant. DEGs were defined as those with a significance (P) value < 0.05, false discovery rate < 5%, and fold change > 1.5, which was determined using the R v.3.1.0 software tag count comparison function [27]. The significance value for functional annotation analysis was P < 0.05. Data were analyzed using R version 3.3.2 (The R Foundation, Vienna, Austria).
Logistic regression analysis was used to assess the impact of selected DEGs on LVRR. Odds ratios (OR) and 95% confidence intervals (CI) were calculated. Confounders with P values < 0.05 in the univariate analysis were entered into the multivariate model. A P value < 0.05 was regarded as statistically significant. Data were analyzed using STATA version 15.1 (StataCorp, College Station, TX, USA).

Patient characteristics
A total of 20 patients were enrolled in this study. Two were excluded due to insufficient RNA amount and heart transplantation in the prior 6 months, respectively. One patient diagnosed with acute fulminant myocarditis was also excluded. Ultimately, 17 patients (3 females and 14 males, mean age: 46.6 ± 11.6 years) without coronary artery disease or other secondary or specific cardiomyopathies were analyzed. LV myocardial samples were obtained by transcatheter biopsy for four patients, by needle biopsy during aortic valve replacement surgery for three patients, or during VAD implantation for 10 patients.
Patients were divided into RR (n = 8) and non-RR (n = 9) groups based on the findings of a follow-up echocardiogram. There were no differences between the two groups in terms of age, sex, disease duration, LVEF on the echocardiogram, brain natriuretic peptide level, cardiac index for right heart catheterization, or ratio of cardioprotective medications (Table 1). On the other hand, proportion of surgical biopsy, inotropes used, and VAD implantation, and serum creatinine concentration were significantly higher in the non-RR than in the RR group (Table 1). Disease duration was more than double in the non-RR group but did not reach statistical significance.

RNA-seq analysis
A total of 22,416 genes were detected by RNA-seq. After excluding miRNAs and genes with low expression (FPKM < 1 in each sample), 14,448 genes were retained for analysis. Of these, 155 were differentially expressed between the RR and non-RR groups, with 150 genes upregulated in the non-RR group ( Table 2). The top three DEGs were LOC100506295, Nicotinamide adenine dinucleotide ubiquinone oxidoreductase subunit (NDUFS)5, and Growth arrest and DNA-damage-inducible protein (GADD)45G (Fig. 1). Genes encoding other NDUF subunits (NDUFB1, NDUFA1, NDUFS6, NDUFA13, NDUFA3, NDUFB7, and NDUFB3) as well as various ribosomal proteins were also upregulated. ROC curve analysis revealed that the three genes (LOC100506295, NDUFS5, and GADD45G) had high predictive capacity for LVRR ( Fig. 2; areas under the curve = 0.99, 0.93, and 0.90, respectively). A functional enrichment analysis showed significant enrichment of genes related to the mitochondrial respiratory chain or ribosome ( Table 3). The PCA showed separate distributions for the RR and non-RR groups (Fig. 3). qPCR analysis qPCR was performed with residual RNA samples. The RR group showed lower expression of NDUFS5 and GADD45G than the non-RR group, consistent with the results obtained by RNA-seq (Fig. 4).

Discussion
In the present study, we performed RNA-seq of myocardium tissue samples from patients with advanced NIDCM and compared the gene expression profiles of those with and without LVRR. Clinically, patients with LVRR showed markedly better renal function. They also needed significantly less inotrope support, surgical biopsy and VAD implantation through shorter disease duration than those without LVRR. The major findings were that (i) several genes including LOC100506295, NDUFS5 along with those encoding other NDUF subunits, GADD45G, and some genes encoding ribosomal protein subunits were upregulated in the non-RR as compared to the RR group; and (ii) genes associated with the mitochondrial respiratory chain and ribosome were enriched in the non-RR relative to the RR group.
Other than LOC100506295 whose function is unknown, NDUFS5 was most highly correlated with LVRR. NDUFS5 encodes a subunit of mitochondrial respiratory chain complex I [28]. Previous studies have shown that disturbance of mitochondrial function can cause heart failure [29][30][31], which is thought to result from increased reactive oxygen species (ROS) production and apoptosis [32]. Decreased complex I activity is also associated with increased ROS production [33][34][35], which can contribute to the progression of heart failure. Another report demonstrated that AF-HF001, a clinical drug candidate for heart failure, reversed the upregulation of NDUFS5 expression in H9c2 rat cardiomyocytes and attenuated ROS production and myocardial cell apoptosis [36]. These observations suggest that NDUFS5 along with other NDUF subunits of mitochondria complex I play a crucial role in the pathophysiology of heart failure. This is supported by the difference in NDUFS5 expression between the RR and non-RR groups as well as the enrichment of genes related to mitochondrial respiratory chain observed in the present study. Thus, NDUFS5 and mitochondrial complex I subunits may serve as predictive biomarkers for LVRR.
Inborn mitochondrial abnormality, i.e. mitochondrial disease (MD), can involve myocardium and induce cardiac dysfunction through disrupted mitochondrial respiratory chain functions [37]. In addition, even isolated cardiac phenotype is found in MD [38]. Thus, inherent but underdiagnosed MD might be included in the study population. This notion is supported by the previous report that 31% of MD patients represented cardiac conditions including cardiac dysfunction before the diagnosis [39]. In the present study, we did not examine mutations of mitochondrial DNA nor mitochondria-related nuclear DNA. However, given the lack of typical clinical signs in the participants of our study, MD with manifest cardiomyopathy such as Mitochondrial encephalomyopathy, lactic acidosis, and stroke-like episodes (MELAS), Chronic progressive external ophthalmoplegia (CPEO), and Myoclonic epilepsy with ragged-red fibers (MERRF) are unlikely. Further study is warranted to determine the prevalence of MD in this population. Apoptosis is another hallmark of mitochondrial function, with mitochondrial permeability transition pore opening to trigger intrinsic apoptosis pathway. Thus, disruption of mitochondrial function can enhance cardiomyocyte apoptosis, leading to progression of heart failure [40,41]. In this study, GADD45G was found to be strongly correlated with LVRR. GADD45G is a member of the GADD45 family of proteins that are involved in p38 mitogen-activated protein kinase-dependent cell death [42]. Recent studies have reported that GADD45G overexpression in mice induced cardiomyocyte apoptosis, fibrosis, LV dysfunction, and heart failure, whereas GADD45G deletion conferred resistance to ischemic injury and cardiomyocyte apoptosis [43]. In our study, GADD45G was more highly expressed in non-RR as compared to RR patients, suggesting that GADD45G upregulation is linked to induction of apoptosis and consequently, a reduced probability of LVRR.
Our study had two novel findings as compared to previous transcriptomic analyses of LVRR. First, we identified markers for predicting LVRR by comprehensive transcriptome analysis of single samples obtained prior to treatment. Some studies showed associations between certain mRNAs or miRNAs and myocardial recovery or subsequent VAD removal after LVAD implantation [13,14,44]. However, these studies analyzed only limited number of for the cardiac samples obtained from the patients with non-ischemic advanced heart failure requiring LVAD implantation, and showed that the expressions of myosin light chain kinase and interleukin-6 genes were significantly higher in the LVAD responders as compared to the non-responders. In this study, however, only 95 genes associated with heart failure were examined [44]. RNAseq in the present study performed a comprehensive transcriptome analysis, thereby demonstrating previously unknown, novel mechanisms of LVRR.
Secondly, cardiac biopsy can reveal not only histological abnormalities but also the potential reversibility of cardiac function, which increases its clinical utility for NIDCM. Cardiac biopsy-particularly by the transcatheter method-is associated with risks such as perforation or idiopathic valvular regurgitation, which has thus far limited the diagnostic yield. However, taking into   consideration information on the reversibility of cardiac function in NIDCM can facilitate clinical decisionmaking based on biopsy samples. Favorable expressions of the specific genes possibly make us advance aggressive medical treatment to avoid VAD implantation, or even to introduce VAD as a bridge to recovery strategy. By contrast, unfavorable expressions would lead us to prompt introduction of VAD as a bridge to transplant strategy. The present study had several limitations. First, this study was carried out at a single institution with a small number of samples. Previous studies examining myocardial samples from severely failing hearts could collect only around 10 samples for each treatment group [14,20,44]. These facts probably reflect the small population of patients whose heart failure status are severe enough to assess pathological specimen of the heart. Multicenter study will overcome both the possible bias from the single-centered study and the less-powered and nonvalidated results due to small sample size. Second, the issue of small sample size also influenced the results of logistic regression analyses. The model was not established with more than three independent variables due to small number of samples. Even with the challenge, GADD45G demonstrated almost significant predictive value for LVRR in the multivariate model, which highlights the strong contribution of the gene. Reconfirmation of the results with more large samples is needed to show clinically definite significance of NDUFS5 and GADD45G. Third, there were unequal number of females and males. Some transcripts might be expressed in a gender specific manner and serve as a bias for the results of the present study. Large-scale and preferably randomized control study is required to resolve this issue. Fourth, adopted treatment strategy differed between the groups. There were more LVAD implantation in the non-RR group, however, this difference was not arbitrary and resulted from the continuing optimal treatment for severe heart failure. Mechanical ventricular support is the most potent cardiac unloading therapy and is thus likely to promote LVRR more efficiently. Nonetheless most patients with LVRR implantation fell upon the non-RR group in the present study, underscoring that the non-RR group actually had little chance to recover cardiac function even after maximal heart failure treatment. Of note, gene expression analyses were performed on the samples prior to LVAD implantation and were not influenced by the difference in the treatment strategy. Fifth, we identified LVRR only by echocardiography in both VAD and non-VAD patients, in contrast to previous studies that defined LV recovery as removal of VAD in VAD-patients [14]. In the present study, there were two VAD patients in the RR group (showing elevations in LVEF from 18 to 39% and from 13 to 37%) who nonetheless required VAD support. However, since these patients showed at least partial recovery of cardiac function, pump flow was successfully weaned without right heart failure, aortic valve insufficiency, or readmission. Sixth, the patients had variable etiology, which included drug-induced cardiomyopathy or valvular heart disease. In fact, the etiology may not be known in some patients diagnosed as idiopathic DCM; additional studies are needed to address this issue. However, we presumed that the observed changes in gene expression are common to heart failure patients irrespective of etiology and can therefore serve as universal biomarkers. Finally, the sampled sites differed among patients-i.e., the apical and posterolateral walls were sampled by surgical and transcatheter biopsy, respectively-yet all patients showed a general decrease in cardiac function, implying homogeneous pathophysiology throughout the heart.

Conclusions
In conclusion, we identified novel genes that are correlated with LVRR by RNA-seq in myocardium tissue samples from patients with advanced NIDCM. NDUFS5 and GADD45G are potential biomarkers for predicting whether patients will develop LVRR, which can inform clinical decision making regarding VAD implantation and/or removal. These findings add value to conventional cardiac biopsy and can potentially improve the treatment of heart failure patients.