Genome-wide meta-analysis of variant-by-diuretic interactions as modulators of lipid traits in persons of European and African ancestry

Hypertension (HTN) is a significant risk factor for cardiovascular morbidity and mortality. Metabolic abnormalities, including adverse cholesterol and triglycerides (TG) profiles, are frequent comorbid findings with HTN and contribute to cardiovascular disease. Diuretics, which are used to treat HTN and heart failure, have been associated with worsening of fasting lipid concentrations. Genome-wide meta-analyses with 39,710 European-ancestry (EA) individuals and 9925 African-ancestry (AA) individuals were performed to identify genetic variants that modify the effect of loop or thiazide diuretic use on blood lipid concentrations. Both longitudinal and cross sectional data were used to compute cohort-specific interaction results, which were then combined through meta-analysis in each ancestry. These ancestry-specific results were further combined through trans-ancestry meta-analysis. Analysis of EA data identified two genome-wide significant (p < 5 × 10−8) loci with single nucleotide variant (SNV)-loop diuretic interaction on TG concentrations (including COL11A1). Analysis of AA data identified one genome-wide significant locus adjacent to BMP2 with SNV-loop diuretic interaction on TG concentrations. Trans-ancestry analysis strengthened evidence of association for SNV-loop diuretic interaction at two loci (KIAA1217 and BAALC). There were few significant SNV-thiazide diuretic interaction associations on TG concentrations and for either diuretic on cholesterol concentrations. Several promising loci were identified that may implicate biologic pathways that contribute to adverse metabolic side effects from diuretic therapy.


Introduction
Hypertension (HTN) is a significant risk factor for cardiovascular morbidity and mortality. However, even after accounting for the beneficial effects of blood pressure reduction on cardiovascular events, cardiovascular risk remains elevated despite intensive lowering of blood pressure [1].
Metabolic abnormalities, including adverse cholesterol and triglycerides (TG) profiles, are frequent comorbid findings with HTN and may contribute to the residual cardiovascular disease (CVD) risk associated with antihypertensive therapy. Thiazide diuretics are recommended as first-line therapy for the treatment of HTN [2]. Loop diuretics are a mainstay in the treatment of heart failure and are frequently used in patients with HTN complicated by renal insufficiency to control volume retention. However, there is little evidence to support a survival benefit with the use of loop diuretics [3]. Both thiazide and loop diuretics have been associated with worsening fasting lipid concentrations [4][5][6], which may dampen the salutary effects of antihypertensive therapy in some patients. While there have been studies identifying genetic sources of inter-individual blood pressure response to primarily thiazide diuretic therapy [7][8][9][10][11], there is little literature regarding the interactions of diuretics with genetic variants on blood lipid concentrations [12].
A meta-analysis of 56 randomized, placebo-controlled trials of diuretic monotherapy for the treatment of HTN, including both European-ancestry (EA) and African-ancestry (AA) individuals, showed that diuretic therapy was associated with average changes of 0.29 mmol/L for total cholesterol (TC), 0.24 mmol/L for low-density lipoprotein cholesterol (LDL), −0.02 mmol/L for high-density lipoprotein cholesterol (HDL), and 0.35 mmol/L for TG [13]. Some long-term studies have suggested that effects of diuretics on lipid concentrations waned after 1 year [14,15]. A retrospective analysis of the Systolic Hypertension in the Elderly Program showed persistent, albeit more modest, increases in lipid concentrations (3-year change in TG: 0.28 ± 0.86 mmol/L with chlorthalidone vs. 0.09 ± 0.71 mmol/L with placebo, P < 0.001) when a cohort was randomized to a thiazide-like diuretic or a placebo [16]; however, these results may have been affected by the use of nondiuretic antihypertensive drugs (atenolol and reserpine). Notably, the larger standard deviations suggest that an underlying genetic component may explain the high variability in individual's lipid response to diuretic therapy.
Although the genetic underpinnings of fasting cholesterol and TG concentrations have been well-described, single nucleotide variants (SNVs) identified through genome-wide association studies (GWAS) explain only~10-12% of heritability [17]. Some studies have shown that the response of lipid concentrations to diuretics is significantly greater in AA than in EA individuals [13], further suggesting the modulating effect of genetic architecture. Unfortunately, most long-term studies were based on either almost exclusively EA subjects [15] or studies including individuals of multiple ancestries where response by ancestry was not characterized [14]. The purpose of this investigation was to determine whether common genetic variants modify the effect of diuretic use on blood lipid concentrations in persons of EA and AA. Exposure to thiazide and loop diuretics are considered separately since their mechanisms action differ; thiazide diuretics inhibit the sodium-chloride cotransporter (SLC12A3) in the renal distal convoluted tubule whereas loop diuretics inhibit the sodium-potassium-chloride cotransporter (SLC17A3) in the thick ascending loop of Henle [18]. These analyses may identify biologic pathways that contribute to adverse metabolic side effects from diuretic therapy.

