Association between nonalcoholic fatty liver disease and incident diabetes mellitus among Japanese: a retrospective cohort study using propensity score matching

Background Previous studies have demonstrated that nonalcoholic fatty liver disease (NAFLD) is a significant risk factor for diabetes mellitus (DM). However, these studies did not completely determine the relationship between NAFLD and DM due to unbalanced confounding factors. The propensity score (PS) is the conditional probability of having a particular exposure, given a set of baseline measured covariates. Propensity score matching (PSM) analysis could minimise the effects of potential confounders. Thus, this study aimed to use PSM analysis to explore the association between NAFLD and DM in a large Japanese cohort. Methods This retrospective PSM cohort study was performed on 14,280 Japanese participants without DM at baseline in Murakami Memorial Hospital between 2004 and 2015. The independent variable was NAFLD at baseline, and the outcome was the incidence of DM during follow-up. One-to-one PSM revealed 1671 participants with and without NAFLD. A doubly robust estimation method was applied to verify the correlation between NAFLD and DM. Results The risk of developing DM in participants with NAFLD increased by 98% according to the PSM analysis (HR = 1.98, 95% confidence interval [CI]: 1.41–2.80, P < 0.0001). The risk of developing DM in the NAFLD participants was 2.33 times that of the non-NAFLD participants in the PSM cohort after adjusting for the demographic and laboratory biochemical variables (HR = 2.33, 95% CI: 1.63–3.32, P < 0.0001). The participants with NAFLD had a 95% increased risk of DM after adjusting for PS (HR = 1.95, 95% CI: 1.39–2.75, P = 0.0001). All potential confounding variables were not significantly associated with NAFLD and DM after PSM in the subgroup analysis. In the sensitivity analysis, the participants with NAFLD had a 2.17-fold higher risk of developing DM in the original cohort (HR = 2.17, 95% CI: 1.63–2.88, P < 0.0001) and were 2.27-fold more likely to develop DM in the weighted cohort (HR = 2.27, 95% CI: 1.91–2.69, P < 0.00001). Conclusions NAFLD was an independent risk factor for the development of DM. The risk of developing DM in the NAFLD participants was 2.33 times that of the non-NAFLD participants in the PSM cohort after adjusting for the demographic and laboratory biochemical variables. The participants with NAFLD had a 95% increased risk of DM after adjusting for PS.


Introduction
Diabetes mellitus (DM) has become a serious global public health problem. According to international epidemiological research on DM, the prevalence of DM in 2019 was 9.3% (approximately 500 million people) [1]. DM and its complications can seriously affect the health of patients and increase medical costs, which can lead to a heavy economic burden on the patients and society [2]. DM is a metabolic disease characterised by hyperglycaemia caused by insufficient insulin secretion or insulin resistance (IR) [3]. Many researchers have explored the pathogenesis and risk factors of DM.
Some prospective cohort studies recently reported that nonalcoholic fatty liver disease (NAFLD) is a significant risk factor for DM [4,5]. NAFLD is often accompanied by DM, obesity, and hyperlipidaemia [6,7]. Additionally, some studies found that NAFLD was an independent risk factor for DM after adjusting for confounding variables [8,9]. A recent meta-analysis of 33 studies involving more than 500,000 individuals showed that participants with NAFLD had a 1.19-fold higher risk of developing DM than participants without NAFLD [10].
The propensity score (PS) is defined as the conditional probability of having a particular exposure (NAFLD versus non-NAFLD), given a set of baseline measured covariates. The propensity score matching (PSM) method is useful in studies in which there are many covariates potentially confounding a rare outcome, there is potential confounding by indication, and there are resource constraints that prevent the conduction of randomized clinical trials. Given the numerous potential confounding variables, a traditional parsimonious regression model could result in bias due to unmeasured or residual confounding, whereas the inclusion of more variables could result in overfitting of the model, potentially preventing identification of the association between the exposure of interest and the outcome [11]. Therefore, PSM analysis was used in this study to explore the actual association between NAFLD and DM in the NAGALA (NAfld in the Gifu Area, Longitudinal Analysis) database of 14,280 Japanese people.

Study design and data source
This was a secondary retrospective study based on NAGALA, sourced from the public DRYAD database (www.Datadryad.org.database). Raw data were provided by Okamura et al. [12]. The original study included 20, 944 participants who underwent medical examinations at Murakami Memorial Hospital from 2004 to 2015. All participants completed a detailed questionnaire on their demographic characteristics and health behaviours. A trained staff member measured the demographic data, such as body weight and waist circumference (WC). Data on laboratory-related biochemical parameters were collected under standardised conditions and processed using a unified process. Since this was a retrospective cohort study, the risk of selection and observation biases was reduced.
The original research was approved by the ethics committee of Murakami Memorial Hospital, and informed consent was obtained from all participants. The authors of the original research handed over all copyrights of these data. Therefore, this study performed a secondary analysis based on their data without prejudice to the authors' rights.

