Development of Biomarkers and Prognosis Model of Mortality Risk in Patients with COVID-19

Introduction

It has been three years since the World Health Organization announced the definition of Coronavirus Disease 2019 (COVID-19) on 30 January 2020. The main clinical manifestations of patients with COVID-19 include fever, cough, pharyngalgia, etc. Some of these patients may develop into acute respiratory distress syndrome (ARDS) and multiple organ dysfunction syndrome (MODS), especially the elders or those patients with underlying diseases which can lead to fatal outcomes. As of 30 April 2023, the COVID-19 pandemic has resulted in over 6.9 million deaths worldwide. In December 2022, the Chinese government issued the new policy for the management of COVID-19. This policy led to a sudden increase in the number of newly infected patients, posing significant challenges to clinical work and experience.¹ The statistics information about COVID-19 can be downloaded at: https://covid19.who.int/. Low lymphocyte (LYM) and white-blood-cell (WBC), high CRP and ferritin were effective in the diagnosis of COVID-19.^2,3 Indicators were found that high ESR, international-normalized-ratio (INR), prothrombin-time (PT), CRP, D-dimer, ferritin and red-cell-distribution width (RDW) values, respectively, were the most effective predictive mortality risk biomarkers of COVID-19.⁴ The accurate and reliable estimation of Oxidant/Antioxidant levels in COVID-19 patients, utilizing biomarkers such as LYM, ferritin, D-dimer, WBC, and CRP, can facilitate the diagnosis and prognosis.⁵

Most studies on COVID-19 were conducted during the initial outbreak, in which virus frequently mutated in the widespread, and genetic recombination occurs when different subtypes infect the human body. Gene mutation or recombination can affect virus’s biological characteristics,^6,7 influencing the whole pathophysiological process. Patients with severe pneumonia caused by COVID-19 infection can be challenging to reverse once the case develops into ARDS, so early identifying cases who may progress to severe are crucial. This study aims to statistically analyze the clinical data of infected cases during the special period, summarize clinical treatment experiences, and identify factors that may affect the prognosis, understand the characteristics of current viral infection process and provide evidence-based guidance for later treatment and early key risk factors.

In recent years, we have observed a significant increase in the utilization of artificial intelligence technologies in various fields, especially in the field of medicine.⁸ Artificial intelligence technologies have started to be used frequently in the diagnosis, prognosis and treatment processes of diseases. The most important reason for this is the machine learning (ML) algorithms have the power to reveal hidden relationship structures between features.^9–11 When the literature is reviewed, there are many attempts using the ML methods to predict the diagnosis and mortality of COVID-19.^12–15Considering that many clinical data have the characteristics of non-linearity, complexity and heterogeneity, this study tries to utilize machine learning to identify essential prognostic factors. Ultimately, the result may remind physicians which patient should be pay attention to guide clinical treatment.

Materials and Methods

Participants

This study collected clinical data of all patients diagnosed with COVID-19 and hospitalized in Xinhua Hospital Affiliated to Shanghai Jiao Tong University from December 15, 2022 to January 15, 2023. The inclusion criteria were (1) the patient had a positive nucleic acid result or a positive antigen test at the first visit (2) hospitalization (3) sufficient clinical data. The exclusion criteria were (1) patients who required hospitalization for reasons other than the respiratory system (2) failure to complete treatment or patients who were discharged voluntarily (3) under 18-year-old. Finally, a total of 1379 cases were included. Flowchart is shown in Figure 1.

Figure 1 (a) Flowchart showing the exclusion and enrolment of COVID- 19 patients; (b) Flowchart of model development and performance by machine learning method. XH hospital: Xinhua Hospital Affiliated to Shanghai Jiao Tong University. Abbreviations: Labs, Laboratory tests; SHAP, SHapley Additive exPlanations.

Data Collection