Materials and methods
Study samples, phenotype, and genotype data Data from 14 EA and 7 AA cohorts was used. Table 1 [20], Minimac [21], IMPUTE2 [22], or BEAGLE [23]. Information on genotype and imputation for each study is presented in Supplementary Table 2. In total, 49,635 subjects over 18 years of age with genotype, phenotype, and covariate information were available in this analysis. Both longitudinal and cross sectional studies were used, and all available data for each subject were used. In longitudinal cohorts that had multiple clinic visits for each subject, those multiple measurements (obtained across clinic visits) were used in the analysis. In cross sectional studies that had only a single visit for each subject, a single measurement was used. Phlebotomy was performed to measure blood lipids after a minimum 8-h fast as per individual study protocols. Three fasting lipid traits were considered for analyses: TG, HDL, and LDL (all mmol/L). LDL was calculated via the Friedewald equation (LDL = TC − HDL − [TG/5]) for those with TG ≤ 4.52 mmol/L [24]. If TG > 4.52 mmol/L, LDL was set to missing unless directly assayed. LDL concentrations were adjusted for statin use as described elsewhere [25]. TG concentrations were natural log-transformed for analysis.
Two patterns of diuretic exposure were used: (1) thiazide and thiazide-like diuretic use (yes/no) without use of loop; and (2) loop diuretic use (yes/no) without use of thiazide or thiazide-like diuretics. For each diuretic exposure, the unexposed group consisted of individuals taking neither diuretic. Drug use (yes/no) was defined as the use of any drug in the drug class, regardless of dosage, and with or without the use of a potassium supplement, as determined by self-report and/or inventory of medication bottles at each clinic visit. Examples of thiazide diuretic, loop diuretic, and statin medications are listed in Supplementary Table 3.