Study sample
In the original study, 5480 participants were excluded from 20,944 Japanese participants based on the following criteria: (1) viral hepatitis (defined by measurements of hepatitis B antigen and hepatitis C antibody at baseline), (2) alcoholic fatty liver disease, (3) DM at baseline, (4) fasting plasma glucose (FPG) level of ≥6.1 mmol/L, (5) use of any medication at baseline, and (6) missing covariate data. Therefore, 15,464 participants were included in the original study. This study further excluded 1184 participants with excessive alcohol consumption (alcohol consumption > 210 g/week in males and > 140 g/week in females [13]). Finally, this study included 14,280 eligible participants. Figure 1 detailed the selection process for all the participants.

Independent variable and covariates
The independent variable was baseline NAFLD, which was diagnosed using abdominal ultrasonography performed by trained technicians [12]. The following covariates were extracted at baseline: age, gender, WC, body mass index (BMI), alcohol consumption, smoking status, regular exerciser, systolic blood pressure (SBP), diastolic blood pressure (DBP), aspartate aminotransferase (AST), alanine aminotransferase (ALT), total cholesterol (TC), gamma-glutamyl transferase (GGT), glycosylated haemoglobin A1c (HbA1c), FPG, high-density lipoprotein cholesterol (HDL-C), and triglycerides (TG). Alcohol consumption was classified into three categories: no or very little alcohol consumption (less than 40 g of alcohol per week), light alcohol consumption (40-140 g of alcohol per week), and moderate alcohol consumption (140-210 g of alcohol per week) [14]. Participants who regularly performed any type of exercise at least once a week were defined as regular exercisers [15]. Visceral fat obesity was defined as WC ≥ 90 cm in males or ≥ 80 cm in females [16].

Outcome measure
The outcome was the incidence of DM. DM was defined as HbA1c ≥ 6.5%, FPG ≥ 7 mmol/L [17], or self-reported during follow-up.