The study collected clinical data under the guidance of a multidisciplinary team of experienced clinicians and informaticists. Data collection involved 45 characteristics, including age, sex, past medical history (Charlson comorbidity index, congestive heart failure, chronic pulmonary disease, rheumatic disease, renal disease, liver disease, diabetes), length of stay, prognosis and laboratory tests including C-reactive protein (CRP), white blood cell (WBC), neutrophil number (ne_num), lymphocyte number (ly_num), monocyte number (mo_num), hemoglobin (Hgb), red blood cell distribution width (rdw), thrombocyte (plt), thrombocyte distribution width (pdw), partial pressure of carbon dioxide (pCO2), apolipoprotein E (APOE), low-density lipoprotein (LDLC), high-density lipoprotein (HDLC), total cholesterol (TCH), triglyceride (TG), troponin I (TNI), N-terminal pro-brain natriuretic peptide (NT-proBNP), international normalized ratio (INR), activated partial thromboplastin time (APTT), fibrinogen (Fb), D-dimer (DD), total protein (tp), albumin (alb), total bilirubin (tbil), alanine transaminase (alt), aspartate transaminase (ast), creatinine (crea), interleukin-8 (IL-8), Interleukin-1β (IL-1β), interleukin-6 (IL-6) tumor necrosis factor (TNF), interleukin-2 receptor (IL-2R), growth stimulation expressed gene (ST2). Biochemical tests were analyzed by the Roche Cobas E702 Fully Automatic Biochemical Analyzer (Roche, Germany). Sysmex XS-1000i Hematology System (Sysmex Corporation, Kobe, Japan) was used to carry out cell blood count. Researchers independently entered and double-checked the data. The first test results at admission were selected as eigenvalues for many examinations to ensure the accuracy of prognosis prediction.

Data Preprocessing

The K-Nearest Neighbor method (R version 4.2.1, package: DMwR2) was employed to fill the missing data.¹⁶ The reference variable “age” was subsequently employed. Age values exceeding three times the standard deviation of the mean squared were considered outliers and excluded from the analysis. In total, 1369 patients were incorporated into the dataset.

Model Development

We divided the patient ‘s final outcome at discharge into two groups that are death or survival. The whole data set was randomly divided into training set (80%) and test set (20%). Because the training set is divided by random values, in order to avoid contingency, we repeated the modeling of each algorithm 10 times, and took the average and variance of the 10 results as the final result. We have built a total of 15 models, their brief introduction is as follows: (1) Light Gradient Boosting Machine (LightGBM) is a gradient boosting framework based on decision trees, utilizing a histogram-based algorithm and a leaf-wise tree growth strategy to reduce memory consumption and enhance computational efficiency. Its features include efficient handling of categorical features and support for large-scale data training through distributed learning. (2) Ridge Regression is a technique that incorporates an L2 norm regularization term into the linear regression loss function to address multicollinearity problems and control model complexity. This method balances the ability to fit data with the need to keep model coefficients small by adjusting the regularization parameter. (3) Logistic Regression is a widely used linear model for binary classification problems, employing a sigmoid function to map the output of the linear model to the (0,1) interval, representing probabilities. The parameters of the logistic regression model are estimated using the Maximum Likelihood Estimation (MLE) method. (4) Adaptive Boosting (AdaBoost) algorithm iteratively trains a series of weak learners, increasing the weight of samples that were misclassified in the previous round, thereby focusing subsequent learners on these difficult-to-classify samples to improve the overall model performance. (5) CatBoost is a gradient boosting decision tree algorithm optimized for categorical feature handling. It reduces overfitting and enhances model accuracy through unique ordered boosting techniques and efficient encoding of categorical features. (6) Linear Discriminant Analysis is a classic linear classification method that finds the best linear classification plane by maximizing the criterion of between-class distance and minimizing the within-class variance. It also serves as a dimensionality reduction technique, representing data through the most discriminative linear combination.

