-
PDF
- Split View
-
Views
-
Cite
Cite
Isabel Witvrouwen, Andreas B Gevaert, Nadine Possemiers, Paul J Beckers, Anne Vorlat, Hein Heidbuchel, Steven J Van Laere, Amaryllis H Van Craenenbroeck, Emeline M Van Craenenbroeck, Circulating microRNA as predictors for exercise response in heart failure with reduced ejection fraction, European Journal of Preventive Cardiology, Volume 28, Issue 15, December 2021, Pages 1673–1681, https://doi.org/10.1093/eurjpc/zwaa142
- Share Icon Share
Abstract
Exercise training is a powerful adjunctive therapy in patients with heart failure with reduced ejection fraction (HFrEF), but ca. 55% of patients fail to improve VO2peak. We hypothesize that circulating microRNAs (miRNAs), as epigenetic determinants of VO2peak, can distinguish exercise responders (ER) from exercise non-responders (ENR).
We analysed 377 miRNAs in 18 male HFrEF patients (9 ER and 9 ENR) prior to 15 weeks of exercise training using a miRNA array. ER and ENR were defined as change in VO2peak of >20% or <6%, respectively. First, unsupervised clustering analysis of the miRNA pattern was performed. Second, differential expression of miRNA in ER and ENR was analysed and related to percent change in VO2peak. Third, a gene set enrichment analysis was conducted to detect targeted genes and pathways. Baseline characteristics and training volume were similar between ER and ENR. Unsupervised clustering analysis of miRNAs distinguished ER from ENR with 83% accuracy. A total of 57 miRNAs were differentially expressed in ENR vs. ER. A panel of seven miRNAs up-regulated in ENR (Let-7b, miR-23a, miR-140, miR-146a, miR-191, miR-210, and miR-339-5p) correlated with %changeVO2peak (all P < 0.05) and predicted ENR with area under the receiver operating characteristic curves ≥0.77. Multiple pathways involved in exercise adaptation processes were identified.
A fingerprint of seven miRNAs involved in exercise adaptation processes is highly correlated with VO2peak trainability in HFrEF, which holds promise for the prediction of training response and patient-targeted exercise prescription.
Introduction
Heart failure (HF) is an increasingly prevalent syndrome with considerable impact on quality of life.1 Exercise training is a powerful adjunctive therapy that improves morbidity and mortality and is therefore recommended to all stable HF patients (Class IA indication).2 Unfortunately, not all HF patients benefit from this approach and ∼55% of HF patients fail to demonstrate a clinically relevant increase in aerobic capacity.3 Peak oxygen consumption (VO2peak) is one of the strongest prognostic factors in HF and a failure to adequately increase VO2peak adds to an adverse prognosis, independent of other risk factors.4
Even with similar training volumes, the variability in VO2peak response is high and underlying mechanisms are not fully explained. Heritability accounts for at least 40–50% of the anticipated effect of exercise training.5 However, at present, it is not possible to predict which patients will show an increase in VO2peak following exercise training.
MicroRNAs (miRNAs) are non-coding RNAs that regulate gene expression at the post-transcriptional level. miRNAs interrupt translation of messenger RNA (mRNA) through base pairing, provoking translational repression or degradation of the mRNA. One miRNA may exert inhibitory effects on several mRNAs.6 miRNAs can be detected in plasma, either packed in exosomes or microparticles, or bound to Ago proteins and HDL cholesterol, making them attractive as biomarkers.7
Plasma levels of miRNAs have been linked to aerobic capacity and they change dynamically following exercise.8 For example, miR-21, miR-210, miR-222 are elevated in healthy individuals with lower VO2peak,9 and miR-1, miR-20a, miR-146a, and miR-486 are higher in endurance athletes with a higher VO2peak.10,11 Exercise training increased plasma levels of miR-21, miR-146a, miR-221, and miR-222 in athletes.10 In HF patients, an acute maximal exercise bout significantly up-regulated circulating miR-21, miR-378, and miR-940 levels immediately after the exercise test,12 and in patients with chronic kidney disease, plasma levels of miR-146a showed a rapid down-regulation after an acute exhaustive exercise test.13
As such, circulating miRNAs are promising as epigenetic markers of physical fitness and exercise-induced cardiovascular adaptation, and they could even play a role in personalized exercise prescription. In this study, we hypothesize that a specific miRNA signature could distinguish exercise non-responders (ENR) from exercise responders (ER), hence allowing to predict which patients will show an increase in VO2peak following exercise training. To this end, we performed an unbiased miRNA screening in HF with reduced ejection fraction (HFrEF) patients and related this to the VO2peak response to exercise training.
Methods
Study design, exercise training, and testing
For this retrospective cohort study, male HFrEF patients [left ventricular ejection fraction (LVEF) <40%] who followed an in-hospital training programme in the setting of a standardized longitudinal study, were screened. To minimize the impact of sex-specific patterns on the miRNA profile, this study was performed in male patients only.14 All patients had to be clinically stable and optimally medically treated for ≥6 weeks. Eighteen patients were included. The programme consisted of 15 weeks moderate aerobic exercise training, three 50-min sessions/week, with a continuous training intensity at 90% of heart rate (HR) at the respiratory compensation point (RCP, n = 17) or four 4-min intervals at 90% of HR at the RCP (n = 1). In addition, all patients performed 10 min of moderate intensity resistance exercise per session. At start and after 15 weeks of training, cardiopulmonary exercise test (CPET) was performed on a bicycle or treadmill ergometer (the same modality at both visits) and VO2peak was determined as the mean VO2peak during the final 30 s of exercise. RCP was identified using the VE/VCO2 curve and the CO2-equivalent. Whole blood was collected at baseline and after 15 weeks after an overnight fast in EDTA tubes, centrifuged within 30 min after collection and stored at −80°C. miRNAs were quantified from baseline blood samples.
This study complied with the Declaration of Helsinki and was approved by the ethical committee of the Antwerp University Hospital. Written informed consent was obtained from each participant.
Definition of exercise responders and exercise non-responders
ER were defined as subjects with an increase of >20% in VO2peak, and ENR as increase of <6% in VO2peak. An improvement of at least 6% in VO2peak has been associated with reduced all-cause mortality and all-cause hospitalization.15 To increase the discriminatory capacity of the miRNA panel, the cut-off for ER was set at >20% VO2peak increase.
miRNA array
miRNAs were profiled and analysed from plasma samples using TaqMan Low Density MicroRNA Array (TLDA) Human Cards A (ThermoFisher), analysing 377 human miRNAs, as previously published.16,17 Briefly, plasma samples were thawed on ice and centrifuged for 10 min (4°C, 16 000 g). Total RNA was isolated using the mirVana Paris kit (ThermoFisher). Reverse transcription and preamplification were performed with MegaPlex primer pools (ThermoFisher) following the manufacturer’s protocol.16 The preamplification product was mixed with TaqMan Universal PCR Master Mix No AmpErase UNG (ThermoFisher) and nuclease-free water before loading to the TLDA card. The arrays were run in a 7900HT Fast Real-Time polymerase chain reaction (PCR) system (ThermoFisher). Raw cycle quantification (Cq) values were calculated in SDS software v.2.4 using automatic baseline and threshold settings. A miRNA was considered non-informative if Cq values were >35 in >80% of samples. As suggested before, geNorm algorithm (NormqPCR package) was used for normalization and relative miRNA levels were expressed as 2−ΔCq.17
Technical validation by miRNA RT-qPCR
Expression of the selected miRNA was repeated by conventional RT-qPCR. A new plasma aliquot of the same individuals was thawed on ice and centrifuged for 10 min (4°C, 16 000 g). RNA enriched for small RNAs was isolated using the mirVana Paris Kit (ThermoFisher). As spike-in control, 20 fmol synthetic Ath-miR-159a (ThermoFisher) was added. Reverse transcription and preamplification were performed using TaqMan miRNA primers (ThermoFisher) and multiplex qPCR was done in a CFX96 thermal cycler (BioRad).13 Raw Cq values were calculated in BioRad CFX manager software v.3.1 using automatic baseline and threshold settings. Data were normalized using Ath-miR-159a and relative miRNA levels were expressed as 2−ΔCq.
Statistical analysis
Data were analysed using R version 3.4.3, SPSS 26.0, and GraphPad Prism 8.3.0. Normality of continuous variables was evaluated using Shapiro–Wilk test. Normally distributed data are expressed as mean ± standard deviation, skewed variables as median and range (1st–3rd quartile). Categorical variables were analysed with Fisher’s exact test, continuous variables with independent samples T-test or Mann–Whitney U test where appropriate.
Unsupervised agglomerative hierarchical clustering was performed using Manhattan distances and Ward linkage, and a heatmap was constructed (cluster and Heatplus packages). The accuracy of this unsupervised model was calculated by superposing unblinded information (responder and non-responder status of the patients) to the clusters.
To identify a miRNA panel that discriminates ENR from ER, miRNA expression between ER and ENR was compared using multiple T-tests and resulting P-values were adjusted for multiple testing using the Benjamini and Hochberg procedure controlling the false discovery rate at 5%. miRNAs were considered significantly differentially expressed when fold change was <0.66 or >1.5, and the adjusted P-value was <0.05, between ER and ENR.
Next, a supervised shrunken centroid model was fitted to predict responder status from baseline miRNA expression (pamr package, threshold 1.0).18 We predefined criteria to further refine the miRNA selection: since we are interested in a biomarker of ENR, miRNAs of interest had to be (i) expressed in all ENR, (ii) up-regulated in ENR, (iii) P < 0.005 for the difference in fold change expression between ENR vs. ER, and (iv) similarly expressed in arrays and RT-qPCR.
Correlation between miRNAs and percent change in VO2peak (%changeVO2peak) was assessed using Spearman correlation analysis. A univariate logistic regression model to predict ENR using relevant baseline characteristics [age, baseline LVEF, and New York Heart Association (NYHA) class]19 was fitted. Area under the receiver operating characteristic curve (AUC) was calculated to predict ENR using miRNAs of interest. Linear mixed models were fitted using time (baseline vs. 15 weeks) and group (ER vs. ENR) as fixed effects and patient ID as random effect, to assess change in CPET characteristics, body mass index, LVEF, strength, and NYHA class. A two-sided P-value <0.05 was considered significant.
Interaction network
To perform a pathway analysis, a three-layer interaction network was constructed:
First, a protein–protein interaction matrix was constructed using STRING database v10.5, specific for heart, blood vessel, skeletal muscle, and kidney—tissues relevant for HF pathophysiology—and annotated with their corresponding genes using GTex database v7. Interactions were relevant when read counts for both genes were >10 in at least one of the relevant tissues.
Second, a miRNA–miRNA interaction matrix was constructed based on a mutual information network after applying the Aracne algorithm, using default settings, on the miRNA expression data.
Third, TargetScan database v7.0 was used to predict human miRNA–gene interactions.
Then, a Markov random walk algorithm was used to smooth expression fold changes for miRNAs with a significantly altered expression over the constructed three-layer network. Three different restart probabilities (0.8, 0.5, and 0.2, respectively limited, average, and extensive smoothing) were executed, thereby simulating the amount of information flow through the network and generating a probability distribution reflecting the likelihood that a gene in the network is perturbed by a miRNA of interest. Next, a gene set enrichment analysis was performed, using the genes ranked in decreasing order according to the perturbation probabilities. Here, 56 gene sets associated with biological processes and signal transduction pathways retrieved from the Kyoto Encyclopaedia of Genes and Genomes (KEGG, retrieved on 20 November 2017) were used. The enrichment score normalized for gene set length and the associated significance level, estimated using 10 000 permutations, were recorded for each of the evaluated gene sets. Resulting P-values were adjusted for multiple testing using the Benjamini and Hochberg procedure controlling the false discovery rate at 5%. Consolidated pathways were those that remained after removal of gene overlap and only positive enriched pathways were selected.
This analysis was performed using data obtained with the different restart probability levels. For network visualization, a minimal spanning tree algorithm was first applied onto the network, to extract the most representative interactions from these dense networks.
Results
Patient characteristics and association with response
From a total of 41 HFrEF patients, nine ER (>20% change in VO2peak) and nine age-matched ENR (<6% change in VO2peak) that underwent an identical training programme, were withheld for further analysis. The evolution of VO2peak is shown in Figure 1. At baseline, ER and ENR were similar with regard to patient demographics, clinical, pharmacological, and CPET characteristics (all P > 0.05, Table 1). Adherence to training was excellent in both groups with 37 completed sessions out of 45 total sessions. ENR trained at 106% and ER at 95% of their target HR during the first 4 weeks (P = 0.040), and ENR at 111% and ER at 96% of their target HR from week 5 to 15 (P = 0.036). Age, baseline LVEF, and NYHA class did not predict exercise response in univariate logistic regression analysis (P > 0.05). There was no difference in change in medical therapy between the groups (P > 0.05, data not shown). Whereas not significantly different between groups (P-value for interaction > 0.05), there was a trend towards improved strength, LVEF and NYHA class in the ER compared to ENR group (Supplementary material online, Table S1).

