Development and Validation of a Machine Learning-Based Nomogram for Predicting Severe Pneumocystis jirovecii Pneumonia in Kidney Transplant Recipients

Introduction

As the first-line treatment for end-stage renal disease, kidney transplantation (KT) greatly improves the survival and quality of life of patients. However, to avoid graft rejection, KTRs must take lifelong immunosuppressive therapy (usually including calcineurin inhibitors, antimetabolites, and glucocorticoids), which simultaneously induce profound and persistent impairment of cellular and humoral immunity [1,2]. This iatrogenic immunosuppressive state makes KTRs prone to many opportunistic infections. Among them, PJP has become one of the most common and fatal post-transplant diseases after KT [3]. PJP has an extremely slow onset, and fever, dry cough, progressive dyspnea, and hypoxemia are the main clinical manifestations. The disease develops rapidly, and critically ill patients are admitted to the ICU for treatment. The mortality rate can reach 50%. Although sulfonamide-based prophylaxis is widely used in clinical practice, some patients still developed breakthrough infections after discontinuation of prophylaxis or because of drug intolerance, which led to a poor prognosis [4–6]. Currently, the lack of reliable and easily available biomarkers for early risk stratification and the heterogeneity of clinical courses of PJP further aggravate the difficulties in treating PJP [7].

In the past few years, researchers have focused on PJP risk factors and constructed machine learning models to predict the risk of PJP onset [8,9]. Additionally, although there have been a few studies that attempted to establish models to predict PJP severity, most of them used only 2 variable selection methods or even only 2 algorithms (eg, logistic regression). These methods have poor robustness in variable selection and inadequate control of multicollinearity [10]. Therefore, based on the clinical data and laboratory data from KTRs with PJP, we collected the baseline demographic data and hematological indicators at hospital admission and constructed a multivariate predictive model by using 3 variable screening methods (Boruta algorithm, RFE, and LASSO). Finally, we visualized the model as a nomogram with the objectives of early identifying high-risk PJP patients who are likely to progress to severe PJP or death, quantifying the risk of adverse events in PJP patients after KT, rationally allocating health care resources, and ultimately improving the prognosis of patients.

Material and Methods

STUDY DESIGN AND PARTICIPANTS:

This was a single-center retrospective cohort study encompassing 169 KTRs diagnosed with PJP that occurred between June 2019 and September 2025 in Shandong Provincial Qianfoshan Hospital. All included patients received regular follow-up at the study center after KT. Moreover, all KTRs received a standardized PJP prophylaxis strategy at the study center (Shandong Provincial Qianfoshan Hospital): low-dose oral trimethoprim-sulfamethoxazole (TMP-SMX) routinely used within 12 months after KT. Prophylactic doses adhered to clinical practice guidelines (typically TMP 80 mg/SMX 400 mg once daily or 3 times weekly), with dosage adjustments based on individual renal function to mitigate the risk of adverse events, including myelosuppression, liver dysfunction, or nephrotoxicity [11,12].

Inclusion criteria were defined as meeting all of the following: 1) kidney transplant recipients; 2) meeting the diagnostic criteria for PJP; 3) age ≥18 years; 4) complete and retrievable clinical data during hospitalization. Exclusion criteria were any of the following: 1) PJP onset prior to kidney transplantation; 2) age <18 years; 3) concomitant malignant tumors (eg, lymphoma, solid tumors) or other immunodeficiency diseases (eg, human immunodeficiency virus infection); 4) concomitant transplantation of other organs (heart, liver, lung, pancreas); 5) incomplete clinical data or loss to follow-up.

The study was approved by the Medical Ethics Committee of Shandong Provincial Qianfoshan Hospital (Ethics Approval Number: 2025[S920]). As this was a retrospective study, and patient personal identifying information (such as name, hospital number) was anonymized, this study followed the regulations set by the ethics committee for waiving informed consent for retrospective studies. Therefore, the requirement for informed consent from the study participants was waived. The study was conducted in accordance with the Declaration of Helsinki and the regulations on medical research ethics in China.

