Discovery and fine-mapping of loci associated with MUFAs through trans-ethnic meta-analysis in Chinese and European populations

MUFAs are unsaturated FAs with one double bond and are derived from endogenous synthesis and dietary intake. Accumulating evidence has suggested that plasma and erythrocyte MUFA levels are associated with cardiometabolic disorders, including CVD, T2D, and metabolic syndrome (MS). Previous genome-wide association studies (GWASs) have identified seven loci for plasma and erythrocyte palmitoleic and oleic acid levels in populations of European origin. To identify additional MUFA-associated loci and the potential functional variant at each locus, we performed ethnic-specific GWAS meta-analyses and trans-ethnic meta-analyses in more than 15,000 participants of Chinese and European ancestry. We identified novel genome-wide significant associations for vaccenic acid at FADS1/2 and PKD2L1 [log10(Bayes factor) ≥ 8.07] and for gondoic acid at FADS1/2 and GCKR [log10(Bayes factor) ≥ 6.22], and also observed improved fine-mapping resolutions at FADS1/2 and GCKR loci. The greatest improvement was observed at GCKR, where the number of variants in the 99% credible set was reduced from 16 (covering 94.8 kb) to 5 (covering 19.6 kb, including a missense variant rs1260326) after trans-ethnic meta-analysis. We also confirmed the previously reported associations of PKD2L1, FADS1/2, GCKR, and HIF1AN with palmitoleic acid and of FADS1/2 and LPCAT3 with oleic acid in the Chinese-specific GWAS and the trans-ethnic meta-analyses. Pathway-based analyses suggested that the identified loci were in unsaturated FA metabolism and signaling pathways. Our findings provide novel insight into the genetic basis relevant to MUFA metabolism and biology.