Statistical analyses
To determine the association of SNV-diuretic interactions on blood cholesterol and TG concentrations, the following regression model was applied. Y is the fasting lipid value, i is the diuretic use (coded 0/1 for the absence/presence of the diuretic), G is the dosage of the imputed genetic variant coded additively (from 0 to 2), and C is the vector of all other covariates, which include age, sex, body mass index (BMI), principal components (to account for population stratification and admixture), and additional cohort-specific covariates (if any). Subjects with missing data for fasting (≥8 h) lipid values, diuretic use, or any covariate were excluded from analysis. Subjects taking both thiazide and loop diuretics were also excluded from analysis.
Each study conducted GWAS analysis and provided the estimated SNV-diuretic interaction effect β GD and standard error. For longitudinal cohorts with repeated measures, generalized estimating equations (GEE) with an independent working correlation were used with an R package boss; for studies with cross sectional data, linear regression was used. To account for relatedness in families, family studies used either GEE treating each family as a cluster using R software or the linear mixed effect model approach with a random polygenic component (for which the covariance matrix depends on the kinship matrix) using ProbABEL (Supplementary Table 2). When minor allele counts of SNVs among participants on the diuretic were small (<10), the standardized interaction effect (β GD /SE) is not normally distributed. Note that, although a normal distribution is often appropriate in a large sample (and/or main effect analysis of GWAS), it is not appropriate in a GxE interaction study. Therefore, following our earlier work [26], cohort-specific P values were computed based on a Student's t-distribution. The degrees of freedom for the t-distribution depend on the number of drug-exposed participants (N_exposed), the SNV imputation quality, and the minor allele frequency. Before meta-analysis, SNVs were excluded if 2 • MAF • N_exposed • imputation quality measure <10 to exclude unstable cohort-specific results that reflect small sample size, low MAF, or low imputation quality measures. As shown in Sitlani et al. [27], we observed that this threshold provides a good balance between over-filtering and under-filtering.
Meta-analysis for each ancestry group was performed to combine these cohort-specific P values with weighted Zstatistics using METAL [28], where weights were based on the number of drug-exposed subjects multiplied by the SNV imputation quality. As cohort-specific P values were computed based on t distribution, inverse-variance weighted meta-analysis was not used. After meta-analysis, SNVs were excluded if heterogeneity P < 1 × 10 −6 or they were available in fewer than three cohorts for EA results or two cohorts for AA results. The ancestry-specific results were BMI body mass index (kg/m 2 ), DBP diastolic blood pressure (mmHg), SBP systolic blood pressure (mmHg), TG triglycerides (mmol/L) further combined to perform trans-ancestry meta-analysis using MANTRA (Meta-analysis of transethnic association studies) [29]. MANTRA accounts for similarity in allelic effects among closely related populations, while allowing for heterogeneity across populations with more diverse ancestries. As MANTRA uses a Bayesian framework, a traditional fixed-effect meta-analysis with weighted Z-statistics was also performed using METAL. Genome-wide significance was defined as P < 5 × 10 -8 from METAL with a fixed-effect meta-analysis or Bayes Factor > 10 6 from MANTRA. Suggestive evidence of association was defined as P < 1 × 10 -6 from METAL or Bayes Factor > 10 5 from MANTRA. To assess type I error due to population stratification and other factors, quantile-quantile (QQ) plots were examined for all cohortspecific GWAS results for each pair of lipid and diuretic exposure. In addition, during meta-analysis, genomic control correction [30] was applied to cohort-specific GWAS results if their genomic control lambda value was greater than 1. The gene locations referenced in the text and tables were obtained from the National Center for Biotechnology Information dbSNP database (reference assembly GRCh38. p2). Functional annotation information was sought using HaploReg [31] and RegulomeDB [32].

Results
The EA group included 39,710 subjects from 14 cohorts; the AA group included 9925 subjects from 7 cohorts ( Table 1). The number of subjects exposed to loop diuretics was 2117 (5.3%) in EA and 784 (7.9%) in AA; the number exposed to thiazide and thiazide-like diuretics was 6878 (17.3%) in EA and 3923 (39.5%) in AA.
The QQ plots (Supplementary Figs. 1-6) showed moderate inflation, in particular for the SNV-loop diuretic interaction terms for TG, LDL, and HDL analyses. The P values based on a t-distribution (red crosses) that we used for our meta-analyses have better genomic control values than the P values based on a standard normal distribution (black circles). Manhattan plots show the adjusted −log 10 (P) values for loop diuretics ( Fig. 1 and Supplementary  Figs. 7, 8) and thiazide diuretics (Supplementary Figs. 9-11) for TG, HDL, and LDL, respectively. Each figure includes plots of EA and AA separately, and for transancestry meta-analysis using METAL and MANTRA. The SNVs reaching genome-wide significance (P < 5 × 10 −8 ) for association with gene-medication interactions on TG concentration are shown in Table 2. Supplementary Table 4 shows the SNVs with suggestive evidence of SNV-diuretic interactions on each lipid trait (with P < 10 −6 from METAL or Bayes Factor > 10 5 from MANTRA).
Trans-ancestry association tests strengthened evidence of association for the KIAA1217 locus, elevating two of the six SNVs in this locus from suggestive in EA to genomewide significance. In particular, at rs10508671, the SNVloop diuretic interaction effect (β GD ) was consistent in both ancestries (−0.0011 vs −0.0010 mmol/L, Supplementary Table 4), and therefore trans-ancestry test provided stronger evidence of association (EA P = 9.43 × 10 −7 , AA P = 0.002, trans-ancestry P = 1.08 × 10 −8 , MANTRA log 10 [BF] = 6.23, Fig. 2c). MANTRA also strengthened the association for a locus on chromosome 8 (most significant rs6985929, EA P = 8.86 × 10 −5 , trans-ancestry P = 1.02 × 10 −5 , MANTRA log 10 [BF] = 6.17) that overlaps an intron in BAALC and that is~32 kb upstream of FZD6. Notably, the effect sizes of SNV-loop diuretic interactions are small (ranging from 0.0010 to 0.0025 mmol/L, Table 2), and were greater for AA than for EA or trans-ancestry analyses.