(7) Decision Trees are a non-parametric supervised learning method used for classification and regression tasks. They recursively partition the dataset to build a tree structure, where each node represents a decision rule based on an attribute, and leaf nodes represent decision outcomes. (8) The Extremely Randomized Trees (Extra-Trees) algorithm introduces randomness in the training process of decision trees by choosing random features and split points to grow trees, thereby increasing model diversity and reducing the risk of overfitting. (9) Random Forest is an ensemble learning method that improves prediction accuracy and stability by constructing multiple decision trees and aggregating their predictions. It introduces bootstrap sampling of the samples and random selection of features to increase the independence of decision trees. (10) Extreme Gradient Boosting (XGBoost) is an efficient gradient boosting algorithm that prevents overfitting through optimized computational resource usage and advanced regularization techniques, such as L1 and L2 regularization. It offers advanced features for parallel processing, tree pruning, and automatic handling of missing values. (11) Multi-Layer Perceptron (MLP) is a type of feedforward artificial neural network, modeling complex nonlinear relationships through multiple layers of nodes (or neurons) and nonlinear activation functions. MLPs consist of at least three layers: an input layer, one or more hidden layers, and an output layer. (12) The K-Nearest Neighbor algorithm is an instance-based learning method that classifies or regresses by measuring distances between different feature values. For a given query point, the algorithm identifies the K nearest neighbors in the training data and predicts the outcome based on these neighbors’ information. (13) Gradient Boosting is an ensemble technique that incrementally adds predictors and optimizes them to reduce the overall model’s loss. It models the residuals at each step to gradually enhance the model’s predictive capability. (14) The naive Bayes classifier is based on Bayes’ theorem and assumes independence among features. Despite this assumption often being unrealistic in practice, naive Bayes can still deliver strong performance in various scenarios, especially in text classification and spam detection. (15) Support Vector Machine (SVM) is a powerful classification technique that finds the optimal hyperplane for separating different categories. SVMs with Radial Basis Function (RBF) kernels address classification problems in nonlinear feature spaces by mapping the input space to a higher-dimensional feature space to find the optimal separating hyperplane.^15,17

Model Evaluation

Accurate value is used to evaluate the accuracy of the model in the prediction task, which represents the ratio of the number of samples correctly predicted by the model to the total number of samples. For example, if a decision tree model correctly predicts 80 samples out of 100 samples, the accuracy is 80%. ROC curve is a visual index to evaluate the performance of the model. Its main function is to judge the recognition ability of a classifier to samples at each threshold. AUC is the area under the ROC curve. The higher the AUC value, the better the classification effect of the machine learning algorithm.

In machine learning and data mining tasks, by ranking the importance of features, we can determine which features play a key role in the prediction of target variables, thereby simplifying the model and reducing computational complexity, so that we can better understand potential risk factors and incorporate them into the decision-making process.

SHapley Additive exPlanations (SHAP) is a method used to interpret the prediction results of machine learning models. It bases on the Shapley value in game theory, each sample is treated as a player in the game, the combination of each eigenvalues is regarded as the cooperation strategy of the participants, and the marginal contribution of features can be measured. We plotted the SHAP value of each feature of each sample to visually show how these variables affect the prognosis of COVID-19 patients.

The confusion matrix allows us to further evaluate the performance of the model. TP (True Positive) is the number of positive examples that the model predicts correctly, FN (False Negative) is the number of positive cases predicted by the model is wrong, FP (False Positive) is the number of negative cases predicted by the model is wrong, TN (True Negative) is the number of negative cases predicted by the model is correct. On this basis, F1 score, sensitivity and specificity are calculated to comprehensively evaluate the model, and the calculation formula is shown:

(1) accuracy = (TP + TN)/(TP + FP + FN + TN), (2) sensitivity/R = TPR = TP/(TP + FN), (3) specificity = TN/(FP + TN), (4) precision/P = TP/(TP + FP), (5) F1 score = 2PR/(P + R).

Statistical Analysis

