Nomogram for Prediction of Diabetic Retinopathy Among Type 2 Diabetes Population in Xinjiang, China

Introduction

The prevalence of diabetes is explosive worldwide, especially in developing countries.¹ In 2019, the Diabetes Atlas of International Diabetes Federation (IDF) (9th Edition) indicated that about 463 million patients exist over the world, which is expected to rise to 700 million in 2045.² In China, the number of diabetics has reached 116.4 million, takes up about one third of the amount of diabetics worldwide.³ Diabetic retinopathy (DR), as one of the complications of diabetes, is considered to be one of the main causes of vision loss and visual impairment in adults, which can lead to blindness in severe cases.^4,5 It’s difficult to detect the symptoms of DR in early stage. Blindness often suddenly occurs. Unfortunately, there is no effective treatment for DR to recover patients’ vision. The incidence rate of DR increases with the increase of DM and the prolongation of lifespan. In recent years, a meta-analysis of 35 epidemiological studies showed that the prevalence of DR in DM patients was 34.6%. It is estimated that the number of DR patients will increase from 127 million in 2010 to 191 million in 2030.⁶ Therefore, the key to reduce diabetic blindness lies in the early detection of asymptomatic DR, so that patients can get effective treatment before blindness occurs.

The current clinical diagnosis of DR is highly subjective description reported as vision loss by patients and regular fundus examination. Once diagnosed, it is difficult to reverse. It is necessary to find the early markers of DR. Many epidemiological studies on DR, including cross-sectional studies^7,8 and cohort studies,^9,10 are carried out around the world to explore the risk factors associated with the disease so as to prevent and treat the disease. Older T2DM patients with longer disease course tend to have higher risk of DR. In addition, diabetic renal complications, poor glycemic control, hyperlipidemia or hypertension have also been reported as risk factors for DR, affecting the progress of DR.¹¹ Early treatment of DR was used to evaluate the effects of these reported factors.¹² It was found that the duration of T2DM, hyperglycemia and hypertension were the most important risk factors for the progression of visual loss.¹³ However, little attention has been paid to DR, related risk factors, as well as eye screening.¹⁴ Therefore, more research are needed to improve the understanding of risk factors.

At present, more and more researchers have made a lot of efforts to build prediction models or improve existing prediction models.^15–17 among which building nomogram is one of the hottest research directions. In the latest research,^15,18,19 the use of nomogram prediction has shown good performance.

In this paper, totally 13,980 T2DM cases were included for exploring the potential relationship between DR and its dangerous factors, as well as establishing a nomogram model benefiting clinical workers for early diagnose of DR patients.

Materials and Methods

Study Design and Participants

The clinical data of 18,904 patients with T2DM in Xinjiang from 2010 to 2019 was retrospectively collected. The exclusion criteria included, (1) Age below 10 years old (a progressive insulin secretory defect on the background of insulin resistance);²⁰ (2) incomplete clinical data; (3) patients with cancer or cardiovascular disease; (4) abnormal data. Through strict data filtering and preprocessing, 13,980 samples were obtained (Figure 1). All subjects were randomly divided into development group (n = 10,485) and validation group (n = 3495) according to the ratio of 3:1. We used retinopathy as the outcome variable. Nomogram was constructed based on development group and was verified using validation group. The study was approved by the Medical Ethics Committee of the First Affiliated Hospital of Xinjiang Medical University, and written informed consent was obtained for each subject. The study was performed in accordance with the Declaration of Helsinki.

Figure 1 Flow diagram of the data screening. Notes: diagnosis unrecognized: Comprehensive eye examinations without mydriasis; incomplete clinical data: Delete samples missing more than 2 variables; DR did not meet the requirements: retinopathy caused by hypertension; abnormal data: data that cannot be included in the evaluation due to severe cataracts and other reasons; 76 patients in incomplete clinical data and DR did not meet the requirements were repeated.

Information Collection

General data, gender, Age, Duration of T2DM, diabetic peripheral neuropathy, Hypertension (The diagnostic criteria for diabetic peripheral neuropathy and Hypertension was based on the Standards of medical care for type 2 diabetes in China (2010)).²¹

DPN, (1) Abnormal temperature perception; (2) Abnormal vibration perception; (3) Nylon filament examination, foot sensation decreased or disappeared; (4) Ankle reflex disappeared; (5) Nerve conduction velocity (NCV). If someone has two or more of the above phenomena, he/she is diagnosed as DPN;

Hypertension: blood pressure ≥140/90 mmHg or on therapy for hypertension.

