Vitamin D and Nonalcoholic Fatty Liver Disease: Bi-directional Mendelian Randomization Analysis.

BACKGROUND: Vitamin D deficiency is associated with nonalcoholic fatty liver disease (NAFLD) in many cross-sectional studies. However, the causality between them has not been established. We used bi-directional mendelian randomization (MR) analysis to explore the causal relationship between 25-hydroxyvitamin D [25(OH)D] and NAFLD.
METHODS: 9182 participants were included from a survey in East China from 2014 to 2016. We calculated weighted genetic risk scores (GRS) for 25(OH)D concentration and NAFLD based on 25(OH)D-related and NAFLD-related single nucleotide polymorphisms. Presence of liver steatosis was assessed using ultrasound. Instrumental variable was used to measure the causal relationship between them.
RESULTS: An SD increase in the 25(OH)D GRS was significantly associated with 25(OH)D (β 1.29, 95%CI -1.54, -1.04, P<0.05) but not with NAFLD (OR 0.97, 95%CI 0.92, 1.01). An SD increase in NAFLD GRS was also strongly associated with NAFLD (OR 1.09, 95%CI 1.04, 1.15, P<0.05) but not with 25(OH)D (β -0.15, 95%CI -0.41, 0.10). Using an instrumental variable estimator, no associations were found for genetically instrumented 25(OH)D with NAFLD and for genetically instrumented NAFLD with 25(OH)D.
CONCLUSION: Our results support the conclusion that there is no causal association between vitamin D and NAFLD using a bi-directional MR approach in a Chinese population.

Introduction

Non-alcoholic fatty liver disease (NAFLD) is a major cause of liver disease worldwide, and carries an increased all-cause mortality (Younossi et al., 2016). In a study including population from 22 countries, the prevalence reached 25.24% globally (Younossi et al., 2016). It has been increasingly recognized as a component of metabolic syndrome (Chalasani et al., 2012), thus the clinical and socioeconomic burden of NAFLD is becoming huge. Despite a known genetic contribution (Macaluso et al., 2015), a sedentary lifestyle and overnutrition have set the stage for the epidemic of NAFLD. However, public health measures and other interventions could be effective in its prevention (Fan et al., 2017).

Vitamin D deficiency is also a pandemic health problem in both developing and developed countries (Palacios and Gonzalez, 2014). Classically, vitamin D exerts a function on calcium/phosphorus homeostasis; however increasing evidence shows it is associated with NAFLD (Kwok et al., 2013). Previous studies found that lower 25-hydroxyvitamin D [25(OH)D] was associated with NAFLD and its severity (Wang et al., 2016c, Zhai et al., 2016), which was consistent with other cross-sectional studies (Wang et al., 2016a, Eliades et al., 2013). Animal studies have also shown active vitamin D could attenuate hepatic steatosis by preventing autophagy and oxidative stress (Li et al., 2017, Zhu et al., 2017). However, the very limited intervention studies testing the effect of vitamin D supplementation on patients with NAFLD have inconsistent findings (Della Corte et al., 2016, Sharifi et al., 2014, Barchetta et al., 2016). Thus, the causality between vitamin D and NAFLD has not been confirmed in human beings.

Mendelian randomization (MR) uses genetic variants in non-experimental data to make causal inferences regarding the effect of an exposure on an outcome (Smith and Ebrahim, 2003). In this study, if low 25(OH)D causally leads to NAFLD, genetic variants associated with lower 25(OH)D should be associated with higher NAFLD risk; conversely, if NAFLD induces low 25(OH)D, then genetic variants associated with higher NAFLD risk should be related to lower 25(OH)D concentrations. These genetic variants are inherited independent of potential confounding factors (Lawlor et al., 2008). Thus, MR could avoid problems in conventional epidemiological studies such as residual confounding and reverse causation (Smith and Ebrahim, 2003).

Based on a large community-based sample of Chinese participants from the SPECT-China study (Survey on Prevalence in East China for metabolic diseases and risk factors), we performed bidirectional MR analyses to explore the causal association between 25(OH)D and NAFLD as detected using ultrasound. Vitamin D and NAFLD genetic risk scores (VD_GRS and NAFLD_GRS) were constructed to represent the genetic susceptibility. We analyzed the causal link between genetically determined 25(OH)D status or NAFLD and risk of NAFLD or low 25(OH)D, respectively.