Continuous variables are presented as median (interquartile ranges [IQR]), while categorical variables as counts and percentages (%). Python version 3.7.16 was used to program. Train_test_split function from sklearn.model_selection module was used to split the data proportionally into test and validation sets. Accuracy_score, roc_auc_score, f1_score, and confusion_matrix function from sklearn.metrics module were used to generate accuracy value, AUC, f1 score, precision, specificity, and sensitivity, respectively.

Results

Characteristics of Participants

The study collected the data of patients who were diagnosed as COVID-19 infection and hospitalized in Xinhua Hospital affiliated to Shanghai Jiaotong University from December 15, 2022 to January 15, 2023. The patients with nucleic acid negative at admission were excluded. Finally, 1379 cases were included in the study. The results showed that the median length of stay was 10 days, with the longest length stay lasting 58 days, and the mortality rate was 11% (Table 1). We additionally extracted data on deceased patients into a separate table (Table 2). It showed the proportion of patients with heart disease increased significantly, and the male-to-female ratio remained basically unchanged.

Table 1 Baseline Characteristics of All Individuals

Table 2 Baseline Characteristics of the Non-Survivors

Model Performance

Fifteen different machine learning algorithms were employed to construct prediction models, and the dataset was split into testing and training sets to evaluate the model. All AUROC for models generated by 15 machine learning algorithms exceeded 90%, and the accuracy of 10 of them also surpassed 90% (Table 3, Figure 2). It shows that these 15 algorithms can accurately predict the prognosis of COVID-19’s patients. The Light Gradient Boosting model emerged as the optimal choice, with accuracy of 0.928 ± 0.0006 and an AUROC of 0.976 ± 0.0028. In addition, naive Bayes (0.946 ± 0.0043) and Extreme Gradient Boosting (0.918 ± 0.0013) also have excellent prediction ability. In order to better judge the performance of the model, we use the confusion matrix to evaluate the model. The results are shown in Figure 3 and Table 4. Boosting, a pivotal category of ensemble learning algorithms, operates on the principle of sequentially constructing a series of “weak learners.” These learners are then aggregated to develop a robust final model. The process involves the successive generation of weak learners, each aiming to address the residuals of the cumulative model by aligning with the negative gradient of its loss function. This approach ensures that the introduction of each new learner reduces the overall model’s loss, enhancing predictive accuracy. XGBoost builds upon the foundational concept of boosting, incorporating several enhancements to improve its efficiency and accuracy. These include the use of second-order derivatives for a more precise loss function, regularization terms to prevent overfitting of trees, and block storage structures that facilitate parallel computations. Similarly, LightGBM’s primary advantage lies in its innovative modifications to the training algorithm, which significantly expedite the process and often lead to the development of more effective models. Notably, the XGBoost incorporate inherent algorithms for variable ranking, whereas both LGBM and the naive Bayes model employ a feature_importance function to determine variable significance.

Table 3 Predictive Performance of Models

Table 4 Confusion Matrix Evaluation of Three Models

Figure 2 The receiver operating characteristic curves of 15 different models. X-axis represents false positive rate, y-axis represents true positive rate.

Figure 3 Confusion matrix of the three models, where true positive is 9,14,4 in naïve Bayes model (a), Light Gradient Boosting model (b) and Extreme Gradient Boosting model (c) respectively.

Analyzing Importance of Features Included in Models

A total of 44 features are included in the construction of machine algorithms. The relative importance rank of all 43 variables is analyzed through the best three machine algorithms, Light Gradient Boosting model, naive Bayes and Extreme Gradient Boosting. In Light Gradient Boosting model, ST2 was most contributing to prognosis (19.3%), the second is IL-8 (17.2%), IL-6 (6.4%), age (6.1%), NT-proBNP (5.1%), IL-2R (5%), TNI (4.6%), congestive heart failure (3.3%), the variables that have little effect on prognosis were tumor, rheumatic disease, diabetes; the result in naïve Bayes model shown ST2 is still most contributing to prognosis (18.4%), after that is IL-8 (17.2%), age (7.8%), IL-2R (6.6%), TNI (5%), NT-proBNP (4.1%), IL-6 (3.8%), d-dimer (3.8%), those with a weight of 0 were tumor, rheumatism disease, chronic pulmonary disease, kidney disease and diabetes; The ranking of feature importance in the Extreme Gradient Boosting model showed that ST2 had the highest weight on prognosis (25.6%), followed by IL-8 (14.3%), age (6.6%), TNI (5.7%), congestive heart failure (4.3%), IL-6 (3.5%), AST (3.1%), NT-proBNP (3%), at the end of the ranking, chronic pulmonary disease, rheumatism and tumor made almost no contribution to prognosis (Figure 4). The results showed that the eigenvalues with high weight in the three algorithms were basically the same, the stable type of the model was excellent, and inflammatory factors and the related indexes of cardiac function had great influence on the prognosis.

