Exploring Targets for Heart Failure Using Mendelian Randomization, Colocalization, Summary Data-based Mendelian Randomization, and Drug Prediction

Summary

This study identified eight potential drug targets for heart failure through genetic analysis and explored candidate therapeutic compounds associated with these identified targets.

Abstract

Heart failure (HF) is a prevalent cardiovascular disease that significantly impairs quality of life in its advanced stages. Despite a variety of current treatment strategies, the disease burden of HF remains significant. Therefore, exploring novel therapeutic targets is urgently needed. This study utilized Mendelian randomization (MR) to evaluate the causal relationships between druggable genes and HF. Colocalization analysis and summary data-based Mendelian randomization (SMR) were conducted to further validate the relationship between them. Finally, 8 potential therapeutic targets for HF were identified, including 4 risk factors (CYP11A1, GALT, KCNH2, and METRN) and 4 protective factors (APOM, CHD4, IL11RA, and LPAR5). The GO enrichment, KEGG enrichment, and protein-protein interaction (PPI) network analysis indicated that these genes were primarily involved in metabolic regulation and cardiac electrophysiological processes. Moreover, potential drugs targeting these targets were identified using drug prediction and molecular docking, including mitotane, Adehl, Benzofurans, 1,3-Di-o-tolylguanidine, 96-69-5, and bromoenol lactone. Furthermore, PCR results also demonstrated that mitotane might be associated with CYP11A1 and KCNH2. Drugs designed based on these potential therapeutic targets may offer higher success rates and improve clinical outcomes.

Introduction

Heart failure (HF) is a prevalent cardiovascular syndrome defined by impaired ventricular filling or myocardial contractility, resulting in inadequate cardiac output to meet metabolic demands^1,2,3. Clinically, HF is characterized by pulmonary and systemic congestion, tissue hypoperfusion, and symptoms such as dyspnea, fatigue, and peripheral edema^1,4,5. Despite therapeutic advances, HF remains a leading cause of morbidity and mortality worldwide, affecting over 64 million individuals according to the Global Burden of Disease (GBD) study^6,7,8. The substantial disease burden underscores the urgent need for novel therapeutic strategies.

Current pharmacologic and device-based interventions for HF are limited by suboptimal efficacy and adverse effects, contributing to persistently high rehospitalization rates and impaired quality of life^9,10,11,12. Identifying new molecular targets with causal links to HF pathogenesis is, therefore, critical for improving clinical outcomes. Integrating genomics into drug discovery represents a highly effective strategy for improving efficiency, as genetically validated therapies demonstrate a significantly higher success rate in clinical trials^13,14,15. Furthermore, proteins encoded by these actionable genes serve as prime therapeutic targets for both small molecules and monoclonal antibodies^16,17.

Recent advances in human genetics have facilitated the prioritization of drug targets through integrative analysis of genome-wide association studies (GWAS) and expression quantitative trait loci (eQTL) data. Mendelian randomization (MR) leverages genetic variants as instrumental variables to estimate causal effects, minimizing confounding and reverse causation inherent in observational studies^13,14. When combined with colocalization and summary-based MR (SMR), this approach enables robust prioritization of genes whose modulation may influence disease risk, offering a genetically anchored strategy for target discovery^18,19,20. Compared with conventional association studies, integrating MR, SMR, and colocalization overcomes confounding and reverse causality to establish robust causal inferences. Furthermore, this framework filters out false-positive signals from linkage disequilibrium by ensuring that gene expression and heart failure risk share the same underlying causal variant, optimizing target prioritization. It is worth noting that while cardiovascular tissues exhibit high specificity, this study primarily utilizes large-scale blood cis-eQTL data to maximize statistical power. Although this approach may not fully capture tissue-specific mechanisms inherent to HF, blood-derived variants serve as accessible biomarkers and provide robust insights into systemic genetic drivers of the disease.

In this study, large-scale GWAS and eQTL datasets were used to evaluate the causal effects of druggable genes on HF risk. Significant associations were further validated through colocalization and SMR analyses. To elucidate the biological functions of prioritized targets, Gene Ontology (GO) and Kyoto Encyclopedia of Genes and Genomes (KEGG) enrichment analyses, along with protein–protein interaction (PPI) network analysis, were performed. Therapeutic potential was subsequently assessed through drug prediction and molecular docking. Finally, the associations between selected druggable genes and candidate drugs were validated through in vitro experiments.