Change in exercise capacity, measured by VO2peak. Linear mixed model using time and group as fixed effects and patient ID as random effect. ENR, exercise non-responder; ER, exercise responder; V1, visit 1; V2, visit 2.
Baseline clinical, pharmacological, cardiopulmonary exercise test, and training characteristics
. | ER (n = 9) . | ENR (n = 9) . | P-value . |
---|---|---|---|
Clinical characteristics | |||
Age (years) | 59.4 (50.7–65.4) | 62.2 (59.8–65.6) | 0.387 |
Male sex | 100% | 100% | 1.0 |
BMI (kg/m2) | 30.5 ± 3.9 | 30.3 ± 5.2 | 0.940 |
Diabetes (%) | 22% | 33% | 1.0 |
Hypertension (%) | 56% | 33% | 0.637 |
NYHA class | II = 89%, III = 11% | II = 44%, III = 56% | 0.131 |
Aetiology of heart failure | 0.793 | ||
Ischaemic cardiomyopathy | 6/9 (67%) | 4/9 (44%) | |
Dilated cardiomyopathy | 2/9 (22%) | 4/9 (44%) | |
Toxic cardiomyopathy (ethyl) | 1/9 (11%) | 1/9 (11%) | |
LV ejection fraction (%) | 26.9 ± 7.7 | 24.6 ± 9.9 | 0.584 |
CRT or ICD | CRT (11%, ICD 22%) | CRT (11%, ICD 56%) | 0.361 |
Pharmacological therapy | CRT (11%), ICD (22%) | CRT (11%), ICD (56%) | |
ACE inhibitor (%) | 89% | 78% | 1.0 |
ARB (%) | 0% | 22% | 0.471 |
Beta blocker (%) | 100% | 100% | 1.0 |
Aldosterone antagonist (%) | 56% | 56% | 1.0 |
Diuretic (%) | 89% | 89% | 1.0 |
Cardiopulmonary exercise test variables | |||
Resting heart rate (b.p.m.) | 68 ± 16.6 | 76 ± 15.0 | 0.280 |
Baseline VO2peak (mL/kg/min) | 17.5 ± 3.4 | 17.2 ± 3.0 | 0.868 |
% Predicted VO2peak (%) | 63.7 ± 11.1 | 64.4 ± 8.6 | 0.882 |
Work economy (W/mL/kg/min) | 6.5 ± 1.6 | 5.8 ± 0.9 | 0.293 |
Peak systolic blood pressure (mmHg) | 124 ± 36.0 | 135 ± 35.6 | 0.532 |
Load (W) | 116.7 ± 46.1 | 101.1 ± 25.7 | 0.390 |
VE/VCO2 slope | 27.7 ± 5.0 | 33.2 ± 7.7 | 0.094 |
Training adherence | |||
Sessions completed (max. 45 sessions) | 37 ± 4.2 | 37 ± 3.4 | 0.719 |
. | ER (n = 9) . | ENR (n = 9) . | P-value . |
---|---|---|---|
Clinical characteristics | |||
Age (years) | 59.4 (50.7–65.4) | 62.2 (59.8–65.6) | 0.387 |
Male sex | 100% | 100% | 1.0 |
BMI (kg/m2) | 30.5 ± 3.9 | 30.3 ± 5.2 | 0.940 |
Diabetes (%) | 22% | 33% | 1.0 |
Hypertension (%) | 56% | 33% | 0.637 |
NYHA class | II = 89%, III = 11% | II = 44%, III = 56% | 0.131 |
Aetiology of heart failure | 0.793 | ||
Ischaemic cardiomyopathy | 6/9 (67%) | 4/9 (44%) | |
Dilated cardiomyopathy | 2/9 (22%) | 4/9 (44%) | |
Toxic cardiomyopathy (ethyl) | 1/9 (11%) | 1/9 (11%) | |
LV ejection fraction (%) | 26.9 ± 7.7 | 24.6 ± 9.9 | 0.584 |
CRT or ICD | CRT (11%, ICD 22%) | CRT (11%, ICD 56%) | 0.361 |
Pharmacological therapy | CRT (11%), ICD (22%) | CRT (11%), ICD (56%) | |
ACE inhibitor (%) | 89% | 78% | 1.0 |
ARB (%) | 0% | 22% | 0.471 |
Beta blocker (%) | 100% | 100% | 1.0 |
Aldosterone antagonist (%) | 56% | 56% | 1.0 |
Diuretic (%) | 89% | 89% | 1.0 |
Cardiopulmonary exercise test variables | |||
Resting heart rate (b.p.m.) | 68 ± 16.6 | 76 ± 15.0 | 0.280 |
Baseline VO2peak (mL/kg/min) | 17.5 ± 3.4 | 17.2 ± 3.0 | 0.868 |
% Predicted VO2peak (%) | 63.7 ± 11.1 | 64.4 ± 8.6 | 0.882 |
Work economy (W/mL/kg/min) | 6.5 ± 1.6 | 5.8 ± 0.9 | 0.293 |
Peak systolic blood pressure (mmHg) | 124 ± 36.0 | 135 ± 35.6 | 0.532 |
Load (W) | 116.7 ± 46.1 | 101.1 ± 25.7 | 0.390 |
VE/VCO2 slope | 27.7 ± 5.0 | 33.2 ± 7.7 | 0.094 |
Training adherence | |||
Sessions completed (max. 45 sessions) | 37 ± 4.2 | 37 ± 3.4 | 0.719 |
ACE, angiotensin-converting enzyme; ARB, angiotensin II receptor blockers; b.p.m., beats per minute; BMI, body mass index; CRT, cardiac resynchronization therapy; ENR, exercise non-responder; ER, exercise responder; ICD, implantable cardioverter-defibrillator; LV, left ventricular; NYHA class, New York Heart Association functional class.
Baseline clinical, pharmacological, cardiopulmonary exercise test, and training characteristics
. | ER (n = 9) . | ENR (n = 9) . | P-value . |
---|---|---|---|
Clinical characteristics | |||
Age (years) | 59.4 (50.7–65.4) | 62.2 (59.8–65.6) | 0.387 |
Male sex | 100% | 100% | 1.0 |
BMI (kg/m2) | 30.5 ± 3.9 | 30.3 ± 5.2 | 0.940 |
Diabetes (%) | 22% | 33% | 1.0 |
Hypertension (%) | 56% | 33% | 0.637 |
NYHA class | II = 89%, III = 11% | II = 44%, III = 56% | 0.131 |
Aetiology of heart failure | 0.793 | ||
Ischaemic cardiomyopathy | 6/9 (67%) | 4/9 (44%) | |
Dilated cardiomyopathy | 2/9 (22%) | 4/9 (44%) | |
Toxic cardiomyopathy (ethyl) | 1/9 (11%) | 1/9 (11%) | |
LV ejection fraction (%) | 26.9 ± 7.7 | 24.6 ± 9.9 | 0.584 |
CRT or ICD | CRT (11%, ICD 22%) | CRT (11%, ICD 56%) | 0.361 |
Pharmacological therapy | CRT (11%), ICD (22%) | CRT (11%), ICD (56%) | |
ACE inhibitor (%) | 89% | 78% | 1.0 |
ARB (%) | 0% | 22% | 0.471 |
Beta blocker (%) | 100% | 100% | 1.0 |
Aldosterone antagonist (%) | 56% | 56% | 1.0 |
Diuretic (%) | 89% | 89% | 1.0 |
Cardiopulmonary exercise test variables | |||
Resting heart rate (b.p.m.) | 68 ± 16.6 | 76 ± 15.0 | 0.280 |
Baseline VO2peak (mL/kg/min) | 17.5 ± 3.4 | 17.2 ± 3.0 | 0.868 |
% Predicted VO2peak (%) | 63.7 ± 11.1 | 64.4 ± 8.6 | 0.882 |
Work economy (W/mL/kg/min) | 6.5 ± 1.6 | 5.8 ± 0.9 | 0.293 |
Peak systolic blood pressure (mmHg) | 124 ± 36.0 | 135 ± 35.6 | 0.532 |
Load (W) | 116.7 ± 46.1 | 101.1 ± 25.7 | 0.390 |
VE/VCO2 slope | 27.7 ± 5.0 | 33.2 ± 7.7 | 0.094 |
Training adherence | |||
Sessions completed (max. 45 sessions) | 37 ± 4.2 | 37 ± 3.4 | 0.719 |
. | ER (n = 9) . | ENR (n = 9) . | P-value . |
---|---|---|---|
Clinical characteristics | |||
Age (years) | 59.4 (50.7–65.4) | 62.2 (59.8–65.6) | 0.387 |
Male sex | 100% | 100% | 1.0 |
BMI (kg/m2) | 30.5 ± 3.9 | 30.3 ± 5.2 | 0.940 |
Diabetes (%) | 22% | 33% | 1.0 |
Hypertension (%) | 56% | 33% | 0.637 |
NYHA class | II = 89%, III = 11% | II = 44%, III = 56% | 0.131 |
Aetiology of heart failure | 0.793 | ||
Ischaemic cardiomyopathy | 6/9 (67%) | 4/9 (44%) | |
Dilated cardiomyopathy | 2/9 (22%) | 4/9 (44%) | |
Toxic cardiomyopathy (ethyl) | 1/9 (11%) | 1/9 (11%) | |
LV ejection fraction (%) | 26.9 ± 7.7 | 24.6 ± 9.9 | 0.584 |
CRT or ICD | CRT (11%, ICD 22%) | CRT (11%, ICD 56%) | 0.361 |
Pharmacological therapy | CRT (11%), ICD (22%) | CRT (11%), ICD (56%) | |
ACE inhibitor (%) | 89% | 78% | 1.0 |
ARB (%) | 0% | 22% | 0.471 |
Beta blocker (%) | 100% | 100% | 1.0 |
Aldosterone antagonist (%) | 56% | 56% | 1.0 |
Diuretic (%) | 89% | 89% | 1.0 |
Cardiopulmonary exercise test variables | |||
Resting heart rate (b.p.m.) | 68 ± 16.6 | 76 ± 15.0 | 0.280 |
Baseline VO2peak (mL/kg/min) | 17.5 ± 3.4 | 17.2 ± 3.0 | 0.868 |
% Predicted VO2peak (%) | 63.7 ± 11.1 | 64.4 ± 8.6 | 0.882 |
Work economy (W/mL/kg/min) | 6.5 ± 1.6 | 5.8 ± 0.9 | 0.293 |
Peak systolic blood pressure (mmHg) | 124 ± 36.0 | 135 ± 35.6 | 0.532 |
Load (W) | 116.7 ± 46.1 | 101.1 ± 25.7 | 0.390 |
VE/VCO2 slope | 27.7 ± 5.0 | 33.2 ± 7.7 | 0.094 |
Training adherence | |||
Sessions completed (max. 45 sessions) | 37 ± 4.2 | 37 ± 3.4 | 0.719 |
ACE, angiotensin-converting enzyme; ARB, angiotensin II receptor blockers; b.p.m., beats per minute; BMI, body mass index; CRT, cardiac resynchronization therapy; ENR, exercise non-responder; ER, exercise responder; ICD, implantable cardioverter-defibrillator; LV, left ventricular; NYHA class, New York Heart Association functional class.
Differential miRNA expression between ENR and ER
Unsupervised analysis of the miRNA array expression data in all patients revealed two separate clusters of patients (Figure 2). These clusters represent ENR and ER with an accuracy of 83%, and a classification error in ER = 22%, in ENR = 11%, which results in an overall classification error of 17% (one ENR misclassified in the ER group, two ER in the ENR group). Clear differences in miRNA expression between the two clusters are illustrated in a heatmap (Figure 2).