Materials and Methods

Participants

The data were from an ongoing SPECT-China study, which is a large cross-sectional study. Recruitment and enrollment have been previously described in detail (Wang et al., 2015, Wang et al., 2017a, Wang et al., 2017b). Chinese citizens ≥ 18 years old who had lived in their current area for ≥ 6 months were selected. We excluded subjects with severe communication problems, acute illness or who were unwilling to participate. From 2014 to 2016, 12666 subjects from 18 to 93 years in age were recruited in the SPECT-China study from 23 sites in Shanghai, Zhejiang, Jiangsu, Anhui and Jiangxi provinces. Among them, genotype information was available for 10,664 participants (84.2%). We excluded participants who had missing information on more than two single nucleotide polymorphism (SNP) genotypes (n = 182), liver ultrasound (n = 399) and 25(OH)D (n = 1). We also excluded those who had a history of excessive consumption (male > 20 g/d, female > 10 g/d) of pure alcohol (n = 706), schistosomal hepatic disease (n = 6), self-reported viral hepatitis (including hepatitis B and hepatitis C virus) (n = 134), or using medications associated with secondary NAFLD (corticosteroids, estrogens, amiodarone, methotrexate) (n = 54). A total of 9182 participants were involved in the final analysis.

The study protocol was approved by the Ethics Committee of Shanghai Ninth People's Hospital, Shanghai JiaoTong University School of Medicine. The study protocol conformed to the ethical guidelines of the 1975 Declaration of Helsinki as reflected in a priori approval by the appropriate institutional review committee. Informed consent was obtained from all participants included in the study.

Measurements

Interviews and collection of biological specimens at each site were undertaken with a single assessment protocol. Blood samples were obtained between 7:00 am and 10:00 am after fasting for at least 8 h. Blood was refrigerated immediately after phlebotomy, and after 2–4 h it was centrifugalized and the serum was aliquoted and frozen in a central laboratory. Glycated hemoglobin (HbA1c) was measured using high-performance liquid chromatography (MQ-2000PT, Medconn, Shanghai, China). Fasting plasma glucose, triglycerides, total cholesterol, high density lipoprotein (HDL) and low density lipoprotein (LDL) were measured using the AU 680 (Beckman Coulter, Brea, USA). The 25(OH)D (ADVIA Centaur XP, Siemens, Germany), were detected using a chemiluminescence assay.

Body mass index (BMI) was calculated as weight in kilograms divided by height in meters squared. Current smoking was defined as having smoked at least 100 cigarettes in one's lifetime and currently smoking cigarettes (Xu et al., 2013).

Definition

As previously described, two experienced ultrasonographers used an ultrasound device (Mindray M7, MINDRAY, Shenzhen, China) to perform an abdominal ultrasonographic examination. The diagnostic criteria for fat accumulation (steatosis) included increased liver echogenicity, stronger echoes in the hepatic parenchyma as compared to in the renal parenchyma, vessel blurring, and narrowing of the lumen of the hepatic veins (Wang et al., 2016b, Saadeh et al., 2002).

Diabetes was determined using a previous diagnosis by health care professionals, fasting plasma glucose level ≥ 7.0 mmol/L or HbA1c ≥ 6.5%. Hypertension was assessed by systolic blood pressure ≥ 140 mm Hg, or diastolic blood pressure ≥ 90 mm Hg, or self-reported previous diagnosis of hypertension by physicians.

Genotyping, Genetic Loci Selection and Genetic Risk Score Construction