Abstract MUFAs are unsaturated FAs with one double bond and are derived from endogenous synthesis and dietary intake. Accumulating evidence has suggested that plasma and erythrocyte MUFA levels are associated with cardiometabolic disorders, including CVD, T2D, and metabolic syndrome (MS). Previous genome-wide association studies (GWASs) have identified seven loci for plasma and erythrocyte palmitoleic and oleic acid levels in populations of European origin. To identify additional MUFA-associated loci and the potential functional variant at each locus, we performed ethnic-specific GWAS meta-analyses and trans-ethnic metaanalyses in more than 15,000 participants of Chinese and European ancestry. We identified novel genome-wide significant associations for vaccenic acid at FADS1/2 and PKD2L1 [log 10  Supplementary key words genetics • fatty acid/desaturases • fatty acid/metabolism • fatty acid/biosynthesis • monounsaturated fatty acid MUFAs are nonessential unsaturated FAs with one double bond in the carbon chain, mainly including palmitoleic acid (16:1n-7), vaccenic acid (18:1n-7), oleic acid (18:1n-9), gondoic acid (20:1n-9), erucic acid (22:1n-9), and nervonic acid (24:1n-9). Oleic acid, the most abundant MUFA in lipids (1), is the predominant dietary MUFA (rich in plant oils, such as olive, canola, hazelnut, almond, and rapeseed, and in animal-derived fats, such as lard, tallow, and butter) (1,2). Oleic and palmitoleic acid can be synthesized through de novo lipogenesis (DNL) by -9 desaturation of stearic and palmitic acid, respectively, in liver and adipose tissue (3,4). Vaccenic, gondoic, erucic, and nervonic acids are elongation products of palmitoleic or oleic acid through endogenous synthesis (1).
MUFAs are important components of cell membranes and serve as energy sources through -oxidation in the mitochondria (for example, in the skeletal muscle during exercise) (5). Previous animal studies and clinical trials have suggested that palmitoleic and oleic acid play a role in lipid and glucose regulation (6)(7)(8)(9). In addition, epidemiological studies indicated that elevated levels of specific plasma and erythrocyte membrane MUFA levels (palmitoleic, vaccenic, erucic, and nervonic acids) were associated with increased risk of T2D (10), metabolic abnormalities (11), and CVD (12)(13)(14) in European populations, and studies in Chinese Hans also observed that higher levels of erythrocyte palmitoleic and oleic acid levels were associated with increased risk of metabolic syndrome (MS) and T2D (15,16). Therefore, these MUFAs are of great importance to cardiometabolic diseases.
Recently, genome-wide association studies (GWASs) have identified seven loci (FADS1/2, PKD2L1, HIF1AN, GCKR, 2p13, LPCAT3, and TRIM58) associated with plasma and erythrocyte palmitoleic acid and/or oleic acid levels in populations of European ancestry (17,18). However, these loci have not been replicated in other ethnic groups with different genetic architecture and dietary intake (19). More importantly, the potential functional variants at these loci and the genetic variants that affect circulating levels of vaccenic, gondoic, erucic, and nervonic acid remain unknown. Therefore, we first performed ethnic-specific GWAS meta-analyses in populations consisting of 3,521 Chinese individuals from two cohorts and of 12,020 European individuals from eight cohorts, respectively, and then combined the GWAS data from these two ethnic groups by trans-ethnic meta-analysis. We further conducted transethnic fine-mapping of the identified loci by construction of 99% credible sets.

Study cohorts
Chinese-specific GWAS meta-analysis included data on 3,521 Chinese-ancestry participants from two cohorts: the Nutrition and Health of Aging Population in China (NHAPC) (20) and the Multi-Ethnic Study of Atherosclerosis (MESA) (21). Europeanspecific GWAS meta-analysis included data on 12,020 Europeanancestry participants from eight cohorts: the Atherosclerosis Risk in Communities study (ARIC) (22), the Coronary Artery Risk Development in Young Adults study (CARDIA) (23), the Cardiovascular Health Study (CHS) (24), the Genetics of Lipid Lowering Drugs and Diet Network (GOLDN) (25), the Nurses' Health Study (NHS) (26), the Health Professionals Follow-Up Study (HPFS) (27), the Invecchiare in Chianti Study (InCHIANTI) (28), and MESA. Summary statistics of the two ethnic-specific GWAS meta-analyses were combined in the trans-ethnic metaanalysis. Written informed consent was obtained from all participants, and each study was approved by local ethics committees. Detailed descriptions of each cohort are provided in the supplemental Methods.

FA measurements
For ARIC, CARDIA, CHS, and MESA, fasting plasma phospholipids were isolated by TLC, and FAs were subsequently quantified by gas chromatography. For InCHIANTI, fasting FAs were measured in total plasma by gas chromatography. For GOLDN, NHS, HPFS, and NHAPC, fasting erythrocyte FAs were measured by gas chromatography or gas-liquid chromatography. FAs were identified for each study, and concentration of each FA was expressed as percentage of total FAs. Detailed methods of MUFA measurement in each cohort are provided in supplemental Methods. Ethnic-specific GWAS meta-analyses were performed for associations with palmitoleic, vaccenic, oleic, gondoic, erucic, and nervonic acids in Chinese populations, and for associations with vaccenic, gondoic, erucic, and nervonic acids in European populations ( Table 1).

Chinese-and European-specific GWAS meta-analyses
For each individual FA, GWAS analysis of approximately 2.2 million genotyped and imputed SNPs was performed separately in each cohort. Linear regression models were applied to test associations between SNPs and individual FA levels under an additive genetic model. All analyses were adjusted for age, sex, site of recruitment (as needed), and principal components to account These results were previously published in European populations from the CHARGE consortium (18).
for population admixture. Genomic control corrections were applied to each study before meta-analysis to further minimize potential confounding by population stratification (34). Cohort-specific association results were then combined by using inverse-variance based meta-analysis in METAL software (http://www.sph.umich. edu/csg/abecasis/metal). In Chinese-and European-specific meta-analyses, SNPs that were present in only one study were excluded. The linkage disequilibrium (LD) measures (r 2 ) were calculated by using the HapMap Phase II data (http://archive. broadinstitute.org/mpg/snap/ldsearchpw.php).

Trans-ethnic meta-analysis
Trans-ethnic meta-analysis was performed to identify additional novel loci for MUFA levels, to narrow the functional regions represented by the identified association signals, and also to test the heterogeneity of the identified loci across ethnic groups. Association statistics in the Chinese and European populations were combined by using Meta-ANalysis of Trans-ethnic Association studies (MANTRA; using a Bayesian framework) (35). MAN-TRA allows for heterogeneity among distinct ethnic populations and has increased power and mapping resolution compared with fixed-and random-effects meta-analysis (35). A log 10 [Bayes factor (BF)] 6 is considered as significant evidence of an association, and SNPs with posterior probability of heterogeneity (Phet) >0.5 were interpreted as having significant heterogeneity. Trans-ethnic meta-analysis using the fixed-effect model by METAL software was also performed. The 95% and 99% credible sets surrounding the most significant SNP at each locus based on the European-ancestry GWAS meta-analysis and the trans-ethnic meta-analysis of European and Chinese populations were calculated, respectively. The 95% and 99% credible sets at each locus was established by 1) defining a ±500 kb region surrounding the most significant SNP; 2) ranking the regional SNPs within this region according to their BF values; and 3) adding SNPs until the cumulative posterior probabilities of the ranked SNPs achieved the 95% and 99% confidence. Variants in the 99% credible sets were examined for sequence overlap with potential regulatory sites (promoter histone marks, enhancer histone marks, and DNase hypersensitivity) by searching the HaploReg V4.1 database (36). The location of each SNP in the fine-mapping section was extracted from the 1000 Genomes Brower (https://www.ncbi.nlm.nih.gov/variation/tools/ 1000genomes/).

cis-Expression quantitative trait loci analysis
In the cis-expression quantitative trait loci (cis-eQTL) analysis, we examined the associations of significant SNPs ( Table 2) with RNA levels of nearby genes in the public Genotype-Tissue Expression database (GTEx; http://gtexportal.org/home/). Adipose tissue, skeletal muscle, and liver were considered as the most important tissues in FA metabolism (4), and these three tissues were selected in the GTEx browser.

Cohort characteristics
The study sample included 3,521 individuals of Chinese origin from the NHAPC and MESA cohorts and 12,020 individuals of European origin from cohorts including ARIC, CARDIA, CHS, GOLDN, NHS, HPFS, InCHIANTI, and MESA (Table 1). Participants comprised mostly middleaged to older individuals (mean age across the cohorts ranged from 45.8 to 75.0 years), and approximately 50% of them were male, except for the NHS (female only) and HPFS (male only) cohorts. FA levels were all expressed as the percentage of total FAs (Table 1).

Fine mapping
Although the SNPs showing the strongest association are not necessarily the functional variants, owing to factors such as sampling variation and LD patterns, it is still reasonable to think that functional variants are among the SNPs interrogated in these loci. By comparing the 99% credible sets calculated on the basis of the Europeanancestry meta-analysis and on the basis of the trans-ethnic meta-analysis, we observed improvements in fine-mapping resolution at FADS1/2 and GCKR loci.
The number of SNPs in the 99% credible set at FADS1/2 locus was reduced from 16 (covering 54.1 kb) to 14 (covering 53.5 kb, from rs174535 located 32.3 kb upstream of FADS2 to rs174577 in the first intron of FADS2), 24 (covering 79.6 kb) to 23 (covering 62.0 kb, from rs174528 located 40.2 kb upstream of FADS2 to rs174578 in the second intron of FADS2), and 15 (covering 54.1 kb) to 12 (covering 53.5 kb, from rs174535 located 32.3 kb upstream of FADS2 to rs174577 in the first intron of FADS2), for association with palmitoleic, vaccenic, and oleic acid, respectively (supplemental Tables S3-S5), and these SNPs show strong LD with each other (r 2  0.590 in CEU and r 2  0.927 in CHB+JPT). After trans-ethnic meta-analysis, the 12 highly correlated SNPs, including rs174535, rs174545, rs174546, rs102275, rs174536, rs174537, rs174550, rs174547, rs174574, rs174576, rs174577, and rs1535 (r 2 = 1 in CEU and r 2  0.927 in CHB+JPT, from rs174535 located 32.3 kb upstream of FADS2 to rs174577 in the first intron of FADS2) in the 99% credible set of oleic acid were shared across the associations with palmitoleic, vaccenic, and oleic acid (supplemental Table S3). We further performed a Fisher's exact test to evaluate whether the SNPs in the trans-ethnic analysis-based 99% credible sets were enriched for potential regulatory sites (promoter histone marks, enhancer histone marks, and DNase hypersensitivity) when compared with the SNPs in the European-ancestry GWAS-based 99% credible sets. However, no significant enrichment was observed after trans-ethnic meta-analysis (P  0.533). These findings suggested that the improved fine-mapping resolution at FADS1/2 locus is likely resulted from the increased sample size after trans-ethnic meta-analysis.
The greatest improvement in fine-mapping resolution was observed at the GCKR locus, where the number of SNPs in the 99% credible set was reduced from 16 (covering 94.8 kb) on the basis of the European-ancestry meta-analysis to five (covering 19.6 kb, from rs1260326 located in the 15th exon of GCKR to rs2911711 located 4.0 kb downstream of GCKR) after trans-ethnic meta-analysis. The five SNPs in the 99% credible, in strong LD with each other (r 2  0.805 in CEU and r 2  0.827 in CHB+JPT), included one missense SNP (rs1260326, p.P446L, located in the 15th exon of GCKR), two intronic SNPs (rs780094 and rs780093 located in the 16th and 17th intron of GCKR, respectively), and two SNPs in the downstream of GCKR (rs1260333 and rs2911711 located 2.1 and 4.0 kb downstream of GCKR, respectively) (supplemental Tables S4 and S5).
The 99% credible sets based on the European-ancestry and trans-ethnic meta-analyses at PKD2L1 locus contained only one SNP (rs603424). Therefore, the posterior probability that rs603424 was functional (or tagged an unobserved functional variant) was greater than 99% for PKD2L1. PKD2L1-rs603424 is located in the second intron of the PKD2L1 gene and overlapped with an enhancer histone mark (H3K4me1, H3K27ac and H3K9ac) in the adipose tissue when we searched the Roadmap Epigenomics database (40) using HaploReg V4.1 (supplemental Tables S4  and S5). The previously reported signal HIF1AN-rs10883511 for association with palmitoleic acid (18) is located 224 kb away from PKD2L1-rs603424, and these two SNPs are independent from each other (r 2 = 0.013 in CEU and r 2 = 0.001 in CHB+JPT). Because PKD2L1-rs603424 was the strongest association signal ( = 0.032, P = 5.32 × 10 14 ) in the ±500 kb genomic region of HIF1AN-rs10883511 ( = 0.023, P = 6.44 × 10 8 ), we could not construct credible sets or perform trans-ethnic fine-mapping at HIF1AN locus.

cis-eQTL analysis
To gain more insight into the potential functional roles of the genome-wide significant loci (Table 2), we performed cis-eQTL analysis by searching the publicly available GTEx database in adipose tissue, skeletal muscle, and liver. Results suggested that the minor allele A of PKD2L1-rs603424 was significantly associated with decreased RNA level of SCD (stearoyl-CoA desaturase; encodes the -9 desaturase in the DNL pathway) in the adipose tissue (P  3.94 × 10 6 ; supplemental Table S6). These results further strengthened our findings from the GWAS meta-analyses that the minor allele A of PKD2L1-rs603424 was significantly associated with decreased levels of palmitoleic and vaccenic acids (P  5.32 × 10 14 ; Table 2).

Gene-and pathway-based analysis
Gene-based analysis improves the statistical power by combining all SNPs in a gene into a gene-based score, which reduces the burden of multiple testing and incorporates multiple independent association signals (41). To identify additional genes and pathways that contribute to circulating MUFA levels and gain insight into the underlying mechanisms, we performed gene-and pathway-based association testing using GWAS summary statistics in Chinese and European populations, respectively. Four genes (FEN1, FADS1, FADS2, and LPCAT3; 27 SNPs) for association with oleic acid in the Chinese populations and six genes (SCD, WNT8B, NDUFB8, FEN1, FADS1, and FADS2; 64 SNPs) for association with palmitoleic, vaccenic, oleic, and/or gondoic acid in the European populations reached gene-based significance (P  1.60 × 10 6 ; supplemental Table S7) using both SPU and GATES methods. Because the association signals at WNT8B and NDUFB8 genes for association with palmitoleic acid in the European populations are in strong LD with the reported SNP HIF1AN-rs10883511 (18) (r 2  0.832), we did not consider them as novel loci. Pathways, including biosynthesis of unsaturated FAs, -linolenic acid metabolism, glycerophospholipid metabolism, and PPAR signaling pathway, were significantly associated with palmitoleic, vaccenic, oleic, and/or gondoic acid levels in the pathway-based analyses (P  2.28 × 10 5 ; supplemental Table S8).

DISCUSSION
In this first trans-ethnic meta-analysis of MUFA levels in Chinese and European populations, we identified novel associations of FADS1/2, PKD2L1, and GCKR with vaccenic and/or gondoic acid and replicated the previously reported associations with palmitoleic and/or oleic acid for loci at FADS1/2, PKD2L1, GCKR, HIF1AN, and LPCAT3 in the Chinese-specific GWAS and the trans-ethnic meta-analyses. We also observed substantial improvement in the finemapping resolution at the GCKR locus after trans-ethnic meta-analysis.
The -9 desaturase, encoded by SCD, plays an important role in MUFA metabolism (42). It catalyzes the desaturation reaction from palmitic and stearic acid to palmitoleic and oleic acid, respectively, in the DNL pathway. The present study and the previous GWAS (18) have identified significant associations of PKD2L1 with palmitoleic and vaccenic acid, which is located near SCD. The results of the cis-eQTL analysis suggested that the most significant SNP, PKD2L1-rs603424, was associated with the RNA level of SCD in the adipose tissue (P  3.94 × 10 6 ; supplemental Table S5). Therefore, PKD2L1-rs603424 may exert its effect on MUFA levels through regulating SCD transcription.
Genetic variants at FADS1/2 also showed significant associations with MUFA levels, including palmitoleic, vaccenic, oleic, and gondoic acid. FADS1/2 encode 5 and 6 desaturases, predominantly involved in the PUFA biosynthesis pathway (43). Recent studies have indicated that 6 desaturase can also catalyze palmitic and stearic acid to produce other unsaturated FAs (1,44). Because palmitic and stearic acid serve as the substrates for MUFA endogenous synthesis, it is possible that FADS1/2 influences specific MUFA levels through substrate regulation. The exact mechanisms underlying the associations between FADS1/2 variants and MUFA levels merit further investigation.
GCKR encodes glucokinase regulator, a protein that inhibits glucokinase (GCK) activity in liver and pancreas (45). In the present study, substantial improvement of finemapping resolution at this locus was observed, and the 99% credible set calculated after trans-ethnic meta-analysis highlighted one missense variant (rs1260326, p.P446L). Functional studies have indicated that GCKR-rs1260326 played a central role in regulation of GCK activity in the liver, which consequently influenced glycolytic flux and DNL (46,47). Therefore, it is likely that rs1260326 is the variant driving the association of GCKR with palmitoleic acid through modifying the DNL pathway.
LPCAT3 encodes lysophosphatidylcholine acyltransferase 3, which is involved in lysophospholipid esterification (48). It was confirmed in the trans-ethnic meta-analysis using MANTRA, but not METAL. The inconsistent trans-ethnic meta-analysis results generated from MANTRA and METAL were mainly due to the different methods implemented in these two software. MANTRA takes account of the expected similarity in allelic effects between the most closely related populations by using a Bayesian partition model and also allows for heterogeneity across ethnic groups (49). In contrast, METAL uses a fixed-effect model that assumes the allelic effect to be the same in all populations. Therefore, MANTRA confers significantly higher power than METAL (50).
In the ethnic-specific GWAS meta-analyses, the magnitude and directions of the identified associations were largely consistent across cohorts (supplemental Table S1 and supplemental Fig. S3). Association analyses further stratified by the measurement methods (erythrocyte and plasma phospholipid FAs) did not materially changed the association results, and no additional erythrocyte-or plasma-specific loci were identified (supplemental Table  S1). Association results after excluding the InCHIANTI cohort, which measured MUFA levels in total plasma, also remained largely unchanged. These results suggested that differences in MUFA measurements did not introduce noise in the association results and were similar to the findings of previously published GWAS of FAs (18,51).
In conclusion, this is the first study to provide evidence for the associations of FADS1/2, PKD2L1, and GCKR with vaccenic and/or gondoic acid levels in the populations of Chinese and European origin. Five previously reported loci (FADS1/2, PKD2L1, GCKR, HIF1AN, and LPCAT3) for palmitoleic and/or oleic acid were also confirmed. Trans-ethnic fine-mapping highlighted a missense variant (rs1260326) at the GCKR locus. Our findings shed light on the genetic basis of MUFA biology and establish the foundation for future genetic and functional investigations.