Statistical analyses
Continuous variables conforming to the normal distribution were presented as mean ± standard deviation (SD), while continuous variables conforming to the skewed distribution were expressed as median and quaternary ranges (25-75th percentile). Categorical variables were expressed as frequencies and percentages. The one-way ANOVA, the Kruskal-Wallis H test and the chi-square test were performed to detect differences between the groups. Missing values of HDL-C were handled by supplementing them with the mean. PSM analysis was used to match the baseline characteristics between the NAFLD and non-NAFLD groups (Table 1), and to form a single group of participants with similar baseline characteristics. The non-parsimonious multivariable logistic regression model was performed to calculate the PS based on NAFLD as the independent variable and 17 baseline variables as covariates. This study used a 1:1 matching protocol without replacement (greedy matching algorithm), and the calliper width was equal to 0.01. The evaluation index of the balance between groups was the standardized differences [18,19]. If the standardized differences were less than 10.0%, the covariates between the two groups were considered to be well balanced [18,19]. Besides, the Kaplan-Meier method was used to assess the incidence of DM in each group, and the log-rank test was conducted to determine significance. P values were calculated for each pair of groups (total three comparisons: Low PS vs. Medium PS, Low PS vs. High PS, Medium PS vs. High PS), with Bonferroni correction [20]. The Cox proportional-hazards regression model was performed to explore the relationship between NAFLD and the incidence of DM by adjusting for covariates in the PSM cohort. The doubly robust estimation method, which combines PS models and the multivariate regression model, was applied to verify the association between NAFLD and the incidence of DM [21,22]. Prespecified subgroup analyses were conducted based on gender, WC, BMI, AST, ALT, TC, GGT, HbA1c, FPG, HDL-C, TG, and PS. Specifically, continuous variables were converted to categorical variables based on the clinical cut-off point or median. Each stratification was adjusted for all factors, except for the stratification factor. In the subgroup analyses, only the corresponding matched pairs in the same subgroup were chosen to maintain the balance of baseline characteristics between the NAFLD and non-NAFLD groups. For example, in the subgroup of participants with BMI < 25 kg/m 2 , only when the matched pairs of the NAFLD and non-NAFLD groups both belonged to the BMI < 25 kg/ m 2 subgroup, these participants could be included in the subgroup analysis. Likelihood ratio tests were used to inspect the modifications and interactions of the subgroups.
For sensitivity analyses, the inverse probability of treatment weights (IPTW) was calculated using the estimated PS. For instance, the weight of NAFLD participants was 1/PS, and the weight of non-NAFLD participants was 1/ (1 -PS). The IPTW model was conducted to create a weighted cohort [22]. In the sensitivity analysis, two association inference models were added to the original and weighted cohorts. A series of sensitivity analysis methods were used to test the robustness of the findings of the study and how conclusions could be affected by applying different association inference models. The effect sizes and P-values were calculated in all models. The results of this study were reported following the STROBE statement [23].
The current research analysis was performed using Empower-Stats (http://www.empowerstats.com, X&Y Solutions, Inc., Boston, MA) and the statistical software package R (http://www.R-project.org, The R Foundation). A two-sided P < 0.05 was considered significant.

Study population
A total of 14,280 participants were eventually included in this study, including 52.10% men and 47.90% women (Fig. 1). Among them, 2515 (17.61%) participants suffered from NAFLD, and 11,765 (82.39%) did not suffer from NAFLD. The average age of the study population was 43.53 ± 8.89 years. During a mean follow-up of 2207.02 ± 1376.51 days, 324 participants developed DM. Some baseline characteristics showed statistically significant differences between the NAFLD and non-NAFLD groups before PSM. Higher levels of age, BMI, WC, SBP, DBP, FPG, HbA1c, AST, ALT, GGT, TC, and TG were observed in the NAFLD group. Participants with NAFLD showed a higher proportion of males, ever smoker, and current smoker. However, participants with non-NAFLD had higher HDL-C levels and higher rate of regular exerciser. In total, 1671 NAFLD patients were matched with 1671 non-NAFLD subjects by using one-to-one PSM. The standardized differences of all covariates were less than 10.0% after PSM, showing a good match. In other words, the differences in baseline characteristics between the two groups were minimal.

The incidence of DM
The results of the Kaplan-Meier analysis revealed that the cumulative incidence of DM among the participants with NAFLD was significantly higher than that among participants without NAFLD before PSM (P < 0.0001; Fig. 2a). This difference still existed in the PSM cohort (P < 0.0001; Fig. 2b). Moreover, the cumulative incidence of DM was significantly higher in participants with higher PS after Bonferroni correction (Fig. 3).

Association between NAFLD and the incidence of DM
The Cox proportional hazards regression model was applied to assess the association between NAFLD and DM risk in the PSM cohort.

Subgroup analysis
Subgroup analysis was applied to discover potential confounding variables that might have affected the association between NAFLD and DM risk. Gender, BMI, WC, Values were n (%) or mean ± SD or median (interquartile range: 25th to 75th percentiles) SD standard deviation, BMI body mass index, WC waist circumference, SBP systolic blood pressure, DBP diastolic blood pressure, FPG fasting plasma glucose, HbA1c glycosylated haemoglobin, ALT alanine aminotransferase, AST aspartate aminotransferase, GGT gamma-glutamyl transferase, TC total cholesterol, TG triglyceride, HDL-C high-density lipoprotein cholesterol a b  Table 3 showed that none of the interactions were observed based on the prior specifications. The results revealed that the variables listed above did not affect the association between NAFLD and DM risk after PSM.

Sensitivity analysis
The estimated PS was used to generate a weighted cohort by establishing an IPTW model. This study evaluated the association between NAFLD and the incidence of DM in both the original and weighted cohorts. Moreover, the unadjusted, partially adjusted, and fully adjusted models were established in both cohorts (

Discussion
The PSM cohort study showed that NAFLD was an independent risk factor for the development of DM. The risk of developing DM in the NAFLD participants was 2.33 times that of the non-NAFLD participants in the PSM cohort after adjusting for the demographic and laboratory biochemical variables. This figure decreased to 95% after adjusting for the PS. In the subgroup analysis, no interaction was observed, indicating that the relationship between NAFLD and DM was robust. The correlation also existed in both the original and weighted cohorts.  NAFLD can develop into liver fibrosis, cirrhosis, and liver cancer and increase the risk of developing diabetes and cardiovascular diseases [24]. Patients with NAFLD have been reported to have a higher prevalence of prediabetes/DM and increased IR [9,25]. The incidence of DM in NAFLD patients was also higher in NAFLD participants than in non-NAFLD patients, even if their plasma glucose levels were within normal ranges [26]. The improvement of NAFLD was related to a decrease in the incidence of DM [27]. It has been reported that NAFLD and DM have some same risk factors and often occur simultaneously in one person [6,10]. Meanwhile, several prospective studies have found that NAFLD strongly increases the incidence of DM [28,29]. In addition, a study explored the relationship between NAFLD and DM using PSM methods [30]. Their findings suggested that NAFLD was a risk factor for DM, which was consistent with the conclusion of this study. However, that study excluded participants with other metabolic diseases (hypertension, dyslipidaemia), and NAFLD was mainly diagnosed by non-invasive scores [30]. Based on these findings, the prevalence of NAFLD might be underestimated. Therefore, the results of the study mentioned above could not be applied to the general population. In this study, NAFLD was diagnosed by abdominal ultrasonography, and the study population was more extensive. These results could better reflect the actual relationship between NAFLD and DM. In contrast, some studies showed different findings. They showed that the association between NAFLD and the risk of DM was not significant after adjusting for confounding factors [31,32]. The possible reasons for these inconsistent findings were as follows: (1) The study population was diverse, including different races, genders, ethnicities, ages, and so on. In this large-scale cohort study, the risk of developing DM in the NAFLD participants was 2.33 times that of the non-NAFLD participants after PSM. There was a difference between this study and previous studies in terms of the risk of DM, which might be related to the fact that this study conducted a PSM analysis and effectively controlled for more confounding variables which were well known to be related to NAFLD and DM, including age, gender, BMI, WC, smoking status, alcohol consumption, regular exerciser, SBP, DBP, ALT, AST, GGT, HbA1c, FPG, TC, TG, and HDL-C [8,33]. Additionally, this study was based on large cohort data (14,280 participants), which further strengthened the statistical power of the results. Exploring the association between NAFLD and DM can help us better guide patients in clinical practice and develop management strategies to reduce DM risk [34,35].. The mechanism by which NAFLD leads to DM remains unclear. A study demonstrated that NAFLD could cause IR, which could further mediate the development of DM [36]. The mechanisms by which NAFLD contributes to IR are as follows: (1) Adipose tissue dysfunction and inflammation promote the secretion of adipokines, increase the secretion of pro-inflammatory factors (such as tumor necrosis factor-α), and increase the release of free fatty acids, resulting in decreased insulin sensitivity. Adipose tissue dysfunction and inflammation interfere with the activation of the pro-inflammatory pathway of insulin signal transmission, leading to decreased insulin sensitivity [37]. (2) Certain incretin related to NAFLD can directly inhibit the production of endogenous glucose through an insulin-dependent mechanism [38]. The reduction of these incretin effects also leads to IR [38].

Study strengths and limitations
This study has the following strengths. The most innovative part of this study is that PSM was used to explore the relationship between NAFLD and the risk of developing DM. In recent years, the PSM method has been widely used in observational research. The acknowledged advantages of the PSM method include a wide range of data requirements, including a reduction of inter-group differences, balancing inter-group confounders, and achieving the effect of "similar randomization". Subgroup analyses were performed to explore other potential risk factors that could affect the association between NAFLD and DM. A series of sensitivity analyses were conducted to ensure the robustness of the results. This study mainly used IPTW to establish a weighted cohort and further explore the association between NAFLD and the incidence of DM in the weighted cohort. More importantly, the sample size of the participants in this study was more extensive than that in most previous retrospective cohort studies. However, the current study has several limitations. First, the population included in this study was Japanese, and therefore, the generalizability of these results to other races requires further validation. Second, the lack of a 2-h oral glucose tolerance test in the original study might have underestimated the incidence of DM. However, it is not feasible to conduct a 2-h oral glucose tolerance test in such a large cohort. Third, the PSM could balance known confounding variables as much as possible, but it could not ensure that all measured baseline characteristics were matched and consider the influence of unknown variables. To reduce the interference of variables on the measurement results, the calliper width was set at 0.01. Fourth, ultrasonography may have some limitations in diagnosing NAFLD. However, some noninvasive scores, such as the FIB4 score, have some advantages. Considering that the original data lack relevant data, such as platelets, FIB4 scores could not be used to diagnose NAFLD. In the future, it would be worthwhile to design studies or collaborate with other researchers to collect as many variables as possible to analyse the actual relationship between the non-invasive score of NAFLD and DM. Fifth, the differences between type 1 and type 2 DM were not considered in the present study. However, type 2 DM is most common, accounting for over 90% of the cases of DM [41]. Therefore, this study aimed to explore the relationship between NAFLD and type 2 DM.

Conclusions
NAFLD was an independent risk factor for the development of DM. After adjusting for the demographic and laboratory biochemical variables, the risk of developing DM in the NAFLD participants was 2.33 times that of the non-NAFLD participants in the PSM cohort. The participants with NAFLD had a 95% increased risk of DM after adjusting for PS.