DNA was extracted from white blood cells using a blood genomic DNA extraction kit (DP603, TIANGEN BIOTECH CO, LTD, Beijing, China) on an automated nucleic acid extraction instrument (YOSE-S32, TIANGEN BIOTECH CO, LTD, Beijing, China). Specific assays were designed using Geneious Pro (v4.8.3) (https://www.geneious.com/). Mass determination was carried out with the JUNO (Fluidigm Corporation, South San Francisco, California, USA) and data acquisition was completed using Fluidigm SNP Genotyping Analysis v4.1.3software (Fluidigm Corporation, South San Francisco, California, USA). Call rates of all SNPs were > 98%.

We selected eight SNPs involved in susceptibility and/or progression of NAFLD, lysophospholipase-like 1 [LYPLAL1]- rs12137855, protein phosphatase 1 regulatory subunit 3b [PPP1R3B]- rs4240624, transmembrane 6 superfamily member 2 [TM6SF2]- rs58542926, patatine-like phospholipase domain containing 3 [PNPLA3]- rs738409, glucokinase regulatory protein [GCKR]- rs780094, sorting and assembly machinery component [SAMM50]- rs738491, parvin beta [PARVB]- rs5764455 and collagen type XIII alpha 1 chain [COL13A1]- rs1227756, for our analysis based on previously published genome-wide association studies for NAFLD-related traits (Wang et al., 2016d, Lin et al., 2014, Kitamoto et al., 2013, Macaluso et al., 2015). The four vitamin D-related SNPs (7-dehydrocholesterol reductase [DHCR7]- rs12785878, cytochrome P450 family 2 subfamily R member 1 [CYP2R1]- rs10741657, vitamin D binding protein [GC]- rs2282679, and cytochrome P450 family 24 subfamily A member 1 [CYP24A1]-rs6013897) were chosen on the basis of a recent genome-wide association study on 25(OH)D (Huang et al., 2016). They all reached a genome-wide significance level (P < 5 × 10− 8) and were not in linkage disequilibrium (r2 = 0, except for 0.61 between rs738491 and rs738409, 0.50 between rs738491 and rs5764455, and 0.51 between rs738409 and rs5764455) (Kitamoto et al., 2013). Because rs4240624 deviated from the Hardy–Weinberg equilibrium with a P < 10− 4, it was eliminated from the GRS.

Statistical Analysis

Data analyses were performed using IBM SPSS Statistics, Version 22 (IBM Corporation, Armonk, NY, USA). All analyses were two-sided. A P value < 0.05 indicated significance. Continuous variables were expressed as the mean ± standard deviation (SD) and categorical variables as a percentage (%), respectively. Serum TG was logarithmically transformed prior to analysis.

The additive genetic model for each SNP (coded as 0–2) was used to construct GRS. For the VD_GRS, we created a weighted score by multiplying each SNP by a weight based on its effect size with 25(OH)D obtained from a large study of an Asian population (Cuellar-Partida et al., 2017). For the NAFLD_GRS, the weights were also from previous Asian population studies (Kitamoto et al., 2013, Wang et al., 2016d, Shang et al., 2015). The characteristics of each SNP in the VD_GRS and NAFLD_GRS are summarized in Supplemental Table 1.

Linear regression analyses were used to determine the association of the two GRSs and present NAFLD with 25(OH)D. Logistic regression models were fitted to analyze the association of the two GRSs and 25(OH)D with NAFLD. Model 1 adjusted for age, sex and BMI. Model 2 adjusted for the terms in model 1, current smoking (yes or no), hypertension (yes or no), diabetes (yes or no), HDL-cholesterol, LDL-cholesterol and triglycerides.

Regarding MR analysis, the weighted VD_GRS and weighted NAFLD_GRS were used as the instrumental variable (IV) estimators to measure the strength of the bidirectional causal relationship between 25(OH)D and NAFLD. The formal MR analyses to estimate the possible causal effect of v25(OH)D on NAFLD (and vice versa) were completed using the IV ratio method (Lawlor et al., 2008). For the causal association of increased risk of NAFLD in relation to lower 25(OH)D, the computational formula was ORIV(VD-NAFLD) = exp (ln (ORVD_GRS-NAFLD)/βVD_GRS-25(OH)D). In the opposite direction, the computational formula was βIV(NAFLD-VD) = βNAFLD_GRS-25(OH)D/ln (OR NAFLD_GRS-NAFLD). The model was the same as the aforementioned model 2, adjusting for age, sex, BMI, current smoking (yes or no), hypertension (yes or no), diabetes (yes or no), HDL-cholesterol, LDL-cholesterol and triglycerides.

To validate the genetic instruments, we assessed the associations between each individual SNP with 25(OH)D and NAFLD, respectively. We also measured the potential pleotropic associations of each individual SNP and the GRSs with BMI, blood pressure, lipid profile and HbA1c.

In the sensitivity analyses, we constructed the unweighted VD_GRS, unweighted NAFLD_GRS and NAFLD_GRS4SNP excluding the three SNPs that may have a pleiotropic effect on NAFLD-related metabolic traits. Then, we used these as the IVs. Moreover, we also further adjusted for seasonal variation where data were categorized into winter (January to March), spring (April to June), summer (July to September) and autumn (October to December).

Results

Association of SNPs With 25(OH)D and NAFLD

We tested the association of seven NAFLD-related and four 25(OH)D-related SNPs with NAFLD and 25(OH)D. Unstandardized coefficients [95% confidence interval (CI)] and odds ratio (OR) (95% CI) of additive linear regression models are shown in Supplemental Figs. 1 and 2. In the seven NAFLD-related SNPs, three SNPs at the GCKR, PNPLA3 and PARVB loci were significantly associated with NAFLD in this study. In the four VD-related SNPs, two SNPs at the GC and DHCR7 loci were significantly associated with 25(OH)D.

Pleiotropic Effects of SNPs

We assessed whether the SNPs showed any association with NAFLD- and 25(OH)D- related major metabolic traits. Therefore, we measured the potential associations of the SNPs with BMI, LDL, HDL, triglycerides and HbA1c using an additive model. The results are summarized in Supplemental Table 2. None of the 25(OH)D- related SNPs had pleiotropic effects. Three of the NAFLD- related SNPs showed a significant association with at least one trait. rs780094 was associated with LDL, HDL and log (triglycerides); rs58542926 was associated with log (triglycerides); and rs738409 was associated with BMI.

Study Characteristics According to VD_GRS and NAFLD_GRS Quartiles

As expected, with increasing VD_GRS, 25(OH)D concentrations significantly decreased and with increasing NAFLD_GRS, the prevalence of NAFLD significantly increased. However, the two GRSs were not consistently associated with any of the investigated biochemical markers such as LDL, HDL, triglycerides, fasting plasma glucose and HbA1c, nor were they with other risk factors such as sex, smoking status, diabetes, hypertension or BMI (Table 1). There was also no trend between VD_GRS and NAFLD and between NAFLD_GRS and 25(OH)D.

Table 1

Characteristics of study participants according to the weighted vitamin D genetic risk score (VD_GRS) and NAFLD_GRS (n = 9182).

Characteristic	Quartile 1	Quartile 2	Quartile 3	Quartile 4	P for trend
VD_GRS
VD_GRS	≤ 0.62	0.63–0.88	0.89–1.16	≥ 1.17
N	2337	2279	2497	2069
Age, yr	54.3 (13.2)	54.2 (13.2)	54.6 (12.9)	54.8 (12.8)	0.18
Men, %	35.1	35.8	37.6	36.1	0.27
25(OH)D, mmol/L	41.8 (12.9)	40.4 (12.3)	39.6 (12.5)	38.7 (11.9)	< 0.001
Body mass index, kg/m2	24.6 (3.6)	24.6 (3.5)	24.6 (3.6)	24.6 (3.7)	0.97
Triglycerides, mmol/L	1.61 (1.36)	1.68 (1.59)	1.68 (1.43)	1.68 (1.32)	0.44
HDL, mmol/L	1.39 (0.32)	1.40 (0.32)	1.41 (0.32)	1.40 (0.32)	0.50
LDL, mmol/L	3.21 (0.83)	3.18 (0.81)	3.17 (0.80)	3.17 (0.79)	0.10
FPG, mmol/L	5.6 (1.4)	5.7 (1.5)	5.6 (1.5)	5.7 (1.5)	0.17
HbA1c, %	5.6 (0.9)	5.6 (1.0)	5.6 (1.0)	5.6 (1.0)	0.33
Current smoker, %	16.4	16.6	18.4	17.5	0.13
NAFLD, %	50.9	50.3	50.8	48.8	0.25
Diabetes, %	13.7	15.2	14.0	14.3	0.85
Hypertension, %	46.6	47.5	46.5	46.2	0.65
NAFLD_GRS
NAFLD_GRS	≤ 0.6596	0.6597–1.4291	1.4292–1.8996	≥ 1.8997
N	2278	2311	2310	2283
Age, yr	54.2 (13.3)	54.2 (13.1)	54.8 (12.7)	54.7 (13.1)	0.07
Men, %	35.5	34.7	37.3	37.3	0.07
25(OH)D, mmol/L	40.2 (12.4)	40.3 (12.7)	40.0 (12.2)	40.0 (12.4)	0.40
Body mass index, kg/m2	24.6 (3.5)	24.7 (3.7)	24.5 (3.4)	24.5 (3.7)	0.12
Triglycerides, mmol/L	1.70 (1.66)	1.60 (1.22)	1.70 (1.48)	1.61 (1.33)	0.67
HDL, mmol/L	1.40 (0.32)	1.41 (0.32)	1.40 (0.33)	1.39 (0.31)	0.40
LDL, mmol/L	3.17 (0.80)	3.19 (0.80)	3.19 (0.82)	3.17 (0.81)	0.94
FPG, mmol/L	5.6 (1.5)	5.6 (1.4)	5.7 (1.5)	5.6 (1.4)	0.53
HbA1c, %	5.6 (1.0)	5.6 (1.0)	5.6 (1.0)	5.6 (1.0)	0.95
Current smoker, %	17.2	16.6	17.4	17.8	0.48
NAFLD, %	48.8	49.4	51.5	51.3	< 0.05
Diabetes, %	14.0	13.6	14.9	14.6	0.32
Hypertension, %	45.8	46.7	47.8	46.5	0.51

The data are summarized as the mean (SD) for continuous variables or as a numerical proportion for categorical variables. P for trend was calculated by ANOVA and chi-square tests. 25(OH)D, 25-hydroxyvitamin D; FPG, fasting plasma glucose; GRS, genetic risk score; HDL, high-density lipoprotein; HbA1c, glycated hemoglobin; LDL, low-density lipoprotein; NAFLD, nonalcoholic fatty liver disease; TG, triglyceride.

Associations of VD_GRS and 25(OH)D With NAFLD

Next, we measured the association of the serum 25(OH)D concentration and VD_GRS with NAFLD. In Table 2, per SD increase in VD_GRS was not significantly associated with a decreased risk of NAFLD after adjustment for age, sex and BMI (OR 0.97, 95% CI 0.92, 1.01) (Model 1). Further adjusting for smoking, hypertension, diabetes and lipid profile did not change the results (OR 0.97, 95% CI 0.92, 1.01) (Model 2). The quartiles of VD_GRS showed no significant association with NAFLD.

Table 2

The associations of VD_GRS and 25(OH)D with NAFLD.

	Model 1	Model 2
Per SD increase of VD_GRS	0.97 (0.92, 1.01)	0.97 (0.92, 1.01)
Quartiles of VD_GRS
Q1	1.12 (0.98, 1.28)	1.14 (0.99, 1.32)
Q2	1.08 (0.94, 1.24)	1.07 (0.93, 1.24)
Q3	1.10 (0.96, 1.26)	1.15 (0.99, 1.32)
Q4	1.00 (Ref.)	1.00 (Ref.)
P for trend	0.16	0.18
Per SD increase of 25(OH)D	0.81 (0.77, 0.85)	0.96 (0.82, 0.91)

Data were presented as odds ratio and 95% confidence interval. 25(OH)D, 25-hydroxyvitamin D; GRS, genetic risk score; NAFLD, nonalcoholic fatty liver disease.

Model 1 adjusted for age, sex and BMI.

Model 2 adjusted for terms in model 1, current smoking (yes or no), hypertension (yes or no), diabetes (yes or no), HDL-cholesterol, LDL-cholesterol and triglycerides.

However, in this cross-sectional study, 1 SD increase of 25(OH)D was associated with a 19% (95% CI 0.77, 0.85) decreased prevalence of NAFLD in model 1. Further adjusting other metabolic profiles attenuated the association but it remained significant (OR 0.96, 95% CI 0.82, 0.91).

Associations of NAFLD_GRS and Present NAFLD With 25(OH)D

As shown in Table 3, each SD increase in NAFLD_GRS was not significantly associated with a decreased level of 25(OH)D in both models. The quartiles of NAFLD_GRS showed similar results. However, present NAFLD was associated with a 2.47 nmol/L decrease in 25(OH)D concentration in model 1. Further adjusting other metabolic profiles attenuated the association but it remained significant (OR 1.68, 95% CI 1.09, 2.27).

Table 3

The associations of NAFLD_GRS and NAFLD with 25(OH)D.

	Model 1	Model 2
Per SD increase of NAFLD_GRS	− 0.16 (− 0.41, 0.10)	− 0.15 (− 0.41, 0.10)
Quartiles of NAFLD _GRS
Q1	0.34 (− 0.38, 1.06)	0.41 (− 0.30, 1.12)
Q2	0.56 (− 0.16, 1.27)	0.50 (− 0.21, 1.21)
Q3	0.01 (− 0.70, 0.72)	0.10 (− 0.61, 0.81)
Q4	0.00 (Ref.)	0.00 (Ref.)
P for trend	0.17	0.16
Present NAFLD (yes)	− 2.47 (− 3.03, − 1.90)	− 1.68 (− 2.27, − 1.09)

Data were presented as unstandardized coefficients and 95% confidence interval. 25(OH)D, 25-hydroxyvitamin D; GRS, genetic risk score; NAFLD, nonalcoholic fatty liver disease.

Model 1 adjusted for age, sex and BMI.

Model 2 adjusted for terms in model 1, current smoking (yes or no), hypertension (yes or no), diabetes (yes or no), HDL-cholesterol, LDL-cholesterol and triglycerides.

25(OH)D and NAFLD: The MR Analysis

Fig. 1 shows the association of genetically determined 25(OH)D with the risk of NAFLD, and conversely, genetically determined NAFLD with 25(OH)D concentrations. In the IV analysis, the causal OR of genetically determined 25(OH)D for risk of NAFLD was 1.03 (95% CI: 0.99, 1.07), and the causal regression coefficient of genetically determined NAFLD for 25(OH)D was − 1.70 (95% CI − 4.63, 1.23). Both directions showed no significant association.

Fig. 1

Bidirectional instrumental variable (IV) estimated association between 25(OH)D and NAFLD by weighted GRSs. Data were presented as regression coefficient (β) or odds ratio (OR) and 95% confidence interval (CI). In this MR framework, the instrumental variable estimators are ORIV(VD-NAFLD) = exp (ln (ORVD_GRS-NAFLD)/βVD_GRS-25(OH)D) and βIV(NAFLD-VD) = βNAFLD_GRS-25(OH)D/ln (OR NAFLD_GRS-NAFLD). Data were adjusted for age, sex, BMI, current smoking, hypertension, diabetes, HDL-cholesterol, LDL-cholesterol and triglycerides. 25(OH)D, 25-hydroxyvitamin D; VD_GRS, vitamin D genetic risk score; NAFLD_GRS, nonalcoholic fatty liver disease genetic risk score; SD, standard deviation.

Sensitivity Analysis

There were three NAFLD-related SNPs having pleiotropic effects, though the weighted NAFLD_GRS was not associated with other metabolic traits (Supplemental Table 3). Thus, we excluded these three SNPs to construct NAFLD_GRS4SNP, which as expected, was not associated with these metabolic traits. The IV estimate for causal relationship from NAFLD to 25(OH)D was β − 1.54 (95% CI − 5.30, 1.22). Using the unweighted GRSs, the results were similar to those using weighted GRSs (See Fig. 2). Both direction showed no significant associations. In addition, when further adjusting for seasonal variation, the IV estimates for causal relationship from NAFLD to 25(OH)D and from 25(OH)D to NAFLD were β − 1.06 (95% CI -4.02, 1.90) and OR 1.03 (95% 0.99, 1.07), respectively.

Fig. 2

Bidirectional instrumental variable (IV) estimated association between 25(OH)D and NAFLD by unweighted GRSs. Data were presented as regression coefficient (β) or odds ratio (OR) and 95% confidence interval (CI). In this MR framework, the instrumental variable estimators are ORIV(VD-NAFLD) = exp (ln (ORVD_GRS-NAFLD)/βVD_GRS-25(OH)D) and βIV(NAFLD-VD) = βNAFLD_GRS-25(OH)D/ln (OR NAFLD_GRS-NAFLD). Data were adjusted for age, sex, BMI, current smoking, hypertension, diabetes, HDL-cholesterol, LDL-cholesterol and triglycerides. 25(OH)D, 25-hydroxyvitamin D; VD_GRS, vitamin D genetic risk score; NAFLD_GRS, nonalcoholic fatty liver disease genetic risk score; SD, standard deviation.

Discussion

Vitamin D status and NAFLD are known to be associated in previous studies but whether it is causal, and if so, its causal direction, is still uncertain. In this cross-sectional survey including nearly 10,000 community-dwelling Chinese adults, we examined whether two GRSs composed of SNPs significantly associated with 25(OH)D and NAFLD were associated with NAFLD and 25(OH)D, respectively. Using an MR approach, our findings indicate a causal role of vitamin D in the development of NAFLD may not exist. Conversely, NAFLD may also not induce lower vitamin D status. Our data provided evidence supporting no association between low vitamin D and NAFLD using the MR approach.

Vitamin D deficiency has silently become increasingly more common (Palacios and Gonzalez, 2014), and simultaneously, during the last decade NAFLD has also become the most common cause of chronic liver disease in China and Western countries (Younossi et al., 2016). It is not unexpected that epidemiological studies, largely cross-sectional designs, point towards an association between low vitamin D and the presence of NAFLD and steatohepatitis, independently of confounders such as obesity and insulin resistance (Eliades et al., 2013). In a meta-analysis, patients with NAFLD were 1.26 times more likely to have vitamin D deficiency (Eliades et al., 2013). However, all studies included in this meta-analysis were cross-sectional and retrospective, and moreover there is extremely limited evidence that vitamin D replacement provides clinical benefit for NAFLD. Barchetta et al. found that cholecalciferol did not improve hepatic steatosis in diabetic patients with NAFLD (Barchetta et al., 2016). Patel et al. reported neither low vitamin D nor VD-related genes expressed in liver related to the presence or histologic severity of NAFLD in patients (Patel et al., 2016). These conclusions were based on a large number of patients with biopsy-proven NAFLD and a well-controlled, non-NAFLD group, and this study was the first to utilize hepatic gene expression to confirm or refute the results of a cross-sectional and case-control analysis (Patel et al., 2016). Furthermore, a more recent meta-analysis showed vitamin D supplementation had no effects on metabolic profiles and liver function in patients with NAFLD (Tabrizi et al., 2017).

In MR, genetic variant(s) are used as IVs to assess the causal effect of the exposure on the outcome. In our study, we used GRS instead of each SNP. GRS is a convenient means of adding multiple genetic variants associated with an exposure. Using a GRS as a single IV helps create stronger instruments. The fundamental conditions that GRS should satisfy to be considered an IV should meet three assumptions (Sheehan et al., 2008). First, the GRS should be associated with the exposure. All SNPs selected in this study have previously been shown to be significantly associated with vitamin D status or NAFLD in large GWAS. In our study, the associations of the two GRSs with the two corresponding exposures were also significant. Second, the GRS should not be associated with any confounder of the exposure-outcome association. In this study, we found the two GRSs were not associated with BMI, lipid profile, diabetes and hypertension, which are common potential confounders of the vitamin D-NAFLD association. We further tested the pleiotropic effects of each SNP on the aforementioned confounders and the results showed three NADLF-related SNPs had pleiotropic effects. Thus, we eliminated these three SNP and constructed the NAFLD_GRS4SNP, but the IV estimate did not significantly change from that of NAFLD_GRS7SNP. Third, the GRS is independent of the outcome, except possibly through its association with the exposure. This means that the only causal route from the genetic variants to the outcome is through exposure (there are no other routes between VD_GRS and NAFLD and between NAFLD_GRS and vitamin D) (Sheehan et al., 2008). Based on this assumption, we analyzed the association of each SNP and GRS with the corresponding outcome (vitamin D or NAFLD). Except rs780094, all of the SNPs and GRSs were not significantly associated with the outcome.

In a cross-sectional setting, 25(OH)D was significantly associated with prevalence of NAFLD. It is interesting to observe this discrepancy in the same data. We are prone to suggest this discrepancy may reflect the flaw (residual confounding) of cross-sectional studies. Though the association was significant in the cross-sectional setting, the remaining association will often still be a biased estimate, due to the existence of unknown or unmeasured confounders (sun exposure, immobilization, physical activity, etc.) or imprecision in measured confounders. This discrepancy can be seen in many previous MR studies (Cuellar-Partida et al., 2017, Dimitrakopoulou et al., 2017, Larsson et al., 2017). On the other side, the possibility that factors strongly modifying vitamin D status other than genetics play an influence on NAFLD risk could not be totally excluded. For example, adipose tissue dysfunction, diabetes and consequent systemic inflammation and oxidative stress might link vitamin D deficiency with NAFLD (Cimini et al., 2017, Zhu et al., 2017). MR study suggests that a higher BMI induces lower vitamin D status, while a lower 25(OH)D is unlikely to increase BMI (Vimaleswaran et al., 2013). About diabetes, there are interesting results. Afzal et al. found the DHCR7 allele score was significantly associated with increased risk of any diabetes (P for trend 0.04), which was quite weak in white Danes (Afzal et al., 2014), whereas Ye et al. found this association was unlikely to be causal in populations of European descent (Ye et al., 2015). Thus, we still emphasize the independent MR studies and large prospective trials using more accurate methods of diagnosing NAFLD are needed to validate our findings.

This study also has some limitations. All participants were of Asian origin. They are not directly applicable, although they can serve to represent other ethnicities with good approximation. Second, 25(OH)D was measured only once at baseline. Hence, we were not able to control for intra-individual variability. Third, ultrasound was used to determine liver steatosis. Liver biopsy is not feasible in such a large sample, because MR analysis commonly needs a large sample size. Thus, many large epidemiological studies use ultrasonography to diagnose fatty liver (Tian et al., 2012, Sinn et al., 2017). Saadeh et al.'s criteria to diagnose fatty liver could provide up to 93% sensitivity with a positive predictive value of 62% for the histological diagnosis of NAFLD. Therefore, ultrasonography is a relatively feasible method with acceptable sensitivity and specificity in MR analysis. Fourth, vitamin D supplementation was not accurately available. In our questionnaire there was not a direct question “are you using vitamin D supplementation?”, but we did have a question “what drugs are you using?”. No one mentioned they were using vitamin D supplementation. Combined with the fact that the prevalence of vitamin D deficiency [25(OH)D < 50 nmol/L] was 81.0% in our data, we suspect that the proportion of vitamin D supplementation was very low, which may not have largely affected the results. Finally, we built up our GRSs only based on common variants, which were considered to represent limited NAFLD and vitamin D heritability. We were unable to assess the potential contribution of rare variants.

In conclusion, although the concentration of 25(OH)D was associated with present NAFLD in both directions, using MR analysis there is evidence that there is no causal relationship between genetically determined NAFLD and 25(OH)D and between genetically determined 25(OH)D and NAFLD. Our study suggested long-term vitamin D deficiency may not affect the development of NAFLD, and conversely, NAFLD also may not induce lower vitamin D status. Independent MR studies and large prospective studies are needed to further validate our findings in other cohorts and ethnicities.

Funding Sources

This study was supported by the National Natural Science Foundation of China (81570726, 81600609); Shanghai JiaoTong University School of Medicine (2014); Science and Technology Commission of Shanghai Municipality (16411971200, 16410723200); Commission of Health and Family Planning of Pudong District (PW2015D-5); the Fourth Round of Three-Year Public Health Action Plan of Shanghai by the Shanghai Municipal Commission of Health and Family Planning (15GWZK0202, 20164Y0079); Municipal Human Resources Development Program for Outstanding Young Talents in Medical and Health Sciences in Shanghai (2017YQ053); and Clinical Research Plan of SHDC (16CR3076B). The funders played no role in the design or conduct of the study, collection, management, analysis, or interpretation of data or in the preparation, review, or approval of the article.

Conflicts of Interest

No potential conflicts of interest relevant to this article are reported.

Author Contributions

Yingli Lu designed the research, contributed to the discussion, reviewed and edited the manuscript, and takes full responsibility for the work as a whole. Ningjian Wang designed the study, performed analysis, wrote the manuscript and contributed to the discussion. Chi Chen, Li Zhao conducted the research, analyzed the data, and reviewed and edited the manuscript. Yi Chen, Bing Han, Fangzhen Xia, Jing Cheng and Qin Li conducted the research and contributed to the discussion.