Biochemical data: the patient’s fasting venous blood was drawn from the morning after admission and sent to the laboratory of our hospital for determination of white blood cell, neutrophilic granulocyte, eosinophil, lymphocyte, total cholesterol, high-density lipoprotein, low-density lipoprotein, aspartate aminotransferase, alanine aminotransferase, glucose, Hemoglobin A1c, Glycosylated Serum Protein, Apolipoprotein A1, apolipoprotein B, triglyceride, Urea.

Fundus examination, all patients were examined by ophthalmologists for vision, intraocular pressure, optometry, ophthalmoscope, and then taken mydriatic fundus photography. For those who were not clearly staged, fundus angiography was performed to determine whether there was retinopathy.²¹

Statistical Analysis

Firstly, we deleted 20 variables which makes no sense in this study (No sense variables included admission time, discharge time, hospitalization days, and pre-hospital inspection indicators (the inspection indicators on the second day of hospitalization were meaningful variables)). Secondly, 4924 samples that did not meet the requirements of this study were excluded. Finally, we used multivariate multiple imputation to impute missing values of 13,980 cases to maximise statistical power (Figure 1).

Statistics description, when describing the baseline characteristics, the count data were described by frequency and proportion, and chi square test was used to compare the two groups. The measurement data with normal distribution was statistically described with mean ± standard deviation, and t-test was used to compared the data between the two groups. The measurement data of skewed distribution was described by median and interquartile distance, and compared between the two groups by non-parametric test.

Establishment of nomogram model, according to the collected clinical data and research purposes, Considering that there may be multicollinearity between variables, which will cause the unstable in the estimation of regression coefficient and intercept coefficient and result in the instability of the model. Therefore, in this study, the minimum absolute contraction and lasso regression are used to select the most useful candidate predictors, as well as improve the performance of the model,^22,23 and the Odds ratio (OR),95% CI and P values of each candidate were calculated to predict the possible diagnosis. Finally, both statistically significant variables and clinically significant variables were included in the multivariate logistic regression analysis, and the nomogram model was constructed using the statistically significant variables selected at last (P < 0.05).^15,24,25

Effectiveness evaluation of nomogram model, the effectiveness evaluation of prediction model mainly includes discrimination ability (discrimination), accuracy (calibration curve) and clinical practicability. Discrimination and calibration are commonly used to evaluate the prediction model. Discrimination refers to the ability to distinguish high-risk and low-risk groups. AUC, also known as C statistic, is used to evaluate the discrimination, whose value is between 0.5 and 1. The closer the AUC value is to 1, the better the prediction model is.²⁶ It is generally believed that AUC greater than 0.8 is of great diagnostic value, and the specific predictive value should be evaluated in combination with clinical practice. Meanwhile, The nomogram was subjected to bootstrapping validation (1000 bootstrap resamples) to calculate an optimism-corrected Harrell’s C statistic (C-index)^19,27 in the development group and the validation group, respectively. Calibration degree is an important index to evaluate the accuracy of disease risk model in predicting the probability of individual future outcome events. Calibration degree, also called consistency, reflects whether the predicted risk is consistent with the actual risk. Akaike information criterion (AIC) and Hosmer-Lemeshow of fit test (H-L test) were also used to evaluate the consistency of the model. The smaller the AIC value is, the closer the p value of H-L test is to 1, indicating that there is no significant difference between the probability of model prediction and the probability of actual observation. Smaller AIC represents a better calibration degree of the model.²⁸ Finally, we conducted a decision curve analysis to evaluate the clinical effectiveness of the model.²⁹ The results showed that the proportion of true positive minus the proportion of false positive, and then weighed the relative risks of false positive and false negative results to obtain the net benefits of making decisions.

Statistical analysis and establishment of nomogram model were performed using R software (version 4.0.2).

Results

Baseline Characteristics of DR

Totally 13,980 T2DM patients were recruited in this study, including 8071 patients without DR and 5909 with DR. Table 1 depicted the baseline information. Compared with T2DM patients without DR, DR patients were older, with longer duration of T2DM and higher HbA1C. Females patients with DR had higher prevalence of DPN and hypertension than those without DR. However, indicators as EOS, Tc, HDL, LDL, GSP presented no differences between DR patients and T2DM patients without DR.

Table 1 Clinical Characteristics of the Study Population (N=13,980)

Characteristics Selection