Unsupervised cluster analysis and heatmap. Each row represents a miRNA and each column a HFrEF patient. Two ER patients and one ENR patient are misclassified, resulting in an overall classification error of 17% and an accuracy of 83%. ENR, exercise non-responder; ER, exercise responder; green, up-regulated miRNAs; Red, down-regulated miRNAs.
Supervised analysis revealed differential expression of 57 miRNAs in the ENR vs. ER: 26 miRNAs were up-regulated and 31 were down-regulated in ENR (Figure 3, P < 0.05). A supervised shrunken centroid model correctly classified 89% of patients as ER or ENR (two ER were misclassified as ENR) using only baseline miRNA expression and responder status as input.

Volcanoplot of relative miRNA expression fold changes. Applying predefined criteria, miRNAs of interest had to be expressed in all ENR, up-regulated in ENR and P < 0.005 for the difference in fold change expression between ENR and ER (red dots). Blue dots, significantly up-regulated miRNAs in ENR and ER (P < 0.05); ENR, exercise non-responder; ER, exercise responder; Green dots, highly significantly up-regulated miRNAs in ER (P < 0.005).
Discriminative miRNA signature for ENR
To develop a discriminative miRNA signature for ENR, predefined selection criteria were applied. miRNAs of interest had to be (i) expressed in all ENR, (ii) up-regulated in ENR, (iii) P < 0.005 for the difference in fold change expression between ENR and ER, and (iv) similarly expressed in arrays and RT-qPCR. This resulted in ENR miRNA signature of 7 miRNAs: miR-23a; miR-339-5p, miR-140, miR-191, miR-210, miR-146a, and Let-7b (Figure 3). All seven miRNAs were significantly correlated with %changeVO2peak (P < 0.05), and all miRNAs of interest had an AUC ≥0.77 for the identification of ENR (Table 2, Supplementary material online, Figure S1).
miRNA . | Amplification . | ER vs. ENR . | Spearman correlation with %changeVO2peak . | AUC to predict ENR . | ||
---|---|---|---|---|---|---|
Fold change . | P-value . | Rho . | P-value . | . | ||
Let-7b | 18/18 | 0.321 | 0.004 | −0.64 | 0.004 | 0.89 |
miR-23a | 13/18 | 0.234 | <0.001 | −0.74 | <0.001 | 0.98 |
miR-140 | 18/18 | 0.443 | 0.001 | −0.86 | <0.001 | 0.94 |
miR-146a | 18/18 | 0.456 | 0.003 | −0.70 | 0.001 | 0.88 |
miR-191 | 18/18 | 0.330 | 0.004 | −0.72 | 0.001 | 0.91 |
miR-210 | 12/18 | 0.165 | 0.004 | −0.64 | 0.004 | 0.91 |
miR-339-5p | 9/18 | 0.076 | <0.001 | −0.52 | 0.027 | 0.77 |
miRNA . | Amplification . | ER vs. ENR . | Spearman correlation with %changeVO2peak . | AUC to predict ENR . | ||
---|---|---|---|---|---|---|
Fold change . | P-value . | Rho . | P-value . | . | ||
Let-7b | 18/18 | 0.321 | 0.004 | −0.64 | 0.004 | 0.89 |
miR-23a | 13/18 | 0.234 | <0.001 | −0.74 | <0.001 | 0.98 |
miR-140 | 18/18 | 0.443 | 0.001 | −0.86 | <0.001 | 0.94 |
miR-146a | 18/18 | 0.456 | 0.003 | −0.70 | 0.001 | 0.88 |
miR-191 | 18/18 | 0.330 | 0.004 | −0.72 | 0.001 | 0.91 |
miR-210 | 12/18 | 0.165 | 0.004 | −0.64 | 0.004 | 0.91 |
miR-339-5p | 9/18 | 0.076 | <0.001 | −0.52 | 0.027 | 0.77 |
Serial two-sample T-tests comparing miRNA fold changes in ER vs. ENR. The association of miRNAs of interest with %changeVO2peak was assessed using Spearman correlation, the association with ENR using area under the curve. Fold change P-values are false-discovery rate-adjusted.
AUC, area under the curve; ENR, exercise non-responder; ER, exercise responder.
miRNA . | Amplification . | ER vs. ENR . | Spearman correlation with %changeVO2peak . | AUC to predict ENR . | ||
---|---|---|---|---|---|---|
Fold change . | P-value . | Rho . | P-value . | . | ||
Let-7b | 18/18 | 0.321 | 0.004 | −0.64 | 0.004 | 0.89 |
miR-23a | 13/18 | 0.234 | <0.001 | −0.74 | <0.001 | 0.98 |
miR-140 | 18/18 | 0.443 | 0.001 | −0.86 | <0.001 | 0.94 |
miR-146a | 18/18 | 0.456 | 0.003 | −0.70 | 0.001 | 0.88 |
miR-191 | 18/18 | 0.330 | 0.004 | −0.72 | 0.001 | 0.91 |
miR-210 | 12/18 | 0.165 | 0.004 | −0.64 | 0.004 | 0.91 |
miR-339-5p | 9/18 | 0.076 | <0.001 | −0.52 | 0.027 | 0.77 |
miRNA . | Amplification . | ER vs. ENR . | Spearman correlation with %changeVO2peak . | AUC to predict ENR . | ||
---|---|---|---|---|---|---|
Fold change . | P-value . | Rho . | P-value . | . | ||
Let-7b | 18/18 | 0.321 | 0.004 | −0.64 | 0.004 | 0.89 |
miR-23a | 13/18 | 0.234 | <0.001 | −0.74 | <0.001 | 0.98 |
miR-140 | 18/18 | 0.443 | 0.001 | −0.86 | <0.001 | 0.94 |
miR-146a | 18/18 | 0.456 | 0.003 | −0.70 | 0.001 | 0.88 |
miR-191 | 18/18 | 0.330 | 0.004 | −0.72 | 0.001 | 0.91 |
miR-210 | 12/18 | 0.165 | 0.004 | −0.64 | 0.004 | 0.91 |
miR-339-5p | 9/18 | 0.076 | <0.001 | −0.52 | 0.027 | 0.77 |
Serial two-sample T-tests comparing miRNA fold changes in ER vs. ENR. The association of miRNAs of interest with %changeVO2peak was assessed using Spearman correlation, the association with ENR using area under the curve. Fold change P-values are false-discovery rate-adjusted.
AUC, area under the curve; ENR, exercise non-responder; ER, exercise responder.
Targets of differentially regulated miRNAs
In order to detect the genes and pathways targeted by the miRNA profile of ENR and ER, a pathway analysis was performed. First, a miRNA–gene interaction network was constructed and a Markov random walk diffusion algorithm was used to rank network genes based on their probability of being perturbed by individual or joint miRNA expression changes. The most important interactions in the miRNA–gene network are shown in Supplementary material online, Figure S2A and B, with key roles for PIK3C2A, DNM2, RAB5C, and HSPA8 genes (targeted by up-regulated miRNAs in ENR), and CRK, EIF4B, and PRKG1 genes (targeted by up-regulated miRNAs in ER). Next, a gene set enrichment analysis was performed, translating the network-based gene perturbation probabilities into pathways. KEGG pathways related to nucleotide-binding oligomerization domain-like receptors (NOD-like receptors), transforming growth factor β (TGF-β), toll-like receptor, adherens junction, apelin signalling, neurotrophin signalling, and miRNAs in cancer were consolidated in the analysis of up-regulated miRNAs in ENR, whereas Notch, mitogen-associated protein kinase (MAPK) and vascular endothelial growth factor (VEGF), epidermal growth factor receptor tyrosine kinase inhibitor resistance, hippo signalling, adherens junction, apelin signalling, neurotrophin signalling pathways and miRNAs in cancer were consolidated in the analysis of up-regulated miRNAs in ER (Figure 4A and B).