CLINICAL DEFINITION:

PJP was defined as satisfying both microbiological confirmation and clinical diagnostic criteria [3,13]. Specifically, microbiological confirmation referred to a positive test result for Pneumocystis jirovecii in respiratory specimens (sputum, bronchoalveolar lavage fluid) or blood samples using metagenomic next-generation sequencing (mNGS) and/or real-time polymerase chain reaction (qPCR) of sputum. Clinical diagnosis was defined as the presence of typical symptoms, including fever, dry cough, chest tightness, or dyspnea, and accompanied by chest computed tomography (CT) features consistent with bilateral diffuse ground-glass opacities or consolidation. Cytomegalovirus (CMV) infection was defined as detection of positive CMV-DNA in peripheral blood and/or urine using qPCR. Similarly, BK virus (BKV) infection was defined as detection of positive BKV-DNA in peripheral blood and/or urine using qPCR. Delayed graft function (DGF) was defined in accordance with the Chinese Scientific Registry of Kidney Transplantation (CSRKT) and the kidney transplantation quality control standards of the National Health Commission: this condition required dialysis within 7 days postoperatively and/or that showed no significant decrease in serum creatinine (at least 10% daily) for 3 consecutive days during the first week after kidney transplantation [14]. Acute graft rejection was defined as either antibody-mediated rejection (ABMR) or T cell-mediated rejection (TCMR) confirmed by renal allograft biopsy and meeting the Banff classification criteria [15,16]. Severe PJP was defined as PJP in patients who required ICU admission or died despite active treatment.

CLINICAL DATA COLLECTION:

All data in this study were extracted from the EMRs of study participants upon the diagnosis of PJP and were completely de-identified before data cleaning and analysis. Data collection included the following aspects: basic demographic characteristics (age, sex, body mass index [BMI], and history of hypertension) and transplant-related clinical characteristics (donor type of kidney, time interval from transplantation to onset of PJP, perioperative induction immunosuppressive regimen, and history of acute rejection). In addition, we systematically collected a comprehensive set of laboratory parameters upon hospital admission: complete blood count (white blood cell count, neutrophil count, lymphocyte count), inflammatory markers (CRP, procalcitonin [PCT]), fungus-specific biomarker ((1→3)-β-D-glucan test [G_test]), key indicators of immune function (CD4+ and CD8+ T cell count), renal function (serum creatinine[Scr]), and PCR detection results for common viruses (CMV/BKV). Finally, to evaluate disease severity and prognosis in a comprehensive manner, we also collected key clinical outcome events: length of hospital stay, ICU admission rate, and cure discharge rate.

DATA PROCESSING:

Initially, a missing value analysis was conducted. Missing data were predominantly observed in laboratory parameters. Features with a missing data rate exceeding 20% were excluded, while those with less than 20% missing data were imputed using the multiple imputation method. Subsequently, data normalization was performed to eliminate discrepancies in data scales, thereby enhancing the suitability of the data for statistical analysis and modeling. The entire cohort was randomly partitioned into a training set and a testing set at a ratio of 7: 3. The training set was used for developing machine learning models, while the testing set served to validate the performance of the established models. All data processing procedures were conducted in the R statistical environment using the “mice” and “caret” packages.

VARIABLE SCREENING, MODEL DEVELOPMENT:

To get the key variables most relevant to the risk of severe PJP in training dataset as accurately as possible, this study applied 3 machine learning-derived variable selection methods, including Boruta algorithm, RFE, and LASSO regression. Specifically, the Boruta algorithm assessed the overall relevance of all predictive variables with respect to the outcome by inducing shadow features and comparing importance values. It randomly permuted original attributes to generate shadow attributes. If an original attribute is unimportant, its importance value should be similar to importance values of shadow attributes. Otherwise, its importance value should be higher than most importance values of shadow attributes. If an original attribute’s importance value was higher than most importance values of shadow attributes, it was considered important, otherwise, it was considered unimportant [17]; RFE selected the optimal feature subset incrementally by recursively building a model on the training dataset and then eliminating the least important training dataset using cross-validation performance as the criterion [18]; LASSO regression introduced an L1 regularization to shrink regression coefficients and the algorithm converged regression coefficients of irrelevant variables to zero so that it could also achieve variable selection and dimensionality reduction [19,20]. To minimize the bias from a single selection method, the variables selected by 3 approaches described above were used, and the intersection of their results was used as candidate variables for subsequent multivariate analysis when severe PJP occurrence was used as the outcome variable, and these candidate variables were included in multivariate logistic regression analysis and then the final predictive model was constructed. Finally, to improve the clinical applicability and visual interpretation of the result for model, the final regression model was visualized as a nomogram to intuitively show the relevance of each predictive variable to the risk of severe PJP and then provide a visual tool for risk evaluation of severe PJP in clinical practice.

MODEL VALIDATION AND CLINICAL UTILITY EVALUATION:

To assess the accuracy and clinical applicability of the developed predictive model in general, the variance inflation factor (VIF) was used to diagnose the multicollinearity of variables included in the model. Multicollinearity is suspected when a VIF is above 5, and variable optimizations are then performed to stabilize the model. For assessment of model performance, the ROC curve was obtained and the AUC was calculated to evaluate the discriminative ability of the model for outcome events. In addition, the Hosmer-Lemeshow goodness-of-fit test and calibration curves were used in conjunction to evaluate the agreement between the model’s predicted probabilities and observed events and thus evaluate its calibration ability. To further visually evaluate its clinical application, the clinical net benefit of the model at different threshold probabilities was calculated, and the decision curve analysis (DCA) curve and clinical impact curve (CIC) were obtained. Finally, the optimal risk threshold was determined according to the nomogram scores, and the patients were divided into different risk groups according to the occurrence of endpoint events in them to evaluate the risk stratification ability of the model.

STATISTICAL ANALYSIS:

Continuous variables were first tested for normality using the Shapiro-Wilk test. Variables with a normal distribution were expressed as mean±standard deviation χ̄±s), with intergroup comparisons performed with the independent-samples t test, those without a normal distribution were expressed as median (interquartile range) [M (IQR)], with intergroup comparisons conducted with the Mann-Whitney U test. Categorical data were expressed as counts (percentages) [n (%)], and intergroup comparisons were performed using the Pearson χ² test or Fisher’s exact test, depending on theoretical frequencies and total sample size. All statistical analyses in this study were performed with SPSS (version 27.0.1) and R software (version 4.5.0). The significance level α was set at 0.05, and a P value <0.05 was considered statistically significant.

Results

BASELINE PATIENT CHARACTERISTICS:

A total of 169 eligible KTRs with post-transplant PJP were enrolled in this study and randomly divided into a training set (n=120) and a testing set (n=49), as shown in Figure 1. The baseline characteristics of the 2 groups are presented in Table 1. In total, males accounted for 81.66% and the prevalence of hypertension was 90.53%. The median time from KT to PJP attack was 10 months (IQR: 5–17). Laboratory investigations revealed a significant inflammatory response (median CRP level was 40.7 mg/L, IQR: 20.05–66.75) and a great decrease in CD4+ T cell count (147 cells/μl, IQR: 94–282.5). Moreover, 33 patients (19.53%) were admitted to ICU because of severe illness. Overall, except for the higher CRP and PCT levels, the 2 groups were well matched in terms of demographic characteristics, clinical indicators, and immune status.

VARIABLE SCREENING RESULTS:

Because the Boruta algorithm calculated the importance of variables by comparing the predictive variables with their shadow variables, the results showed that there were 6 variables whose mean importance values were significantly higher than those of other variables: PCT, G_test, CRP, N (neutrophil count), Time KT to PJP (time from kidney transplantation to PJP onset), and WBC (white blood cell count) (Figure 2A). In addition, when RFE analysis was performed (Figure 2B), we found that with the increase in variable count from 1 to 5, the cross-validated RMSE curve rapidly decreased to its minimum. With further increases in variable count, there was no significant improvement in cross-validated RMSE, and the RMSE fluctuated slightly. This clearly indicated that the model with 5 variables had the best predictive performance. Moreover, the coefficients of irrelevant variables were compressed to zero by L1 regularization in LASSO regression (Figure 2C, 2D), and the results of the 3 methods were integrated in Figure 2E. The results of LASSO regression further verified that PCT, G_test, CRP, and Time KT to PJP were the established variables significantly associated with PJP outcomes. Multicollinearity analysis was then performed on these 4 established variables. The results showed that there was no multicollinearity between these established variables, and these variables were included in the predictive model. The coefficient path plot showed that as the penalty parameter λ increased, most variables were eliminated, and only a few had non-zero coefficients. PCT, G_test, CRP, and Time KT to PJP were retained throughout the regularization path (Figure 2C).

PREDICTIVE MODEL DEVELOPMENT AND PERFORMANCE EVALUATION:

Through the multivariate logistic regression analysis, we found the independent predictors of severe PJP in KTRs and displayed the results in forest plot (Figure 3). The OR for Time KT to PJP was 0.87 (95% confidence interval [95% CI]: 0.74–0.96, P=0.041). In other words, a longer post-transplantation time was associated with a lower possibility of severe PJP. Although the P value for PCT in the logistic regression model was marginally significant (P=0.087), the 95% CI (1.02–9.13) indicated a potential risk association. It is noteworthy that the CIs were estimated using the profile likelihood method, which accounts for asymmetry in the likelihood function better than the Wald approximation used for P value calculation. Given that PCT was rigorously screened using the Boruta, RFE, and LASSO algorithms, and the profile likelihood CI suggests a positive association, we consider PCT to be a clinically and potentially prognostically relevant predictor in the final model. The ORs for G_test and CRP were 2.94 (95% CI: 1.42–7.75, P=0.014) and 3.08 (95% CI: 1.12–9.85, P=0.042), respectively, and both showed statistically significant association with increased risk of severe PJP.

Based on the regression coefficients of these variables, we developed a nomogram model for the individualized prediction of severe PJP risk in KTRs (Figure 4). This model allows for the intuitive calculation of the corresponding probability of disease occurrence by summing the scores of each variable and mapping them to the total score axis. The discriminative ability of the model was evaluated using the ROC curve. In the training set, the AUC was 0.935 (Figure 5), indicating excellent discriminative ability of the model. In the testing set, the AUC was 0.886 (Figure 5), which verified the good generalization performance of the model. The calibration of the model was validated using the Hosmer-Lemeshow test and bootstrap calibration curves. For the training set, the P value of the Hosmer-Lemeshow test was 0.652, suggesting that the model had a good fit, with no significant bias. Additionally, calibration curves based on 1000 bootstrap samples demonstrated that the calibration effects of the training and testing groups were close to the ideal scenario, confirming the reliable predictive performance of the model (Figure 6A, 6B). To assess the clinical utility of the model, we performed DCA and plotted CIC. The results showed that within a certain range of high-risk thresholds, the nomogram provided more net benefits compared with the “treat all” and “treat none” strategies (Figure 7A, 7B). The CIC further confirmed that at a threshold of 0.177, the model could effectively identify high-risk individuals, and the number of predicted high-risk individuals matched well with the number of actual events (Figure 8A, 8B).

RISK STRATIFICATION AND SURVIVAL ANALYSIS BASED ON THE NOMOGRAM MODEL:

Based on the nomogram scoring model constructed in this study, we first calculated the total score of all KTRs with PJP. Through ROC curve analysis and Youden index calculation, −0.442 was determined as the optimal cutoff value, which allowed for the stratification of patients into the low-risk group (total score ≤−0.442) and high-risk group (total score >−0.442). Subsequently, the Kaplan-Meier method was used to plot the survival curves of the 2 groups to enable intuitive comparison of prognostic differences (Figure 9). The results demonstrated that as hospital stay duration increased, the survival probability of patients in the high-risk group showed a more significant downward trend. In contrast, the survival probability of patients in the low-risk group remained relatively stable throughout the observation period, and it was significantly higher than that in the high-risk group at all time points. Furthermore, the log-rank test was performed to quantify the survival difference between the 2 groups, and the results showed that the survival curves of the 2 groups were statistically significantly distinct (P<0.0001). These results suggest that the risk stratification developed from the nomogram model could also discriminate the prognostic heterogeneity among KTRs with PJP. As the method of calculation was simple and the cutoff value was clear, this model has important clinical predictive value for in-hospital survival prognosis of KTRs with PJP and could be a reference for clinicians to formulate individualized intervention strategies.

Discussion

Risk assessment after KT for severe conversion from mild PJP is an important clinical issue in daily clinical practice. In particular, when patients are iatrogenically immunosuppressed, they often exhibit subtle early manifestations that frequently result in severe pneumonia. Clinically valid and reliable biomarkers that can be easily applied to predict which patients are at high risk are still lacking [7,21,22]. This clinical gap highlights the need and clinical value of constructing a predictive model for severe PJP, which was the primary goal of this study.

We retrospectively analyzed clinical information of 169 patients with laboratory-confirmed PJP after KT in a single center. This study applied the Boruta algorithm, RFE, and LASSO regression in conjunction to screen core variables 3 times. We constructed a nomogram predictive model and internally validated it. The model was composed of 4 predictors: Time KT to PJP, PCT, G_test, and CRP. The model demonstrated excellent discriminative performance, as indicated by area under the AUC values of 0.935 in the training set and 0.886 in the testing set, respectively, suggesting strong discriminatory performance. In addition, it also demonstrated good calibration and significant clinical net benefit, thereby providing a practical tool for achieving individualized management of high-risk patients.

To further enhance the robustness and clinical interpretability of variable screening, we incorporated 3 variable screening algorithms combined with cross-validation into this study design. The 4 final variables reflect the pathophysiological process of PJP from 3 perspectives. Specifically, PCT and CRP reflect the intensity of the host’s systemic inflammatory response to infection [23,24], the G_test directly reflects the fungal load [25,26], and Time KT to PJP acts as a surrogate for the cumulative exposure to immunosuppressive therapy [27]. Compared with most previous studies that merely used 1 variable screening method, our strategy effectively reduced selection bias caused by the limitations inherent to each method. The core variable set obtained by our strategy is not only small but also free from multicollinearity and fully follows the “Parsimonious Model” principle proposed by Steyerberg [20,28].

Detailed analysis of the direction and strength of effects of each predictor showed agreement with previous pathophysiological knowledge and recent breakthroughs in research, as well as several disagreements needing further study. Regarding classical inflammatory biomarkers, numerous previous studies have demonstrated a positive association between elevated CRP and PCT levels and the severity of PJP [23,24,29]. The amounts of these 2 markers are directly correlated with the strength of the host’s inflammatory response to Pneumocystis jirovecii infection. It has been further confirmed in our study that they play a key role in the prognosis of severe PJP. Although PCT did not achieve the traditionally accepted level of statistical significance (P>0.05) in this study, there was still a very strong directional trend of association with severe risk. The conclusion that elevation of PCT is not significant in isolated PJP infection in some previous literature may differ from our study [30,31]. We speculate that this may be because the studied population was composed of KTRs receiving immunosuppressive therapy, who are more likely to be affected by other infections or to present with uncontrolled inflammatory response, a feature that makes the prognostic significance of PCT in this studied population more apparent [32,33]. This small difference makes it evident that the interpretation of biomarkers in an immunosuppressed host requires extreme caution and must be accompanied by the clinical situation to allow for a comprehensive judgment [34].