Figure 4 Rank of importance to prognosis of COVID-19 patients in Light Gradient Boosting model (a), naive Bayes model (b), Extreme Gradient Boosting model (c), x-axis represent proportion of weight, y-axis represent features, st2 is most important feature accounting for 25.6%, 19.3%, 18.4% in the three models respectively. IL-8 is second one accounting for 14.3%, 17.2%, 17.2% in the three models respectively; SHAP value of each feature of every sample in Light Gradient Boosting model (d), naive Bayes model (e), Boosting model (f), red dot represents the positive effect and the blue dot represents the negative effect. Abbreviations: SHAP, SHapley Additive exPlanations; crp, c-reactive protein; WBC, white blood cell; ne_num, neutrophil number; ly_num, lymphocyte number; mo_num, monocyte number; hgb, hemoglobin; rdw, red blood cell distribution width; plt, thrombocyte; pdw, thrombocyte distribution width; pco2, partial pressure of carbon oxide; apoe, apolipoprotein E; ldlc, low-density lipoprotein; hdlc, high-density lipoprotein; tch, total cholesterol; tg, triglyceride; tni, troponin I; nt-probnp, n-terminal pro-brain natriuretic peptide; INR, international normalized ratio; APTT, activated partial thromboplastin time; fb, fibrinogen; dd, d-dimer; tp, total protein; alb, albumin; tbil, total bilirubin; alt, alanine transaminase; ast, aspartate transaminase; crea: creacreatinine; il-8, interleukin-8; il-1β, interleukin-1β; il-6, interleukin-6; tnf, tumor necrosis factor; il-2r, interleukin-2 receptor; st2, growth stimulation expressed gene.

The honeycomb diagram showed 20 characteristic variables in our models, in which the red dot represents the positive effect and the blue dot represents the negative effect. It can be seen that ST2, IL-8, IL-6, heart disease, TNI, and NT-proBNP have significant effects on the poor prognosis of COVID-19 patients.

Discussion

The study collected the data during a major COVID-19 outbreak in Shanghai, a special period, which has better timeliness than the cases in Wuhan. We improved accuracy and stability by using machine learning algorithms rather than traditional retrospective analysis. The clinical information extracted from this study is easy to obtain, making it suitable for broader applications and highly practical. Most blood test can be easily finished when patients are admitted to hospital even in the emergency department. Physicians can collect all the data and predict patients’ outcomes and intervene the treatment plan upon hospital admission. A limitation of the study is that all cases were from a single-center data cohort (Xinhua Hospital affiliated to Shanghai Jiaotong University), the machine algorithm adopted internal verification, lacked external verification, and had certain limitations and regionality. If conditions permit, we would try to collect multi-center data for external verification in the future.