Positively enriched signalling pathways. These result from a gene set enrichment analysis, grouping genes in KEGG pathways. (A) KEGG pathways targeted by up-regulated miRNAs in ENR. (B) KEGG pathways targeted by up-regulated miRNAs in ER. NES is plotted for each pathway and for three different restart probabilities (respectively 0.2 = extensive smoothing, 0.5 = average smoothing and 0.8 = limited smoothing) of the Markov random walk diffusion. Blue dot, consolidated pathways; KEGG, Kyoto Encyclopaedia of Genes and Genomes; NES, normalized enrichment score.
Discussion
The present study describes for the first time a plasma miRNA profile that distinguishes HFrEF patients with a favourable response to exercise training from ENR, despite excellent training adherence and similar patient characteristics in both groups. The panel of seven miRNAs (Let-7b, miR-23a, miR-140, miR-146a, miR-191, miR-210, and miR-339-5p), up-regulated in ENR, is involved in exercise adaptation processes such as angiogenesis, skeletal muscle function, and inflammation. In silico gene set enrichment analysis revealed several pathways involved in the regulation of exercise response, which may serve as novel therapeutic targets.
Clinical variables to predict exercise response
Early identification of ENR is of great importance, as 55% of HFrEF patients show impaired VO2peak response (<6% increase in VO2peak) after standardized cardiac rehabilitation.3 Indeed, patient-tailored training modifications, either through changing the type or through increasing the duration, intensity or frequency of the training programme, might result in a more favourable and cost-effective VO2peak response.19–21 In the past decade, several studies have investigated different clinical and training-related determinants of exercise response, but a clear discriminative marker is still lacking. In the SMARTEX-HF study, investigating the effect of exercise intensity in 215 HFrEF patients, lower NYHA class, younger age, higher LVEF, and high-intensity interval training (HIIT) or moderate continuous training (MCT) significantly increased the odds for being a VO2peak responder.19 In a large study on VO2peak trainability, including 677 participants of both healthy, elderly and clinical populations (coronary artery disease, type-2 diabetes, and metabolic syndrome), age, sex, exercise volume, population group, and the average between pre- and post-training VO2peak explained only 17% of the variance in VO2peak trainability.20 In patients with coronary artery disease, older age, history of elective percutaneous coronary intervention, and higher baseline VO2peak significantly predicted ENR.22
In contrast with the SMARTEX-HF findings,19 a significant association between age, baseline LVEF, NYHA class, and exercise response was not found in the present study, possibly due to its smaller sample size.
Epigenetics to predict exercise response
In the past two decades, extensive research on (epi)genetic biomarkers of aerobic capacity has been performed. The HERITAGE Family Study attributed 47% of the variability in VO2peak response to genetic factors, and Karvinen et al.23,24 confirmed that the genome plays a significant role in exercise participation with a narrow-sense heritability of physical activity estimated at 53%. Furthermore, the latter also suggested that genetic pleiotropy might partly explain the association between high physical activity, cardiorespiratory fitness and survival.24 As epigenetic regulators of exercise response, several miRNAs have been investigated in healthy individuals or well-trained athletes. To the best of our knowledge, we are the first to describe a miRNA signature to identify non-responders to exercise training prior to exercise prescription. Plasma-derived miRNAs are stable and relatively easy detectable with conventional RT-qPCR techniques, and circulating miRNAs are withheld as promising diagnostic or prognostic biomarkers in different pathologies such as cancer and HF. Moreover, their role as therapeutic targets is currently being investigated in various clinical conditions.25,26
The proposed miRNA panel, consisting of Let-7b, miR-23a, miR-140, miR-146a, miR-191, miR-210, and miR-339-5p, is up-regulated in patients with an unfavourable VO2peak response. All miRNAs correlated with the change in VO2peak and showed good performance for the prediction of non-response (AUCs ≥ 77%). Circulating miR-210 has been inversely related to aerobic capacity in healthy subjects,9 and plasma miR-146a has been positively correlated with VO2peak in endurance athletes.10 Both acute and chronic exercise also up-regulate plasma levels of miR-146a in patients with chronic kidney disease and in healthy athletes.10,13
Epigenetics to identify underlying mechanisms of exercise response
Since the expression profiles of the proposed miRNA panel were significantly different between ER and ENR, miRNAs could also aid in unravelling the underlying mechanisms of the lack of response to exercise training in HFrEF patients. miRNAs can be found in the circulation in resting conditions, but also upon release by skeletal muscle or endothelial cells during exercise, thereby mediating exercise adaptation pathways. Previous research demonstrated that an increase in plasma volume with endurance training contributes to the increase in VO2peak.27 To unravel epigenetic mechanisms of exercise response, we performed an in silico analysis that revealed several pathways involved in exercise adaptation processes such as VEGF, Notch, apelin, MAPK, NOD-like and toll-like receptor pathways, providing potential therapeutic targets which may contribute to achieving higher VO2peak response rates in the future.
Endothelial function and angiogenesis
Both miR-146a and miR-210 stimulate VEGF expression by decreasing the expression of neurofibromin 2 and increasing Ras-related C3 botulinum toxin substrate 1 and p21-activated kinase-1, and by decreasing runt-related transcription factor-3 respectively.28,29 Furthermore, the proposed miRNA panel is involved in angiogenesis, with clear roles for miR-23a (inhibits Semaphorin-6A and -6D, and sprouty homolog 2),30 miR-146a (increases fibroblast growth factor binding protein 1 expression31), and miR-191 [targets nuclear factor kappa B (NF-kB) signalling and vascular endothelial zinc finger 1].32 Since these miRNA appeared to be up-regulated in ENR, either the stimulation of endothelial function and angiogenesis is not sufficient or downstream pathways are impaired, which could contribute to ENR.
Skeletal muscle mass and function
Some of the miRNAs have been related to skeletal muscle; Let-7b and miR-191-5p regulate myogenesis through respectively paired box protein-7 inhibition33 and MAPK, interleukin-6 signalling, and serine/threonine-protein phosphatase PP1-β catalytic subunit, and signal transducer and activator of transcription 3 pathway.34 Furthermore, miR-140 and miR-23 protect against skeletal muscle atrophy through inhibiting Wnt family member 11 expression35 and the ubiquitin-proteasome pathway respectively,36 and miR-146a-5p has been related to TGF-β, which is one of the key pathways in cardiac remodelling and fibrosis and in skeletal muscle repair after exercise.37 Additionally, miR-23a and miR-140-3p mediate cardiac hypertrophy, through targeting the ubiquitin-proteasome pathway and GATA binding protein 4, respectively.38,39
Inflammation
Several of the proposed miRNAs have been related to inflammatory processes; overexpression of miR-23a could down-regulate heat shock protein 90 and NF-kB protein in inflammatory macrophages and foam cells.40 miR-146a suppresses inflammation through increasing erb-b2 receptor tyrosine kinase 4 expression and decreasing tumour necrosis factor receptor-associated factor 6, interleukin-1 receptor-associated kinase 1, NF-kB, early growth response factor 1 expression and toll-like receptor 4 activation,41 and miR-339- targets fibroblast growth factor receptor substrate 2.42
Future perspectives
The proposed miRNA signature for identification of ENR is promising. Whether the high expression of the 7-miRNA panel in ENR reflects an increased release in the circulation or a compensatory rise to the lack of downstream response, remains to be elucidated before miRNA can assist in patient-tailored exercise prescription or can lead towards new therapeutic targets. Therefore, these findings deserve validation in a larger independent cohort.
Limitations
The findings of this retrospective study should be validated in a large prospective trial. A second limitation relates to the analysis of the miRNA expression data in relation to the constructed interaction networks in order to define biological themes, as the number of genes and interactions were limited to only those relevant in tissues of interest. Hence, we emphasize that the presented computational analysis is hypothesis-generating.
Conclusion
In HFrEF patients participating in exercise training, we investigated whether circulating miRNAs are able to predict VO2peak response. In this discovery cohort, we found a significantly different expression pattern of baseline plasma miRNA levels between ENR and ER, which could distinguish ENR from ER with 83% accuracy in an unsupervised analysis and 89% in a supervised analysis. A fingerprint of seven miRNAs was strongly correlated with %changeVO2peak and showed AUCs of ≥0.77 for predicting ENR. Pathway analysis revealed several targets involved in exercise adaptation processes such as angiogenesis, skeletal muscle adaptation, and inflammation. Therefore, the proposed miRNA panel can be an asset in optimizing personalized exercise prescription and could even open new therapeutic avenues in HF. Therefore, these findings deserve prospective validation in a large cohort of HFrEF patients.
Supplementary material
Supplementary material is available at European Journal of Preventive Cardiology online.
Poster presentation at the Heart Failure: Crossing the Translational Divide, Keystone, CO, USA, 14–18 January 2018.
Poster presentation at the Heart Failure Association Winter Meeting, Les Diablerets, Switzerland, 25 January 2018.
Poster presentation at the Heart Failure 2018—Annual congress of the ESC Heart Failure Association, Vienna, Austria, 26–29 May 2018.
Acknowledgements
The authors would like to thank all the participants in this study and the staff of the Cardiac Rehabilitation Centre of the Antwerp University Hospital.
Funding
This work was supported by the Flanders Research Foundation [Predoctoral mandate to I.W. 1194918N, senior clinical investigator grant to E.M.V.C. 1804320N] and the King Baudouin Foundation.
Conflict of interest: none declared.
Data availability
We declare that the materials described in the manuscript, including all the relevant raw data, will be available to any scientist wishing to use them for non-commercial purposes, without breaching participant confidentiality. All data are available from the corresponding author upon request.
Comments