The G_test was one of the most stable predictors in our study. From the pathophysiologic standpoint, the cell wall of Pneumocystis jirovecii is rich in (1→3)-β-D-glucan (BG). When the fungus proliferates massively in alveoli and breaks the lung tissue barrier, BG is released into the circulation, causing increased serum BG levels [35]. Many studies have confirmed the clinical application of the G_test in the diagnosis of PJP. However, studies on the association between the G_test and disease severity or prognosis are scarce and the conclusions are controversial [25,36]. Some studies have reported that BG levels can reflect the fungal load in lungs indirectly, and patients with PJP who have increased BG levels have a significantly higher rate of poor outcomes [36]. However, other studies have reported that PJP patients do not show a trend of decreasing G_test results after antifungal treatment, and the trend cannot be used to evaluate the prognosis of invasive fungal diseases [37,38]. In contrast, in the present study, we found that there was a clear and direct relationship between G_test levels and the risk of severe PJP. This relationship not only reflects the positive correlation between fungal load and the degree of lung injury, but also provides direct clinical evidence that the G_test can be used to stratify the severity of PJP [39,40]. Some studies have reported that the false-positive rate of the G_test in KTRs is about 15%, which is caused by factors such as existing bacterial infection, intestinal mucosal barrier damage, or the use of medical devices containing glucan [25]. However, our study effectively excluded the interference of false-positives through use of rigorous inclusion criteria, satisfying both the standard of microbiological confirmation and clinical diagnostic criteria. Therefore, the results are reliable [25,38]. In addition, the G_test is rapid and easy to operate. It is particularly suitable for quick clinical risk stratification.

Importantly, there was a significant negative correlation between Time KT to PJP and the risk of severe disease, which basically reflects the “immunosuppression intensity-time dependence” effect. Patients who developed PJP in the early post-transplant period were more likely to present with severe disease due to higher immunosuppression intensity or incomplete immune reconstitution [27,39]. This was also verified visually by risk-stratified Kaplan-Meier survival analysis. In addition, our study quantified this time-dependent association in which for every additional month after KT, the risk of severe PJP will decrease by 13%; OR=0.87 [27], and clearly shows the limitation of the current traditional “unified 6–12 month PJP prophylaxis regimen”. Specifically, for patients with slow immune reconstitution (eg, persistent CD4+ T-cell count <150 cells/μL), discontinuing PJP prophylaxis at 12 months may be premature and increase the risk of severe PJP. Conversely, low-risk patients with a CD4+ T-cell count >200 cells/μL at >12 months after transplant can safely discontinue prophylactic medications. Therefore, in contrast to rigid time-based thresholds, a risk-stratified individualized management strategy is essential for achieving precision in prevention of PJP [3,4,6,40].

This study has several limitations. First, as a single-center retrospective study, despite the use of multivariate analysis to control for confounding factors as much as possible, there remains a potential risk of selection bias. For instance, our center routinely uses a 12-month post-transplant TMP-SMX PJP prophylaxis regimen, while some centers extend the PJP prophylaxis period to 24 months or use nebulized pentamidine as an alternative prophylactic agent. Such variations in PJP prophylaxis strategies may affect the predictive value of Time KT to PJP. Second, the sample size of this study was relatively small. Third, despite random partitioning, statistically significant imbalances in baseline CRP and PCT levels were observed between the training and testing sets, with the training cohort exhibiting higher inflammatory markers. Such distributional shifts, while occasionally unavoidable in limited sample sizes, can introduce estimation bias during internal validation. However, the nomogram maintained robust discriminative power (AUC 0.886) in the testing set despite these lower baseline inflammatory indices, suggesting a degree of model resilience. Nevertheless, this baseline discrepancy further underscores the necessity for future external validation to rigorously test the model’s stability across populations with varying spectra of disease severity.

Conclusions

Based on the clinical data of patients with PJP after kidney transplantation, this study employed 3 machine learning algorithms to screen variables to develop and validate a prognostic prediction model. This model exhibited excellent discriminative ability and calibration performance, which can effectively stratify the risk of severe PJP and demonstrates potential in individualized risk management. Future multicenter studies are required to further validate its generalizability.