Key findings of the study included the significance of “ST2”, “IL-8”, “TNI”, “d-dimer”, “age”, “congestive heart failure”, “NT-proBNP”, “IL-6”, “IL-2R” in predicting COVID-19 patients’ prognosis. During a COVID-19 infection, the immune system’s response can sometimes become dysregulated, leading to an excessive release of inflammatory cytokines. This overactive immune response, often referred to as a “cytokine storm”, can result in widespread inflammation in the body and lead to severe respiratory distress and other complications.^18,19 IL-8, IL-6 and IL-2R are interleukins, a kind of cytokines, were identified as significant factors in predicting COVID-19 patient outcomes, which was proved in our study. These molecules exert profound biological effects even at low concentrations and play a pivotal role in modulating immune responses and the growth of immune cells.^20,21 Diane Marie found that the elevated serum concentration level of IL-6 and IL-8 serve as independent predictors of survival in 1484 patients with COVID-19;²² Bo Diao found that COVID-19 patients requiring ICU treatment exhibited lower IL-6 and IL-10 levels;²³ Aiping Ma‘s retrospective study further proved that IL-8 and IL-2R were related to the duration and severity of COVID-19.²⁴ Based on these studies, IL-8, IL-6 and IL-2R are important risk factors which physicians should pay attention to. Our study underscores the importance of monitoring interleukin levels in clinical work, with a particular emphasis on the pronounced significance of IL-8. IL-8, also known as CXC8, can produce chemotaxis to neutrophils by binding to G protein-coupled receptors CXCR1 and CXCR2.²⁵ Researchers observed that neutrophil extracellular traps (NETs) in the interstitium and peribronchiolar epithelium were positively correlated with IL-8 mRNA level in lung tissues of dead patients with COVID-19, significantly higher than non-COVID-19 related ARDS patients.²⁶ Furthermore, abnormal regulation of CXCL8-axis formed by IL-8 (CXCL8) and its receptor CXCR1/CXCR2 was associated with respiratory diseases in other studies. The CXCL8-axis can enhance airway permeability to inflammatory responses by recruiting neutrophils and stimulating airway epithelium²⁷ and also can increase airway hyperresponsiveness.²⁸ We believe that the role played by IL-8 in COVID-19 deserves further exploration. IL-6 has garnered significant attention because it is proved to be associated with the hyperinflammatory response that can occur in severe cases of COVID-19.^29–31 Tocilizumab has been employed in the management of COVID-19 patients by disrupting the interaction between IL-6 to its receptor. However, due to the difficulty in obtaining tocilizumab, we only able to administer it in two cases regrettably without achieving a positive reversal in their fatal outcomes. Further clinical studies are needed to ascertain the appropriate timing and actual effect of tocilizumab.

In this study, the prognostic significance of congestive heart failure in COVID-19 patients surpasses that of other factors. Congestive heart failure, in conjunction with crucial biomarkers such as ST2, NT-proBNP, and TNI, acts as pivotal indicators of cardiac function. These parameters suggest that patients with underlying heart disease or when the virus has caused an impact on the heart are at a heightened risk of experiencing adverse outcomes. These indexes have more important clinical significance compared with other organ functions, which strongly resonates with our experiences during the course of clinical treatment. This perspective is further fortified by a lot of preexisting research,^32,33 serving to underscore its validity and importance. Rather than only focusing on the heart disease history, we believe that physicians should not ignore the special damage caused by the virus itself to the cardiac system. Physicians must try to reduce the risk factors associated with death during the course of the disease. Currently, the precise mechanisms through which the virus triggers or exacerbates myocardial injury remain unclear. Some researchers believe that the collective basic energy consumption caused by virus infection and the hypoxia environment increases the heart burden, which may be related to angiotensin-converting enzyme 2 (ACE2).^34,35 Moreover, electrolyte imbalances and adverse drug reactions may also promote the process of heart failure. Therefore, in clinical treatment, physicians should ensure adequate oxygen supply for patients, use sedatives to reduce energy consumption when necessary, and avoid using negative myodynamia drugs, so as to protect cardiac function and avoid poor prognosis.

Conclusions

The machine model algorithm can accurately predict the prognosis of COVID-19 patients. The CatBoost Classifier model was identified as the optimal choice, with “ST-2”, “IL-8”, “TNI”, “d-dimer”, “Age” and “congestive heart failure” were significant predictors of in-hospital mortality among COVID-19 patients.