In this study, the least absolute shrinkage and selection operator (Lasso) regression was applied to the development group. With the change of the penalty coefficient λ, the variables included in the model were gradually reduced. Finally, the 10-fold cross-validation error was selected as the minimum λ+1 (lambda. 1se=0.012) as the optimal value of the model, at this time we get 8 variables (Figure 2). Then the 8 variables obtained from Lasso regression were included in the logistic multiple regression model, including DPN, Age, WBC, NE, HDL, HbA1C, Duration of T2DM and GSP (details in Table 2).

Table 2 Baseline Characteristics Were Analyzed by Multivariate Logistic Regression in the Training Group. (N= 10,485)

Figure 2 Demographic and clinical feature selection using the LASSO binary logistic regression model (figure was created by R software, “glmnet” package, version 2.0–18, https://CRAN.R-project.org/package=glmnet). Notes: (A) Optimal candidate (lambda) selection in the LASSO model used 5-fold cross validation via minimum criteria. The area under the receiver operation characteristic curve was plotted versus log (lambda). Dotted vertical lines were drawn at the optimal values by using the minimum criteria and the 1 SE of the minimum criteria; (B) LASSO coefficient pro-files of the 21 candidates. A coefficient profile plot was produced against the log (lambda) sequence. Vertical line was drawn at the value selected using 5-fold cross validation, where optimal lambda resulted in 8 candidates with nonzero coefficients.

Candidate variables were selected by Lasso regression. Finally, the candidate predictive variables were included in multivariate logistic regression analysis. Then OR, 95% CI and P value of each candidate were calculated. Except WBC ((p>0.05), DPN, Age, NE, HDL, HbA1C, Duration of T2DM and GSP were all statistically significant (Table 2). There was no multiple collinearity among the risk factors included in the model, and the maximum variance expansion factor (VIF) was 1.395 while the lowest was 1.060. Therefore, we construct the clinical risk prediction model with these 7 statistically significant variables.

Development of a DR-Predicting Nomogram

Nomogram was established based on multiple regression (Figure 3). Each sub-type in these variables is assigned a score, DPN was 61, Age was 79, NE was 68, HDL was 57, HbA1C was 77, Duration of T2DM was 100, GSP was 27. The total score ranged from 153 to 255, and the corresponding risk rate ranged from 0.1 to 0.9. The higher the total score, the higher the risk of DR.

Figure 3 Nomogram to predict the risk of DR. Notes: To use the nomogram, an individual participants value is located on each variable axis, and a line is drawn upward to determine the number of points received for each variable value. The sum of these numbers is located on the total points axis to determine the risk of DR.

Take an example of nomogram usage, a sample without DPN was randomly selected,47-year-old T2DM patients, the NE, HDL, HbA1C, Duration of T2DM and GSP were 3.34 mmol/L, 1.13 mmol/L,9 mmol/L, 3 years and 3.4 mmol/L, respectively. A total score of this sample is calculated to 272 (Uncorrected). Hence, the corresponding probability of DR was 12% (Figure 4A).

Figure 4 According to Nomogram’s estimate, the probability of DR risk in patient no.1020 (A) and no.2089 (B) was 0.12 and 0.816, respectively, P < 0.001. Notes: Figure a and b are the line segment dynamic nomogram, which is different from the line segment static nomogram in the paper. After optimization corrected, each variable starts with 0 (nomogram). In the line segment dynamic nomogram, it is uncorrected and reflects the score corresponding to each variable in the actual situation. Therefore, the score corresponding to each variable starts from the specific value rather than 0. Although the scores of the two are different, the prediction probability is the same.

Another sample suffering DPN was randomly selected, 50-year-old T2DM patients, the NE, HDL, HbA1C, Duration of T2DM and GSP were 4.2 mmol/L,0.6 mmol/L,11 mmol/L,10 years and 3.96 mmol/L, respectively. A total score of this sample is 332 (Uncorrected). The corresponding probability of DR was 81.6% (Figure 4B).

Validation of Nomogram for T2DM

The model quality was evaluated from three aspects, discrimination, consistency (calibration) and clinical effectiveness. The ROC curve was drawn for the prediction effect of nomogram (Figure 5). In the development group, the AUC was 0.882 (95% CI, 0.875–0.888); according to the Youden-index (sensitivity + specificity − 1), the maximum of which was selected as the theoretical optimal diagnostic cut-off point, which was 0.528, (specificity, 0.879, sensitivity, 0.779). In the validation group, the AUC was 0.870 (95% CI, 0.856–0.881) and the theoretical optimal cut-off point was 0.487 when the Youden-index was the largest (specificity, 0.848, sensitivity: 0.804). In addition, the AUC value of the validation group was only 0.012 lower than that of the development group, indicating that the prediction model has good prediction discrimination in both groups. Meanwhile, the nomogram were subjected to bootstrap-ping validation (1000 bootstrap resamples) to calculate an optimism-corrected C-index were 0.878 and 0.867 in the development group and the validation group, respectively. In addition, in the development group and the validation group, we conducted ROC analysis to evaluate the performance of seven risk factors in DR risk prediction. As shown in Figure 5, in the prediction capacity of DR risk, the combination model is superior to the single risk factor model.

Figure 5 The ROC curves of the nomogram for DR risk (left, development group, right, validation group. including ROC curves of single risk factor model). Notes: In the development group, the AUCs were 0.882 (95% CI: 0.875–0.888); In the validation group, the AUCs were 0.870 (95% CI: 0.856–0.881), respectively. ROC, receiver operating characteristics curve, AUC, area under curve.

Figure 6 depicted that in both development and validation groups, the calibration curve of nomogram model showed great consistence between the prediction values and actual values, while Hosmer-Lemeshow test showed the model was in line with observed data (P >0.05) as well.

Figure 6 Calibration curves for the validation and development group models (left, development group, right, validation group). Notes: The diagonal dotted line represents a perfect prediction by an ideal model. The solid line represents the performance of the nomogram, of which a closer fit to the diagonal dotted line represents a better prediction.

Results of DCA curve in development and validation groups were shown in Figure 7. The gray line represented the model. The blue line represented the net benefit when all participants had DR, while the red line represented the net benefit when no participant had DR. The area between the “red line” and “blue line” of the model curve meant the clinical validity of the model. If the gray line is above the horizontal red line and the left oblique blue line, we can think that the gray line value of this section can benefit. The DCA curve showed that using the nomogram model to predict the occurrence of DR risk had a higher net income. The nomogram model is a good assessment tool because of available selectivity threshold probability.

Figure 7 The decision curve analysis of the nomogram for DR risk (left, development group, right, validation group). Notes: The red line represents the net benefit when no participant was considered to exhibit DR, while the blue line represents the net benefit when all participants were considered to suffer from DR. The area among the model curve, “treat none line” (red line) and “treat all line” (blue line), represents the clinical usefulness of the model. The farther the model curve is to the blue and red lines, the better clinical value the nomogram holds.

Discussion

This study is a retrospective research based on T2DM patients in Xinjiang Uygur Autonomous Region, aiming to predict the risk of DR among whom. DPN, Age, NE, HDL, HbA1C, Duration of T2DM and GSP were independent risk factors of DR. DPN, Age, HDL, HbA1C, Duration of T2DM and GSP were positively correlated with DR, having the strongest correlation with DPN, Duration of T2DM, HbA1C, while these factors were negatively correlated with NE. Original data was randomly divided into development (n=10,485) and validation (n=3495) groups. The verification results showed both of which had good risk prediction ability. Calibration chart showed that the nomogram was accurate in predicting the risk of T2DM; similarly, our decision curve analysis also proved the clinical application value of the nomogram.

DPN is one of the most common chronic complications of T2DM, which often coexists with DR in patients with T2DM. The prevalence of DPN increases with the prevalence of DR, both in parallel,³⁰ which are consistent with the results of this study. DR and DPN are two manifestations of microvascular complications in the eye and peripheral nerves, respectively, and the pathological basis of both is related to microangiopathy, microcirculatory disorders and glucose metabolism. Both pathologies are related to microangiopathy, microcirculatory disorders and glucose metabolism disorders. The gradual thickening of the capillary basement membrane after various oxidative stresses is the common morphological basis of both.³¹ The strong correlation between them suggests that we should pay higher attention to DR patients in our clinical work, and we should not only pay attention to their ocular lesions, but also carefully inquire about neurological symptoms to exclude DPN, and also perform detailed ophthalmologic examinations in time for patients with autonomic dysfunction, neuralgia, limb sensory numbness and mobility disorders, so that they can be detected early and treated in time.³²

Among the independent risk factors in this DR prediction model, the duration of T2DM is the factor with the greatest risk weight. Duration of T2DM is also recognized as the most important risk factor for diabetes complications. It is closely associated with the development of various macrovascular and microvascular complications of diabetes. Liu et al³³ found in a 16-year cohort study of patients with T2DM that the prevalence of DR was 25% at 10 years of diabetes and increases to 50.0% at 15 years of disease progression. In the 35 epidemiological studies reported by Yau et al showed a significant increase in the prevalence of DR with increasing disease duration, with a 21.1% prevalence of DR in patients with less than 10 years of diabetes and a 76.3% prevalence of DR when the disease duration exceeded 20 years.³⁴ The results of this study showed that the risk of developing DR in patients with diabetes mellitus of more than 15 years duration was 3.45 times higher than that in patients with disease duration less than or equal to 5 years (Table 1). As the duration of diabetes increases, the risk of DR increases significantly, which is consistent with the domestic and foreign literature.^33–35 According to Calderon et al³⁶ and Solomon,³⁷ the blood glucose index is closely related to the severity of DR, mainly including HbA1c, HDL and GSP. The higher the HbA1C level in T2DM patients is, the higher the incidence of DR is. When HbA1C level decreased, the risk of DR also decreased, consistent with the results of the present study. In addition to glycemic indicators, some clinical studies found that earlier diabetes onset (at younger age) was associated with the occurrence of DR.^38,39 Previous epidemiological studies showed an effect of age on DR. Unlike Cheung et al and Sasongko’s et al,^40,41 who found a negative association between age and DR, we demonstrated that age as an independent risk factor was positively associated with DR. Older patients showed higher likelihood of DR. This difference might be explained by the fact that Cheung et al conducted their study in a Healthy population, whereas Sasongko’s study included only type 1 diabetic patients aged 12–20 years. In addition, these studies were conducted in a different ethnic group compared to our patients. Earlier it has been demonstrated⁴² that there were differences in retinal vascular curvature values between races.

There are several risk prediction models for DR based on different ethnic populations, and the predictors in the models vary. The main factors include age, duration of T2DM, blood pressure, urine protein, and blood glucose etc. Referring to the model developed by Liu Y et al¹² the age of patients, course of disease, SBP, HbA1C, PBG and Triglyceride were incorporated into the DR risk prediction model. Results showed that the AUC was 0.661 and 0.658 for the development and validation groups, respectively (Supplementary Figure S1). Soto-Pedre et al⁴³ developed a model based on a Spanish population to evaluate whether retinopathy endangers patients’ vision, with predictors including gender, type of diabetes, duration of diabetes, glycosylated hemoglobin, blood pressure, and grading of DR, with AUC of 0.74. Mendoza-Herrera et al⁴⁴ carried out a study among a Mexican population to predict DR, with model predictors including the duration of T2DM, exercise condition, blood glucose level, systolic blood pressure, with AUC of 0.78. Cichosz et al⁴⁵ constructed a risk prediction model for DR based on a white population that included glycated hemoglobin, BMI, waist circumference, age of patients, systolic blood pressure, urine protein, and urine creatinine, with AUC of 0.74. Because of the heterogeneity between races, risk prediction models based on different ethnic groups often could not be applied to the whole population and have limited predictive validity.

The capillaries of diabetic patients are more sensitive to blood glucose changes in a short period of time,⁴⁶ and GSP can just reflect the average blood glucose concentration in the past 1 to 3 weeks. Studies⁴⁷ have confirmed that NE played an important role in retinopathy, but there were few studies about predicting the risk of DR with GSP and NE involved. Previous studies mainly included systemic factors of diabetes and related indicators of hypertension and renal lesions, such as age, duration of T2DM, HDL, HbA1C et al. In this study, the two detection indicators of GSP and NE were incorporated into the DR risk prediction model simultaneously, finding that the area under the ROC curve of the training group was 0.882 which far exceeded the model established by Liu Y et al¹² (Supplementary Figure S1). Therefore, we thought that NE and GSP played a certain role in predicting the risk of DR. However, although the data used to validate the prediction model in this study were independent of the modeling data, case selection bias might exist since the data was from the same Medical center. In the follow-up study, we plan to collaborate with other centers and use their clinical data to conduct a more extensive and in-depth external validation of the prediction model. Overall, much further work was need to build the prediction model to clinical work.

Conclusion

Diabetic retinopathy is one of the most common and serious microvascular complications of diabetes. While the clinical diagnosis of DR mainly relies on regular fundus examinations and subjective patient descriptions, once diagnosed, it is difficult to reverse. Therefore, it is necessary to search for risk factors for diabetic retinopathy. In this study, we present a nomogram model applicable to Xinjiang region in western China. With the visual visualization of nomogram model, it can help clinical workers to provide reference for the screening and early diagnosis of DR by including patient information into the risk prediction model and calculating the risk score to assess the risk of DR in patients.