Thiazide diuretics with blood TG, HDL, and LDL
For analyses of SNV-thiazide diuretic interactions in EA and AA, no loci reached genome-wide significance for TG, HDL, or LDL ( Supplementary Figs. 9-11). In EA, suggestive association was noted with two loci including MYO16 (rs1926511, P = 1.01 × 10 −7 ) for TG and one locus for LDL. In AA, three loci yielded suggestive associations with HDL including a seven-SNV locus on chromosome 14 in TRIP11 (most significant rs10083510, P = 3.20 × 10 −7 ). SNVs in this locus are eQTLs for FBLN5.

Discussion
Genome-wide meta-analyses with 39,710 EA individuals and 9925 AA individuals were performed to identify genetic variants that modify the effect of loop or thiazide diuretic use on blood lipid concentrations. As the effects of loop and thiazide diuretic classes on lipids may differ [5,6], each with different effects on renal sodium excretion [18], analyses were performed on each of diuretic classes separately. Ancestry-specific analyses (of EA and AA data) identified genome-wide significant variants in three loci with SNV-loop diuretic interaction on TG concentrations; trans-ancestry analysis strengthened evidence of association at two additional loci.
The most significant associations with TG concentrations in EA were observed at intronic SNVs clustered in COL11A1 (collagen, type XI, alpha 1) on chromosome 1. COL11A1 produces a structural/adhesion protein secreted by primary rat adipocytes [33]. The expression of this gene is upregulated during human adipogenesis [34], suggesting a role in the modulation of lipid traits. Although it remains unclear how diuretics modulate the effect of COL11A1 variants on lipid traits, a COL11A1 mutation has been shown to be associated with nephrogenic diabetes insipidus [35], which is characterized by an inability of the kidney to concentrate urine, thus providing a possible mechanistic link between COL11A1, diuretics, and lipids.
Trans-ancestry analysis combining European and African results provided evidence of association at intronic SNVs clustered in KIAA1217 (an uncharacterized long coding DNA) which are just 5′ to an intronic microRNA (MIR603) on chromosome 10. KIAA1217 variants have been associated with obesity [36]. MIR603 has been associated with benign and malignant tumors [37][38][39][40]. Of potential interest, SNVs in both KIAA1217 and COL11A1 have been associated with lumbar disk herniation which, in turn, has been associated with blood lipid concentrations [41] and cardiovascular risk traits and disease [42].
Loci with genome-wide significant associations found in intragenic regions may impact the transcription of genes tens or hundreds of kilobases away from their target genes [43]. The cluster of SNVs on chromosome 20 reaching genome-wide significance for TG in AA lies between two genes. FERMT1 (fermitin family member 1) is a membrane-associated protein that links intracellular structural proteins to the extracellular matrix. While mutations in this gene have been implicated in a heritable skin disorder [44], connections to cardiovascular traits are sparse. Perhaps more relevant is this locus' proximity to BMP2 (bone morphogenetic protein 2), the expression of which is regulated by angiotensin II [45], a key regulator of blood pressure and fluid/electrolyte balance. Furthermore, BMP2, a member of the transforming growth factor-β superfamily of cytokines, is expressed in endothelial cells, and has been shown to stimulate adjacent smooth muscle and endothelial cells to proliferate, differentiate, and deposit extracellular matrix, thus affecting blood pressure and contributing to atherosclerosis by processes that also includes Wnt signaling [46]. BMP2 also participates in the differentiation of mesenchymal stem cells to mature adipocytes [47], thus potentially affecting blood lipid concentrations. Variants in the BMP2 receptor gene (BMPR2) have also been implicated in familial forms of pulmonary arterial HTN [48].
A cluster of SNVs spanning the 3′ end and intragenic region adjacent to BAALC and 5′ of FZD6 was nominally associated with TG response to loop diuretic therapy in both races; however, trans-ancestry meta-analysis by METAL and MANTRA strengthened the association. There is little evidence to link BAALC (brain and acute leukemia, cytoplasmic) with TG concentrations. However, both FZD6 (frizzled class receptor 6), which is~32 kb downstream from this locus on chromosome 8, and DKK1 (dickkopf WNT signaling pathway inhibitor 1), included in a suggestive locus on chromosome 10, are negative regulars of the β-catenin/Wnt signaling cascade [49,50] that not only alters the expression of genes which regulate renal sodium, chloride and potassium handling [51], but is inhibited by the loop diuretic ethacrynic acid [52]. Both FZD6 and DKK1 have also been associated with renal development [53,54], glomerular damage, and proteinuria [55,56].
Trans-ancestry analyses identified SNVs in TEAD1 (TEA domain family member 1 [SV40 transcriptional enhancer factor]) as significantly associated with HDL concentrations. This transcription factor is another β-catenin/Wnt pathway member [57] that has been shown to modulate adipocyte differentiation and proliferation [58] and play a role in regulating insulin sensitivity via skeletal muscle fiber type switching [59,60].
In the analyses of SNV-loop diuretic interactions in AA data with LDL, an SNV reaching genome-wide significance adjacent to PHIP (pleckstrin homology domain interacting protein) is within 270 kb of another SNV (rs16890334) previously identified as having a genome-wide significant association with blood pressure response to a high-sodium diet intervention in Han Chinese [61].
No loci from analyses of thiazide and thiazide-like diuretics with lipid traits reached genome-wide significance in EA or AA. However, trans-ancestry analyses identified a MACROD2 (MACRO domain containing 2) with genomewide significance for a diuretic interaction effect on HDL. This gene was found to have a significant association with coronary artery disease and HTN using a two-marker testing approach in a reanalysis of data from the Wellcome Trust Case Control Consortium [62], and a suggestive association with brain infarcts in another meta-analysis [63]; however, the role of this gene in determining cardiovascular traits remain uncertain.