In conclusion, this study offers important insights for the discovery of new therapeutic targets for HF. By employing MR, SMR, colocalization analysis, gene enrichment analysis, PPI network construction, drug prediction, molecular docking, and in vitro experiments, these findings may contribute to the development of more effective therapeutic strategies for HF.

Protocol

All datasets used in this study were obtained from previously published studies with ethical approval. The reagents and the equipment used are listed in the Table of Materials.

1. Study design

The workflow of this study is illustrated in Figure 1. First, exposure data were obtained from the cis-eQTLs of druggable genes. Second, MR analysis was conducted to investigate the causal relationship between druggable genes and HF. Candidate therapeutic targets for HF were subsequently identified by colocalization and SMR analyses. Gene enrichment analysis and PPI network construction were then performed. Finally, drug prediction and molecular docking were conducted to explore the affinity of drugs for genes.

2. Data source

A total of 6,888 druggable genes were identified from two sources. One dataset included 4,463 genes from a previously published study that integrated multiple data sources¹⁷. The other dataset contained 5,012 genes from the Drug-Gene Interaction Database (DGIdb)²¹. Detailed information was illustrated on the website (https://www.dgidb.org/). After merging and removing duplicates, 6,888 unique druggable genes were retained for analysis.

Cis-eQTL data were obtained from the eQTLGen Consortium (https://eqtlgen.org/cis-eqtls.html), which includes expression profiles from 31,684 blood samples and covers 16,987 genes²².

HF summary statistics were obtained from the R12 release of the FinnGen database²³, which included 37,653 cases and 462,695 controls. The diagnosis of HF was based on the International Classification of Diseases (ICD), including ICD-10—I11.0, I13.0, I13.2, I50, ICD-9—4029B|428, and ICD-8—42700|42710|428|7824. All participants were of European ancestry. Detailed information is available on the FinnGen website and the original article (https://www.finngen.fi/en/access_results).

3. Mendelian randomization analysis

Instrumental variables (IVs) were selected based on three core MR assumptions²⁴. Firstly, IVs are directly associated with exposure factors. Secondly, IVs are not associated with confounding factors. Finally, IVs affect outcome through exposure only. Given that proximal eQTLs exert a more immediate and robust regulatory control over target genes, MR exclusively utilized SNPs located within 100 kb regions of the genes to maximize biological proximity. Strict criteria were employed for SNP selection. All selected SNPs were located within 100 kb upstream of the transcription start site and 100 kb downstream of the transcription end site of each druggable gene. Only SNPs with genome-wide significance (p < 5×10^-8) were included²⁵. To minimize linkage disequilibrium (LD), LD clumping was performed using an r² threshold of 0.001 and a clumping distance of 10,000 kb²⁶. SNPs with F-statistics >10 were retained, calculated using the formula: F = [R²(N−k−1)]/[k(1−R²)]²⁷. SNP information is provided in Supplementary Table 1.

The primary analytical methods were inverse-variance weighted (IVW) and Wald ratio. A p-value less than 0.05 was considered statistically significant^28,29. False discovery rate (FDR) was applied to adjust the p-value to minimize false positive results³⁰. When only one SNP was available, the Wald ratio method was used; otherwise, IVW was applied. Analyses included MR-Egger, weighted median, weighted mode, and simple mode when enough SNPs were present.

The Cochran's Q test and MR-Egger intercept analysis were employed to test heterogeneity and pleiotropy. A Q-p-value greater than 0.05 indicated the absence of heterogeneity³¹. A p-value for the MR-Egger intercept less than 0.05 indicated the presence of directional pleiotropy¹⁸. Finally, the leave-one-out test was conducted to test if the results could be influenced by a single SNP. MR analyses were performed using the TwoSampleMR package (version 0.6.11) in R (version 4.4.2).

4. Colocalization analysis

Colocalization analysis was conducted using the coloc R package (version 6.0.1) to identify whether shared genetic variants influenced both gene expression and HF risk³². Prior probabilities for an SNP being associated with either trait and both traits were set to 1 × 10^-4 and 1 × 10^-5, respectively. The posterior probabilities of the five hypotheses can be calculated: (1) PPH0: SNPs are not linked with two traits. (2) PPH1/PPH2: SNPs are associated with gene expression or HF. (3) PPH3: SNPs are linked with HF and gene expression, but they are driven by different SNPs. (4) PPH4: SNPs are linked with HF and gene expression, and they are driven by the common SNPs. A PPH4 > 0.8 was considered evidence for colocalization.

5. SMR analysis

The SMR analysis was carried out to explore the relationship between the gene and HF using the SMR software (version 1.3.1)³³. The heterogeneity in dependent instruments (HEIDI) test was carried out to determine whether observed associations were due to linkage. A significant SMR p-value (< 0.05) combined with a HEIDI p-value > 0.05 was interpreted as evidence supporting a true causal relationship rather than a linkage scenario-caused association. The detailed operational protocols and computational settings can be accessed on the official website.

6. Enrichment analysis

GO and KEGG enrichment analysis was employed to characterize the biological functions and pathways of the screened potential therapeutic target genes³⁴. GO enrichment included biological processes (BP), molecular functions (MF), and cellular components (CC). The R package clusterProfiler (version 4.14.6) and Pathview (version 1.46.0) were utilized to conduct the enrichment analysis^35,36.

7. Candidate drug prediction

To identify potential therapeutic compounds, drug prediction analysis was performed using the Drug Signatures Database (DSigDB) (http://dsigdb.tanlab.org/DSigD Bv1.0/)³⁷. 22,527 gene sets are included in the DSigDB, including 17,389 unique compounds covering 19,531 genes. Associations between candidate genes and compounds were analyzed, followed by enrichment analysis to identify potentially relevant drugs targeting HF-associated genes.

8. Protein interaction network construction

Protein-protein interaction (PPI) networks were constructed using GeneMANIA (https://genem ania.org/)³⁸. GeneMANIA could help to generate gene function hypotheses, analyze gene lists, and prioritize genes for functional testing.

9. Molecular docking

Molecular docking was performed to validate the druggability of the genes and their relationship to the drug candidates. Docking simulations help assess binding affinity and interaction patterns between compounds and target proteins, thereby informing candidate prioritization and optimization. Protein structures were retrieved from the Protein Data Bank (PDB) (http://www.rcsb.org/), and compound structures were obtained from PubChem (https://pubchem.ncbi.nlm.nih.gov/)³⁹. Docking was performed using CB-Dock2 (https://cadd.labshare.cn/cb-dock2/index.php), which uses cavity detection to guide blind docking^40,41. Detailed operation procedures and calculation settings information could be obtained on the official website.

10. Cell culture

Human immortalized cardiomyocyte cells (AC16) were cultured in AC16 cell-specific culture medium under sterile conditions. Cells were maintained at 37 °C in a humidified incubator containing 5% CO2. The culture medium was replaced every 2–3 days, and cells were passaged at approximately 70%–80% confluence to maintain optimal growth conditions. Cell morphology and confluence were monitored routinely using light microscopy before each experiment. All cell culture procedures were performed in a biosafety cabinet to minimize contamination risk. For treatment experiments, AC16 cells were seeded into 6-well plates and allowed to adhere overnight before treatment administration.

11. Quantitative reverse transcription PCR (RT-PCR)

AC16 cells were treated with 200 µM palmitic acid (PA) or 20 µM mitotane for 24 h after reaching approximately 60%–70% confluence. PA was dissolved in bovine serum albumin (BSA), whereas mitotane was dissolved in dimethyl sulfoxide (DMSO). All potentially hazardous reagents were handled in accordance with standard laboratory safety procedures using appropriate personal protective equipment, including gloves and laboratory coats.

Following treatment, total RNA was extracted using an RNA extraction kit according to the manufacturer’s instructions under RNase-free conditions. RNA extraction procedures were performed on ice whenever possible to minimize RNA degradation. RNA concentration and purity were assessed before downstream experiments, and only samples with an A260/A280 ratio between 1.8 and 2.0 were used for subsequent analyses.

Complementary DNA (cDNA) was synthesized using a reverse transcription kit according to the manufacturer’s protocol. Quantitative RT-PCR was subsequently performed using SYBR Green Master Mix on a real-time PCR detection system. Each reaction was prepared in a final volume of 10 µL containing 5 µL SYBR Green Master Mix, 0.4 µL forward primer, 0.4 µL reverse primer, 1 µL cDNA template, and 3.2 µL nuclease-free water. Amplification was performed under the following cycling conditions: initial denaturation at 95 °C for 30 s, followed by 40 cycles of denaturation at 95 °C for 5 s and annealing/extension at 60 °C for 30 s. Each experiment was performed with three biological replicates, and all reactions were conducted in triplicate.

GAPDH was used as the internal reference gene, and the relative mRNA expression levels of CYP11A1 and KCNH2 were calculated using the 2^−ΔΔCt method. A single peak in melt curve analysis was considered indicative of specific amplification. The primer sequences used for RT-PCR are listed in Supplementary Table 2.

Results

Mendelian randomization
The MR results evaluating the causal effects between druggable gene expression and HF were presented in Supplementary Table 3. An FDR-adjusted P value less than 0.05 was considered statistically significant. A total of 11 genes with significant causal associations with HF were identified (Figure 2). Among them, six genes were associated with an increased risk of HF, including CDKN1A (P_FDR: 0.03, OR: 1.293), CYP11A1 (P_FDR: 0.034, OR: 1.494), GALT (P_FDR: 0.042, OR: 1.074), KCNH2 (P_FDR: 0.018, OR: 1.163), LST1 (P_FDR: 0.03, OR: 1.112) and METRN (P_FDR: 0.018, OR: 1.115). The remaining five genes were associated with a decreased risk of HF, including APOM (P_FDR: 0.03, OR: 0.773), CAMK2G (P_FDR: 0.018, OR: 0.699), CHD4 (P_FDR: 0.018, OR: 0.811), IL11RA (P_FDR: 0.018, OR: 0.947), and LPAR5 (P_FDR: 0.021, OR: 0.915).

Sensitive analysis was performed in this study, and the results showed that no significant heterogeneity or horizontal pleiotropy was detected (Table 1), and the MR estimates were not driven by any single SNP (Supplementary Figure 1).

Colocalization analysis
Colocalization analysis was performed on the 11 identified genes (Figure 3 and Table 2). Nine genes demonstrated strong evidence of sharing causal variants with HF, including APOM (PP.H4: 0.923), CHD4 (PP.H4: 0.942), CYP11A1 (PP.H4: 0.91), GALT (PP.H4: 0.967), IL11RA (PP.H4: 0.967), KCNH2 (PP.H4: 0.815), LPAR5 (PP.H4: 0.938), LST1 (PP.H4: 0.913) and METRN (PP.H4: 0.907). There were two genes that were excluded from further analysis because the PP.H4 was less than 0.8, namely, CDKN1A and CAMK2G.

SMR analysis
To further validate gene-disease associations, SMR analysis was conducted for the nine colocalized genes (Table 2). The SMR analysis manifested that eight genes remained significantly associated with HF. There were four genes negatively associated with HF, including APOM (P_SMR: 0.018, P_HEIDI: 0.115, OR: 0.838), CHD4 (P_SMR: 5.45E-05, P_HEIDI: 0.458, OR: 0.732), IL11RA (P_SMR: 1.00E-04, P_HEIDI: 0.388, OR: 0.938), and LPAR5 (P_SMR: 3.71E-05, P_HEIDI: 0.168, OR: 0.889). And four were positively associated with HF, including CYP11A1 (P_SMR: 5.68E-04, P_HEIDI: 0.056, OR: 1.81), GALT (P_SMR: 1.02E-04, P_HEIDI: 0.726, OR: 1.085), KCNH2 (P_SMR: 2.76E-04, P_HEIDI: 0.162, OR: 1.236) and METRN (P_SMR: 1.55E-05, P_HEIDI: 0.068, OR: 1.103). Although LST1 (P_SMR: 4.61E-05, P_HEIDI: 0.025, OR: 1.135) showed an association with increased risk of HF, the HEIDI test indicated that the association might be caused by linkage disequilibrium. Therefore, LST1 was excluded.

Enrichment analysis
The eight genes identified through MR, colocalization, and SMR analyses were subjected to GO enrichment and KEGG enrichment analysis (Figure 4A–E). In class BP, genes were enriched in disassembly of metabolic processes (cortisol metabolic process, galactose metabolic process, glucocorticoid biosynthetic process, and reverse cholesterol transport), lipoprotein particle-related processes (high-density lipoprotein particle assembly, high-density lipoprotein particle remodeling, and high-density lipoprotein particle clearance), and membrane repolarization during ventricular cardiac muscle cell action potential. In class CC, genes were enriched in lipoprotein particles (low-density lipoprotein particle, very-low-density lipoprotein particle, triglyceride-rich plasma lipoprotein particle, high-density lipoprotein particle, plasma lipoprotein particle, and lipoprotein particle) and protein complexes (NuRD complex, CHD-type complex, and protein-lipid complex). In class MF, genes were enriched in electrophysiological processes (voltage-gated potassium channel activity involved in ventricular cardiac muscle cell action potential repolarization, voltage-gated potassium channel activity involved in cardiac muscle cell action potential repolarization, delayed rectifier potassium channel activity, and inward rectifier potassium channel activity) and nuclear structure binding (nucleosomal DNA binding, scaffold protein binding, and nucleosome binding). As shown in Figure 4D,E, the pathways analyzed by KEGG enrichment were metabolic processes (galactose metabolism, biosynthesis of nucleotide sugars, amino sugar and nucleotide sugar metabolism, steroid hormone biosynthesis and cortisol synthesis, and secretion) and the prolactin signaling pathway. The GO enrichment and KEGG enrichment network analysis indicated that these genes were primarily involved in metabolic regulation and cardiac electrophysiological processes.

Candidate drug prediction
Drug prediction of druggable genes was performed through the DSigDB database. Based on the adjusted p-values, the top ten potential drug candidates were shown in Table 3 and Figure 4F,G. The results showed that the top ten drug candidates primarily targeted three genes. Mitotane was the most significant drug associated with KCNH2 and CYP11A1. In addition, CYP11A1 also interacted with mitotane, vanadic sulfate, adehl, 2,4-DIBROMOPHENOL, 3,3',4,4'-Tetrabromobiphenyl, 3,3',4,4',5,5'-Hexabromobiphenyl, and benzofurans. Moreover, GALT only interacted with vanadic sulfate.

Protein interaction network construction
A protein-protein interaction (PPI) network was constructed based on the eight druggable genes (Figure 5). The PPI network included 20 additional interacting proteins. The interactions were primarily based on co-expression (71.28%), colocalization (28.26%), predicted (5.37%), and genetic interactions (0.46%). The network highlighted hormone-related metabolic processes as the most enriched functional categories, including glucocorticoid metabolic process, steroid hormone biosynthetic process, mineralocorticoid metabolic process, hormone biosynthetic process, cellular hormone metabolic process, steroid biosynthetic process, and hormone metabolic process.

Molecular docking
Molecular docking was employed to assess the binding affinities of the top ten predicted drugs to their target proteins. A binding energy less than -5 kcal/mol suggested that a stable bond may exist between the drug and the target. The information on the potential drugs is shown in Supplementary Table 4. The results showed six drugs with potentially effective docking profiles (Table 4 and Figure 6). Drug candidates that formed hydrogen bonds and exhibit electrostatic interactions with their molecular targets had the potential to effectively fit into the binding sites of these targets. The lowest binding energy was observed between CYP11A1 and adehl (−6.3 kcal/mol), and between KCNH2 and bromoenol lactone (−6.3 kcal/mol), suggesting stable binding configurations.

Quantitative reverse transcription PCR (RT-PCR)
The result has shown that mitotane was the most significant potential drug associated with KCNH2 and CYP11A1. The in vitro experiment was utilized to test the relationship. The PCR results showed that both CYP11A1 and KCNH2 were downregulated after mitotane treatment. Due to the absence of a cellular model for HF, AC16 cells were stimulated with PA based on the pathophysiological mechanisms of HF, and both CYP11A1 and KCNH2 were found to be upregulated upon PA stimulation. Meanwhile, when mitotane was added concurrently with PA stimulation, both CYP11A1 and KCNH2 were downregulated (Figure 7).

DATA AVAILABILITY:
The corresponding source links for the raw data and accession numbers are provided in Supplementary File 1.

Figure 1: Schematic overview of the study design. Abbreviation: eQTLs, expression quantitative trait loci; HF, heart failure; MR, Mendelian randomization; SMR, summary-data-based Mendelian randomization; HEIDI, heterogeneity in dependent instruments; GO, gene ontology; KEGG, Kyoto Encyclopedia of genes and genomes; PPI, protein-protein interaction. Please click here to view a larger version of this figure.

Figure 2: Forest plot illustrating the MR results. Abbreviation: nsnp, the number of snp; OR, odds ratio. Please click here to view a larger version of this figure.

Figure 3: Results of the colocalization analysis. (A) APOM, (B) CHD4, (C) CYP11A1, (D) GALT, (E) IL11RA, (F) KCNH2, (G) LPAR5, (H) LST1, and (I) METRN. Please click here to view a larger version of this figure.

Figure 4: Functional enrichment and candidate drug prediction results. (A) Circos plot of GO enrichment analysis; (B) Bar plot of GO terms; (C) Bubble plot of GO terms; (D) Bar plot of KEGG pathways; (E) Bubble plot of KEGG pathways; (F) Bar plot of candidate drug prediction results; (G) Bubble plot of candidate drug prediction results. Please click here to view a larger version of this figure.

Figure 5: PPI network constructed using GeneMANIA. Each node represents a gene, and color indicates the biological process associated with each gene. Please click here to view a larger version of this figure.

Figure 6: Molecular docking results of drug-target interactions. (A) CYP11A1 with mitotane; (B) CYP11A1 with Adehl; (C) CYP11A1 with benzofurans; (D) KCNH2 with mitotane; (E) KCNH2 with Adehl; (F) KCNH2 with 1,3-di-o-tolylguanidine; (G) KCNH2 with 96-69-5; (H) KCNH2 with bromoenol lactone. Please click here to view a larger version of this figure.

Figure 7: PCR results of CYP11A1 and KCNH2 following mitotane treatment. (A,B) The mRNA expression of CYP11A1 and KCNH2 treated with mitotane for 24 h; (C,D) The mRNA expression of CYP11A1 and KCNH2 treated with PA for 24 h; (E,F) The mRNA expression of CYP11A1 and KCNH2 treated with mitotane plus PA for 24 h (N = 3). The data are presented as the mean ± standard deviation. *p < 0.05; **p < 0.01; ***p < 0.001. Please click here to view a larger version of this figure.

Exposure	nSNP	Method	Q	Q_pval	Egger_intercept	Pval
CDKN1A	2	IVW	4.466	0.035	NA	NA
GALT	2	IVW	0.036	0.849	NA	NA
IL11RA	5	IVW	1.742	0.783	-0.034	0.371
KCNH2	2	IVW	0.408	0.523	NA	NA
LPAR5	2	IVW	0.647	0.421	NA	NA
LST1	4	IVW	3.155	0.368	-0.127	0.369

Table 1: The result of sensitive analysis. Abbreviation: nSNP: the number of SNPs; IVW: inverse variance weighted; NA: not available.

Exposure	PP.H4	β_SMR	OR_SMR	p_SMR	p_HEIDI
APOM	0.923	-0.176	0.838	0.018	0.115
CAMK2G	0.517	-0.527	0.59	NA	NA
CDKN1A	0.172	0.427	1.533	NA	NA
CHD4	0.942	-0.311	0.732	5.45E-05	0.458
CYP11A1	0.91	0.593	1.81	5.68E-04	0.056
GALT	0.967	0.081	1.085	1.02E-04	0.726
IL11RA	0.967	-0.064	0.938	1.00E-04	0.388
KCNH2	0.815	0.212	1.236	2.76E-04	0.162
LPAR5	0.938	-0.117	0.889	3.71E-05	0.168
LST1	0.913	0.126	1.135	4.61E-05	0.025
METRN	0.907	0.098	1.103	1.55E-05	0.068

Table 2: The result of colocalization analysis and SMR by using the cis-eQTL dataset. Abbreviation: SMR: summary-data-based mendelian randomization; HEIDI: heterogeneity in dependent instruments; OR: odds ratio, NA: not available.

Drug	gene	Pvalue	P.adjust
mitotane	KCNH2/CYP11A1	<0.001	0.028
Vanadic sulfate	CYP11A1/GALT	0.002	0.028
Adehl	KCNH2/CYP11A1	0.002	0.028
1,3-Di-o-tolylguanidine	KCNH2	0.004	0.028
2,4-DIBROMOPHENOL	CYP11A1	0.004	0.028
3,3',4,4'-Tetrabromobiphenyl	CYP11A1	0.004	0.028
3,3',4,4',5,5'-Hexabromobiphenyl	CYP11A1	0.004	0.028
96-69-5	KCNH2	0.004	0.028
Benzofurans	CYP11A1	0.004	0.028
bromoenol lactone	KCNH2	0.004	0.028

Table 3: The result of candidate drug prediction.

Drug	Target	Binding energy
mitotane	CYP11A1	-5
Vanadic sulfate	CYP11A1	-3.2
Adehl	CYP11A1	-6.3
2,4-DIBROMOPHENOL	CYP11A1	-4
3,3',4,4'-Tetrabromobiphenyl	CYP11A1	/
3,3',4,4',5,5'-Hexabromobiphenyl	CYP11A1	/
Benzofurans	CYP11A1	-6.2
Vanadic sulfate	GALT	-4.2
mitotane	KCNH2	-5.1
Adehl	KCNH2	-6
1,3-Di-o-tolylguanidine	KCNH2	-6.1
96-69-5	KCNH2	-5.6
bromoenol lactone	KCNH2	-6.3

Table 4: The result of molecular docking.

Supplementary Figure 1: The plot of the leave-one-out test.Please click here to download this file.

Supplementary Table 1: Summary information of the SNPs used in the Mendelian randomization analysis.Please click here to download this file.

Supplementary Table 2: Primer sequences used for RT-qPCR.Please click here to download this file.

Supplementary Table 3: Mendelian randomization results for 3,721 druggable genes in relation to heart failure.Please click here to download this file.

Supplementary Table 4: The information on the potential drugs.Please click here to download this file.

Supplementary File 1: Detailed list of data sources in this study.Please click here to download this file.

Discussion

This research systematically evaluated the causal effects of druggable genes on HF risk by using MR, identifying eight potential therapeutic targets. These included four protective factors (APOM, CHD4, IL11RA, and LPAR5) and four risk factors (CYP11A1, GALT, KCNH2, and METRN). These associations were further supported by colocalization analysis and SMR. Functional annotation through GO, KEGG enrichment, and PPI network analyses suggested that these genes may influence HF primarily through pathways related to metabolism and electrophysiology. Moreover, drug prediction and molecular docking studies highlighted the translational potential of these targets.

APOM was identified as a protective factor for HF. APOM binds sphingosine-1-phosphate (S1P) and exerts anti-inflammatory effects by reducing TNF-α expression, limiting monocyte-endothelial adhesion, and maintaining endothelial barrier integrity, thereby mitigating inflammation-induced myocardial injury^42,43,44. Additionally, APOM supports cardiomyocyte survival and proliferation while reducing apoptosis, preserving cardiac structure and function. Its deficiency has been linked to adverse cardiac outcomes, consistent with this finding^42,45. However, the underlying mechanism linking APOM to HF requires further investigation.

CHD4, another protective gene identified, is a core component of the nucleosome remodeling and deacetylase (NuRD) complex and plays an essential role in cardiac development and myocardial growth^46,47,48. Loss of CHD4 leads to incorporation of inappropriate myofibrillar proteins into cardiac sarcomeres, resulting in structural and functional abnormalities, including chamber dilation and impaired ventricular development. This result reinforces the critical role of CHD4 in maintaining normal cardiac function.

It was also revealed that IL11RA could decrease the risk of HF. IL11RA is the primary ligand-binding subunit of the IL-11 receptor complex. While IL-11 has been shown to drive fibrotic responses in the heart, the specific role of IL11RA remains unclear⁴⁹. The protective effect observed in this study suggests that IL11RA may exert context-dependent functions or interact with other signaling components to limit fibrosis under certain conditions⁵⁰. This apparent discrepancy highlights the need for further mechanistic studies to clarify the functional role of IL11RA in HF pathogenesis.

This study showed that LPAR5 was associated with a reduced risk of HF. LPAR5 is a receptor for lysophosphatidic acid (LPA). Although LPA is generally considered pro-inflammatory and has been implicated in adverse cardiac remodeling, individual LPA receptors may have divergent roles^51,52. For instance, LPA signaling modulates cytoskeletal stress fiber formation through Gα12/13 and Gαq/11 pathways, playing an important role in maintaining cardiomyocyte morphology, contractile function, and signal transduction^53,54. Nevertheless, the precise role of LPAR5 in HF pathophysiology warrants further investigation.

This investigation indicated that CYP11A1 was a risk factor for HF. CYP11A1 encodes the rate-limiting enzyme in steroid hormone biosynthesis^55,56. Elevated production of steroid hormones, particularly aldosterone, contributes to sodium retention and increased cardiac workload, key features of HF progression^57,58. Additionally, CYP11A1 generates reactive oxygen species during its metabolic activity, which may exacerbate myocardial damage⁵⁹. A real correlation between CYP11A1 and mitotane was observed in vitro, providing evidence that mitotane may have potential therapeutic value for the treatment of HF. However, further experiments are needed to explore the relationship between mitotane and HF.

GALT could increase the risk of HF in this research. GALT is a key enzyme in galactose metabolism and shows a positive association with HF risk^60,61. Elevated GALT expression has been proposed as a biomarker for cardiac fibrosis and remodeling, with higher levels linked to worse clinical outcomes, possibly due to its association with profibrotic metabolic pathways^62,63. This finding suggests that targeting GALT may improve HF prognosis; however, whether this is mediated through attenuation of fibrosis requires further investigation.

This research demonstrated that KCNH2 was a risk factor for HF. KCNH2 encodes a subunit of the voltage-gated potassium channel, which is involved in cardiac repolarization⁶⁴. Mutations in KCNH2 are known to cause arrhythmias, increasing myocardial stress and HF risk^65,66. Emerging evidence also suggests that KCNH2 may influence mitochondrial function. Given the established role of mitochondrial dysfunction in HF pathogenesis, the potential effects of KCNH2 on mitochondrial function in the context of HF represent an important direction for future research^67,68. The result of RT-PCR indicated that KCNH2 was associated with mitotane; further investigation was needed to explore the link between mitotane and HF.

Finally, METRN was identified as a risk factor, although its role in HF remains controversial. METRN is a secretory protein highly expressed in adipose tissue and mucosa, implicated in fat browning, insulin sensitivity, and inflammation^69,70. While some studies have reported protective effects through AMPK and PPARδ pathways, others have associated elevated plasma METRN levels with worse outcomes in HF patients^71,72. This finding supports the latter view, underscoring the need for further studies to clarify METRN’s function in HF.

This study has several strengths. First, MR was constructed to infer causal relationships between gene expression and HF, reducing potential confounding. FDR correction was applied to account for multiple testing. Colocalization and SMR analyses further strengthened the credibility of these findings. Functional enrichment and PPI network analyses provided mechanistic insights, while drug prediction and molecular docking highlighted the translational potential of identified targets. Furthermore, these experiments further suggested a relationship between CYP11A1, KCNH2, and mitotane, which supported the reliability of this study.

However, several limitations should be acknowledged. Firstly, all the analyses were based on genetic data from individuals of European ancestry, which may limit generalizability to other populations. Secondly, eQTL data were derived solely from blood samples, which may not fully capture gene expression patterns in cardiac tissue. Thirdly, smaller sample sizes in some datasets may have reduced statistical power, underscoring the need for larger studies. Fourthly, the cellular model used in this study was PA-stimulated AC16, which did not fully recapitulate the complex pathophysiological features of HF in vivo, necessitating further experimental investigation into the underlying mechanisms.

This study identified eight potential drug targets for HF by using MR, SMR and colocalization analyses. Drug prediction and molecular docking further identified some potential drugs for these candidate targets. These findings offer new promising candidates for future therapeutic development. However, further research is needed to investigate the feasibility and mechanism of these targets.

Materials

List of materials used in this article

Name	Company	Catalog Number	Comments
2× Universal SYBR Green Fast qPCR Mix	Abclonal	RK21203
AC16	HyCyte	TCH-C119
AC16 cell-specific culture medium	HyCyte	TCH-G119
cDNA reverse kit	Abclonal	RK20400
CFX96™ Real-Time system	Bio-Rad Laboratories
Mitotane	TargetMol	T1199
Palmitic acid	TargetMol	T2908
RNA extraction kit	Abclonal	RK30120