Limitations
Although this study identified interesting SNV-loop diuretic interactions as modulators of TG concentrations, there were limitations in this study. First, the study design assumed therapeutic class effects although individual drugs may have different effects on lipids [13,64,65]. Information regarding the dosage, duration, or adequacy of diuretic therapy may have also affected the traits but were not available in most studies. Second, because potassium-sparing diuretics were frequently coadministered with thiazide and thiazide-like diuretics, it was not feasible to pursue the sole contributions of potassium-sparing diuretics in this study. Third, the effects of other classes of antihypertensive medications that may affect lipid levels, such as calcium, β-adrenergic, and α-adrenergic receptor blockers [66], were not assessed in this study. Fourth, the potential contribution of comorbidities (e.g., diabetes, renal dysfunction) and lifestyle choices (e.g., alcohol consumption, physical activity) on lipid traits was not considered. Fifth, the relatively smaller number of AA subjects exposed to loop diuretics may have reduced the power to detect significant association or increased the risk for spurious findings in this ancestry group. Finally, given the number of tests performed and the relatively small effect sizes for the interactions, this effort should be construed as a discovery analysis that warrants replication.

Conclusions
This genome-wide SNV-diuretic interaction meta-analysis on blood lipid concentrations used data from 39,710 EA individuals and 9925 AA individuals. The results of the present study suggest stronger interaction effects for loop versus thiazide diuretics, identifying several genome-wide significant loci of small effect sizes. The results of this study are typical of most GWAS studies wherein loci of small-tomodest effects are identified [67][68][69][70][71][72][73][74][75][76][77]. Although these loci are likely not action able in terms of genotype-guided therapy, some may nonetheless be biologically and clinically relevant [78], perhaps by identifying biologic pathways that contribute to adverse metabolic side effects from diuretic therapy. These findings may also explain, at least in part, some of the residual CVD risk following reduction of blood pressure in patients treated with antihypertensive medications.