Accurate and timely assessment of organ dysfunction remains fundamental to mortality prediction and clinical decision-making in the intensive care unit (ICU). The Sequential Organ Failure Assessment (SOFA) score represents a cornerstone tool in critical care medicine, providing clinicians with a quantitative framework to assess organ dysfunction and predict short-term mortality [1]. Originally developed in 1996 [2] (Table S1 in Multimedia Appendix 1) and updated in 2025 [3] (Table S2 in Multimedia Appendix 1), the SOFA score evaluates 6 organ systems—respiratory, cardiovascular, hepatic, coagulation, renal, and neurological—assigning points based on objective clinical and laboratory parameters. This scoring system has demonstrated remarkable utility in stratifying patient risk and informing clinical decision-making across diverse critical care settings [4].

Beyond its predictive value for mortality, the SOFA score serves as a dynamic monitoring tool that captures the trajectory of organ dysfunction over time. The score has been validated in numerous patient populations, including patients with sepsis, acute respiratory distress syndrome, or trauma and patients after cardiac arrest, demonstrating its versatility across heterogeneous ICU cohorts [5]. Furthermore, the SOFA score has been incorporated into sepsis definitions (Sepsis-3 criteria) and quality improvement initiatives, solidifying its role as an essential metric in critical care practice [6].

Despite its clinical utility, the SOFA score was originally developed through expert consensus, with physiological thresholds and point increments defined heuristically rather than derived from outcome-driven data [2]. Each organ domain is discretized into fixed intervals scored from 0 to 4, implicitly assuming equal contributions across score increments and stability of threshold relevance over time. These assumptions may no longer hold in contemporary ICU populations.

Changes in case mix, baseline comorbidity burden, and advances in organ support therapies have altered both the distribution and prognostic impact of key physiological variables [7,8]. Several studies have highlighted limitations of SOFA in modern settings, including imperfect calibration, reduced discrimination in some subgroups, and challenges in accurately capturing neurological dysfunction and chronic organ impairment [5,9,10]. Moreover, SOFA does not distinguish acute organ failure from chronic baseline dysfunction, potentially conflating prognostic signals in patients with preexisting disease [9].

In response to these limitations, there have been calls to revise or update SOFA using contemporary data. Recent large-scale analyses combining expert input with outcome-based modeling have demonstrated that data-driven recalibration of organ failure thresholds is feasible and clinically meaningful [3]. However, such updates remain infrequent, and the core additive structure of SOFA has remained largely unchanged for decades.

Recent advances in AI and machine learning (ML) have opened new avenues for enhancing traditional clinical scoring systems. Multiple studies have demonstrated that AI algorithms can identify complex, nonlinear relationships within clinical data that conventional linear models may overlook [11,12]. ML models have demonstrated strong predictive performance for clinically relevant outcomes in critical care [13]. These AI-enhanced approaches leverage large datasets to recalibrate the weighting of individual organ dysfunction components, potentially capturing temporal patterns and patient-specific trajectories that static scoring systems cannot [14].

However, many high-performing ML models function as black boxes, relying on complex feature interactions that are difficult to reconcile with clinical reasoning. This lack of transparency has been identified as a major barrier to clinical trust and real-world adoption in high-stakes medical settings [15]. Although post hoc explanation techniques such as Shapley Additive Explanations [16] or Local Interpretable Model-Agnostic Explanations [17] are commonly used to interrogate black box models, they operate outside the model’s core logic and may provide unstable or clinically unintuitive explanations [18]. In contrast, additive point-based models remain preferred by many clinicians because they are transparent, verifiable, and aligned with established clinical intuition. This trade-off between accuracy and interpretability has motivated growing interest in hybrid modeling strategies that integrate ML with clinically meaningful structures.

Recent advances in interpretable ML for health care have led to the development of data-driven frameworks for constructing transparent clinical risk scores that combine statistical learning with clinically meaningful representations. Score generation frameworks such as AutoScore [19] have demonstrated that data can produce new additive integer scores with performance comparable to regression approaches. These methods address the problem of constructing scores de novo. This work addresses a complementary problem: updating the parameters of an existing, widely adopted score without altering its structure.

Complementary approaches have further expanded the space of interpretable scoring models. FasterRisk, for instance, formulates score construction as a constrained optimization problem and uses efficient beam search strategies with rounding to produce a set of near-optimal sparse scoring systems [20]. Other methods, such as SymScore, use symbolic regression with monotonicity constraints to optimally group responses and assign integer weights, yielding interpretable additive score tables; it has been applied to shorten existing clinical questionnaires while approximating their total score [21]. Collectively, these approaches demonstrate a growing methodological toolkit for deriving interpretable, data-driven risk scores that align with clinicians’ preference for simple additive models.

In parallel, hybrid modeling strategies have been proposed to reconcile predictive performance with interpretability by explicitly incorporating clinical structure into ML architectures. One recent example is a tree-structured hybrid model for ICU mortality prediction, in which patient features are grouped into clinically motivated modules corresponding to different organ systems or disease processes. Each module is learned from data, while the overall structure reflects established pathophysiological knowledge of multiorgan dysfunction. The outputs of these modules are subsequently aggregated to generate a final risk estimate, allowing interpretation at the level of organ systems or clinical domains. While such hybrid frameworks improve transparency relative to black box models, they typically rely on multiple interacting submodels and do not yield a single, simple scoring rule that can be readily compared to existing clinical scores [22].

Despite these advances, relatively few efforts have focused on recalibrating established clinical scoring systems that are already deeply embedded in routine practice. Widely used severity and risk scores, such as Acute Physiology and Chronic Health Evaluation (APACHE), Simplified Acute Physiology Score (SAPS), and SOFA, were largely developed decades ago using expert consensus and limited datasets. Although numerous ML models have since demonstrated superior predictive performance, these legacy scores remain in widespread use due to their familiarity, interpretability, and extensive validation across clinical settings.

Given the need to modernize aging scoring systems, large-scale collaborative initiatives have occasionally undertaken comprehensive revisions. A notable example is the recent development of SOFA 2.0, an updated organ dysfunction score informed by contemporary ICU practice and large international datasets [1]. This effort combined expert consensus methods with extensive empirical analysis of more than 1 million ICU admissions to revise variable definitions, introduce additional measurements, and adjust thresholds across organ systems. Importantly, SOFA 2.0 preserved the original additive structure and the familiar 0 to 4 per-organ scale, thereby maintaining continuity with the original score while improving its relevance to modern critical care. However, the development of SOFA 2.0 required substantial expert involvement, bespoke data curation, and a resource-intensive consensus process, limiting its generalizability as a methodological template.

More broadly, there remains a clear methodological gap in the literature: the absence of widely adopted, general, automated approaches for outcome-driven recalibration of legacy clinical scores that preserve an additive structure grounded in the clinical domains of the original score. Traditional statistical recalibration techniques typically focus on adjusting model intercepts or recalibrating predicted probabilities in new populations, without modifying the internal composition of the score itself. While such approaches can improve calibration, they do not offer a principled way to update variable thresholds or point allocations in longstanding scoring systems as clinical practice and patient populations evolve [23]. Consequently, when established scores exhibit signs of performance degradation or miscalibration, clinicians are often faced with a binary choice between continued use of an outdated tool or adoption of an entirely new model that sacrifices interpretability and historical continuity. Even in the case of SOFA 2.0, maintaining comparability with the original score has required explicit mapping strategies and periods of dual scoring, underscoring the challenges of updating entrenched clinical instruments.

The widespread deployment of electronic health records (EHRs) in intensive care settings creates an opportunity to leverage routinely collected clinical data for continuous model refinement. Large-scale, publicly available critical care databases—such as Medical Information Mart for Intensive Care IV (MIMIC-IV) and eICU—now provide multicenter data spanning diverse patient populations and institutional practices, enabling robust development and external validation of recalibrated scoring systems. Unlike static scoring systems, data-driven recalibration frameworks can be systematically updated as patient populations evolve, supporting periodic parameter updating of decision support tools that maintain clinical validity over time. Crucially, preserving the interpretable, additive structure of established scores facilitates integration into existing EHR-based workflows and supports clinician trust in automated risk estimates.

To our knowledge, few existing frameworks provide automated, outcome-driven updating of a widely used legacy-score structure while retaining an interpretable, additive, point-based format. Addressing this gap would enable systematic, data-driven modernization of established clinical scores without abandoning the additive structures and clinical familiarity that underpin their sustained adoption.

In this study, we address this gap by developing a data-driven recalibration of a SOFA-derived score that preserves an additive, point-based structure and full clinical interpretability. We introduce a modeling framework that learns optimal threshold cutpoints and interval-specific weights for each SOFA component directly from outcome data. By imposing explicit structural constraints, the approach ensures that each organ system contributes through transparent, clinically interpretable intervals rather than opaque transformations.

By retaining the organ-domain organization and additive, point-based logic of SOFA, the recalibrated score preserves clinical familiarity while allowing thresholds and weights to be updated from data. Importantly, the recalibrated SOFA-derived score is not proposed as a replacement for the canonical SOFA assessment: it addresses mortality risk estimation while preserving the additive structure of SOFA, but its total is not expressed on the conventional 0 to 24 SOFA scale. Through this work, we aim to demonstrate a practical pathway toward interpretable ML in critical care, one that balances predictive performance, robustness, and clinical credibility.

The main contributions of this study include (1) the development of a constrained optimization framework that updates a SOFA-derived score for contemporary ICU populations while preserving its additive structure and clinical transparency; (2) the incorporation of data-driven discretization and integer-weight optimization techniques that maintain full interpretability of the scoring system; and (3) internal validation on a large multicenter ICU dataset (208 US hospitals), demonstrating enhanced discrimination relative to the SOFA-derived reference score. Table 1 summarizes how this study differs from existing approaches. The aim of this study was to develop and internally validate a data-driven recalibration of a SOFA-derived score that preserves its additive, point-based structure and to test the hypothesis that outcome-driven optimization of its cutpoints and integer weights improves discrimination for ICU mortality relative to a prespecified SOFA-derived reference score, without loss of interpretability.

This retrospective multicenter cohort study used the eICU Collaborative Research Database (eICU-CRD), a publicly available repository of deidentified EHR data from adult ICUs across the United States. The database contains 200,859 ICU stays for 139,367 unique patients admitted between 2014 and 2015 to 335 units at 208 hospitals participating in the Philips eICU telehealth program. Available data include high-resolution bedside monitor vital signs, nurse-validated observations, laboratory measurements, medication and infusion records, procedures, APACHE IV/IVa severity of illness scores and diagnostic codes, care plan documentation, and detailed administrative and outcome information, all stored as time-stamped structured data [24].

The eICU-CRD was selected for its large sample size, multicenter design, and high temporal granularity of physiological and treatment variables, making it well suited for developing and validating ML models in critical care. The database is freely accessible to qualified investigators via PhysioNet after completion of required ethics training and acceptance of a data use agreement.

We included adult patients (aged ≥18 y) with their first ICU admission in eICU-CRD for whom all SOFA-related physiological and laboratory measurements were available during the first 24 hours of ICU admission. After applying inclusion criteria and the preprocessing described in the following section, 14,848 patients were included in the final analysis. The detailed selection process is presented in Figure 1.

All SOFA-related variables were obtained from 3 eICU-CRD tables: apacheApsVar, infusionDrug, and lab. From apacheApsVar, we extracted variables directly corresponding to SOFA components: mechanical ventilation status, Glasgow Coma Scale (GCS) total score, and mean arterial pressure (MAP).

Continuous infusions of the 4 SOFA vasoactive agents were extracted from infusionDrug. As dosing units could not be reliably harmonized across participating hospitals, each agent was encoded as a binary indicator of administration during the first ICU day.

Laboratory measurements (serum creatinine, total bilirubin, partial pressure of arterial oxygen [PaO₂], fraction of inspired oxygen [FiO₂], and platelet count) were obtained from the lab table, restricting records to samples drawn within the first 24 hours after ICU admission. For each patient and analyte, the most extreme (worst) value in that window was retained, following the convention used in APACHE for day 1 physiology. Age, sex, and ICU length of stay were additionally extracted from demographic tables to characterize the study population. The patient-level SOFA dataset was constructed by merging these derived variables using the unique patientUnitStayId identifier.

For the primary analysis, encounters were retained only when all continuous determinations required for the common analytic input set were available: paired PaO_₂ and FiO_₂ values required to calculate PaO_₂/FiO_₂, MAP, bilirubin, creatinine, and platelet count. Missing GCS components were handled according to a prespecified operational rule rather than used as an exclusion criterion: encounters receiving mechanical ventilation were assigned GCS a score of 7, whereas nonventilated encounters were assigned a score of 15. No other imputation was performed in the primary analysis.

As vasoactive-agent dose units could not be reliably harmonized across participating hospitals and urine output was incompletely available, the canonical original SOFA score could not be reconstructed exactly. We therefore prespecified a SOFA-derived reference score. Respiratory, neurological, coagulation, hepatic, and creatinine-based renal components followed the original SOFA thresholds for PaO₂/FiO₂ (mm Hg, with mechanical ventilation incorporated in the respiratory criteria), GCS, platelet count (×10⁹/L), total bilirubin (mg/dL), and serum creatinine (mg/dL), respectively. Urine output was not used for renal scoring. The cardiovascular component used MAP (mm Hg) and binary indicators of vasoactive-agent administration: MAP >70 mm Hg was assigned 0 points; for MAP ≤70 mm Hg, absence of all 4 vasoactive agents was assigned 1 point; dopamine, dobutamine, or epinephrine without norepinephrine was assigned 2 points; norepinephrine alone was assigned 3 points; and norepinephrine together with any of dopamine, dobutamine, or epinephrine was assigned 4 points. Laboratory variables used the most extreme value recorded within the first 24 hours, physiological variables obtained from apacheApsVar used the day 1 value selected according to the APACHE extraction logic described earlier, and vasoactive exposure denoted any administration during the first ICU day. Missing GCS components and the complete-case requirements were handled according to the prespecified rules described earlier. This SOFA-derived score served as the reference comparator in all paired performance analyses.

All data extraction and preprocessing were implemented in Python (version 3.10; Python Software Foundation) using pandas (The pandas development team), NumPy (NumPy Developers), and pyarrow (Apache Software Foundation).

To construct an interpretable mortality prediction model, all continuous predictors were transformed into data-driven categorical intervals. The framework proceeds in 2 sequential stages: supervised discretization of each variable, followed by constrained optimization of integer interval weights.

Baseline characteristics were summarized overall and stratified by ICU survival status and by prespecified ICU length of stay (LOS) strata (2, 7, and 30 d). Continuous variables were reported as median with IQR; categorical variables were reported as counts and percentages. Between-group differences were assessed using the Mann-Whitney U test for continuous variables and the chi-square test for categorical variables. All tests were 2 sided, with P <.05 considered statistically significant.

Model performance was estimated using nested cross-validation grouped by hospital. Hospitals were assigned to 10 outer folds, with approximate stratification by ICU mortality and LOS distribution, ensuring that all encounters from a given hospital remained within the same fold and that no center contributed to both development and evaluation. Within each outer development set, hospitals were further partitioned into internal training and validation subsets, with approximately 25% of development encounters assigned to validation by whole hospitals. Candidate discretization schemes and model hyperparameters were fitted using the internal training data and selected according to validation log loss; the selected configuration was then refitted on the complete outer development set and evaluated once on the corresponding held-out hospitals. Each patient therefore contributed exactly one out-of-fold prediction, and pooled predictions were used to estimate overall performance. All comparator models were evaluated using the identical outer folds and equivalent hospital-grouped tuning procedures. Model selection and convergence criteria were prespecified before the revised analysis: configurations were selected by validation log loss, and runs with validation area under the receiver operating characteristic curve (AUROC) below 0.70 or with more than 3 predictors receiving identical weights across all their intervals were treated as nonconvergent. A single random seed was prespecified and held fixed throughout. Following completion of the hospital-grouped cross-validation, the final published score was refitted using the complete analytic cohort with the same prespecified development procedure. This full-cohort model was used only to derive the final thresholds and integer weights; it was not used to estimate predictive performance.

Discrimination was quantified using the AUROC and the area under the precision-recall curve (AUPRC), with 95% bias-corrected and accelerated bootstrap CIs from 5000 stratified resamples, computed on the pooled out-of-fold predictions. Paired comparisons between the recalibrated SOFA-derived score and the SOFA-derived reference score were performed on identical patients within each fold and summarized as the median and IQR of the paired difference across folds. Calibration was assessed using calibration curves comparing observed and predicted ICU mortality across risk strata, together with calibration slope, calibration-in-the-large, the observed-to-expected ratio (O:E), and the Brier score. The logistic calibration layer was fitted on unweighted development data within each fold, so that predicted probabilities correspond to the natural prevalence of the cohort. All analyses were performed in Python (version 3.10).

To assess the sensitivity of the findings to complete-case selection, we repeated the complete hospital-grouped nested cross-validation procedure after expanding the cohort to encounters missing no more than 1 of the 5 required continuous inputs. Under a rule specified before the sensitivity analysis was performed, a single missing bilirubin, creatinine, platelet, or MAP value was randomly assigned from a predefined range entirely contained within the corresponding zero-point category of the original SOFA: bilirubin, 0.2 to 1.1 mg/dL; creatinine, 0.4 to 1.1 mg/dL; platelets, 150 to 300×10⁹/L; and MAP, 70 to 95 mm Hg. Values were sampled from uniform distributions using a fixed random seed. This stochastic normal-range assignment avoided placing all recovered encounters at an identical arbitrary value while preserving the same prespecified assumption of no dysfunction in the affected SOFA domain. For nonventilated encounters with missing PaO_₂/FiO_₂, a ratio between 400 and 500 mm Hg was assigned under the same rule. Ventilated encounters with missing PaO_₂/FiO_₂ (n=5796) were not imputed because assuming normal oxygenation in the presence of mechanical ventilation was not considered clinically defensible. Encounters missing 2 or more required inputs were excluded. GCS was handled according to the same prespecified operational rule used in the primary analysis. The imputation procedure used no outcome information, and the complete model development and hospital-grouped validation procedure was then repeated in the expanded cohort.

Ethical approval from an institutional review board (IRB) or ethics committee was not required for this study. The eICU-CRD is a publicly available, deidentified critical care database that complies with the Health Insurance Portability and Accountability Act Safe Harbor provision for deidentification. Access to the database was granted through PhysioNet after the corresponding author completed the Collaborative Institutional Training Initiative program in human subjects research and signed a data use agreement. As this study involved secondary analysis of preexisting, fully anonymized data with no possibility of patient reidentification, it was deemed exempt from IRB review and individual informed consent requirements.

A total of 14,848 ICU stays (153 hospitals) were included in the analysis, of which 2529 (17.03%) resulted in death and 12,319 (82.97%) in survival (Figure 1). Of 148,532 encounters with outcome information, 133,684 (90%) were excluded from the primary complete-case analysis because at least one required determination was unavailable. Missingness was concentrated in PaO_₂/FiO_₂ data and bilirubin, which were unavailable in 76.5% (113,595/148,532) and 64% (95,114/148,532) of encounters, respectively, compared with 19.6% (29,164/148,532) for creatinine, 8.8% (13,017/148,532) for platelet count, 1.54% (2290/148,532) for GCS, and 0.3% (481/148,532) for MAP. PaO_₂/FiO_₂ data and bilirubin were the sole missing determination in 27.4% (36,584/133,684) and 11.6% (15,477/133,684) of excluded encounters, respectively (Table S3 in Multimedia Appendix 1). Overall, 52,417 (39.2%) excluded encounters lacked exactly one required determination, whereas 81,267 (60.8%) lacked 2 or more. GCS was imputed according to the prespecified operational rule in 437 (2.9%) included encounters, of which 353 were mechanically ventilated. Urine output was available in 51.5% (76,447/148,532) of encounters with outcome information and was missing in 52% (7727/14,848) of the primary analytic cohort; requiring urine output would therefore have reduced the analytic cohort from 14,848 to 7121 encounters. The renal domain was consequently operationalized using creatinine.

Included and excluded encounters differed systematically in several measures of illness severity and treatment intensity (Table S3 in Multimedia Appendix 1). Mechanical ventilation was used in 64.2% (9531/14,848) of included compared with 20.2% (26,978/133,684) of excluded encounters (standardized mean difference [SMD] 0.995), norepinephrine in 28.6% (4243/14,848) versus 8.3% (11,153/133,684; SMD 0.540), and ICU mortality was 17% (2529/14,848) versus 4.5% (6066/133,684; SMD 0.411). Age was comparable between groups (absolute SMDs <0.03), and sex distributions were similar. Among included encounters, ICU mortality was 25.4% (111/437) among those with imputed GCS components and 16.8% (2,418/14,411) among those with recorded components.

Baseline characteristics, stratified by survival status and ICU LOS categories (≤2, ≤7, ≤30, and >30 d), are summarized in Table 2. Most deaths occurred early: 44.40% of nonsurvivors died within 2 days and 81.53% within 7 days. Among survivors, 27.94% were discharged within 2 days and 76.86% within 7 days, with only 118 (0.96%) remaining beyond 30 days.

Across all LOS strata, nonsurvivors were consistently older than survivors (median age 67‐68 vs 61‐64 y; P<.001) and exhibited more severe organ dysfunction. Respiratory status was worse in nonsurvivors, with higher FiO_₂ requirements and lower PaO_₂ (all P<.001). GCS scores were markedly lower in nonsurvivors (P<.001). Hemodynamic and laboratory markers were also more severely abnormal: nonsurvivors had lower MAP, higher serum creatinine and bilirubin, lower platelet counts, and higher heart rate (all P<.001). Invasive mechanical ventilation and vasoactive support—particularly norepinephrine and epinephrine—were significantly more frequent among nonsurvivors (all P<.001). Sex distribution was similar between survivors and nonsurvivors in the ≤2-day stratum but differed in the ≤7- and ≤30-day strata (P=.04 and P=.03, respectively).

After completion of the hospital-grouped cross-validation, the prespecified model development procedure was applied to the complete analytic cohort (n=14,848) to derive the final scoring instrument. The number of bins and model hyperparameters were selected using the prespecified internal validation procedure, after which the selected configuration was refitted on the full cohort. This final refit was used solely to define the thresholds and integer weights of the published score; all predictive performance estimates reported in this study were derived exclusively from out-of-fold predictions from the hospital-grouped cross-validation and not from the full-cohort fit.

The recalibrated SOFA-derived model learned data-driven discretization thresholds and interval-specific weights for all components. The number of intervals per variable ranged from 2 to 5, reflecting variable-specific granularity selected automatically by the algorithm. The model preserved a strictly additive and transparent structure in which each variable contributes through explicitly defined clinical ranges. Table 4 presents the complete threshold values and interval weights, along with those of the original SOFA-derived score. A direct comparison between the SOFA-derived reference and recalibrated discretization thresholds is provided in Figure S2 in Multimedia Appendix 1.

The learned weights are not on the SOFA scale: the recalibrated SOFA-derived score ranges from 0 to 71 rather than 0 to 24, so its weights represent relative contributions within the fitted prognostic score rather than revised 0 to 4 organ-domain points. Maximum contributions were largest for the cardiovascular and respiratory predictor groups and smaller for the coagulation, hepatic, and renal domains. This partly reflects model structure because MAP, mechanical ventilation, and the 4 vasoactive agents enter as separate predictors rather than being combined within composite SOFA domains.

The direction of each contribution was imposed through the prespecified monotonicity constraints, whereas the cutpoint locations and relative magnitudes were estimated from the data. PaO_₂/FiO_₂, platelet count, and bilirubin showed ordered adverse-direction contributions, with the largest single interval contribution assigned to PaO_₂/FiO_₂ ≤49.9 mm Hg (Table 4). Creatinine contributions increased across thresholds of 1.085, 1.425, and 1.665 mg/dL. These outcome-derived cutpoints should not be interpreted as revised physiological definitions of respiratory, coagulation, hepatic, or renal dysfunction.

Other components departed more substantially from conventional SOFA grading. The neurological component retained a single cutpoint, assigning 9 points at GCS ≤4.5 and no additional points above this threshold. Thus, the fitted model concentrated neurological prognostic contribution at the most severe end of the GCS distribution rather than reproducing the graded neurological categories of SOFA. MAP remained bidirectional for the reason described in the Methods section: the APACHE-derived source field captures the first-day value associated with the greatest APACHE abnormality and may therefore represent either a hypotensive or hypertensive extreme. Among vasoactive therapies, norepinephrine and epinephrine each contributed 6 points, dopamine 2 points, and dobutamine 0 points. The zero dobutamine contribution should not be interpreted as evidence of clinical neutrality; it indicates that the fitted additive model assigned no additional prognostic contribution to this binary exposure variable. As vasoactive therapies were available only as binary indicators rather than harmonized dose measurements, these weights should be interpreted as cohort-specific prognostic contributions rather than dose-equivalent cardiovascular severity categories.

Overall, these cutpoints and weights characterize prognostic structure within this cohort under the prespecified model constraints and are not proposed as new physiological definitions of organ dysfunction. The recalibrated SOFA-derived score should therefore be interpreted as a mortality risk model rather than as a replacement for the descriptive organ dysfunction grading provided by conventional SOFA.

The recalibrated SOFA-derived model retained full interpretability: each prediction decomposes into interval-specific contributions that sum to an additive score, which is then mapped to mortality risk through logistic calibration, ensuring a monotonic relationship between score and predicted mortality.

Several variables (PaO_₂/FiO_₂, GCS, and platelets) showed clear separation between extreme intervals, with adverse ranges receiving large positive contributions. Creatinine and bilirubin also displayed ordered increases in contribution across adverse ranges, whereas MAP showed a less ordered interval-specific profile. Vasoactive therapies, encoded as binary indicators, showed variable patterns: norepinephrine and epinephrine each contributed 6 points when present; dopamine contributed 2 points; and dobutamine contributed 0 points. Mechanical ventilation contributed 3 points when present.

The model remains fully transparent: changes in predicted risk can be attributed directly to specific clinical ranges and their additive contributions. From a clinical decision-making perspective, decision curve analysis was performed using the out-of-fold predictions from the hospital-grouped cross-validation and is shown stratified by LOS and overall cohort (Figure S3 in Multimedia Appendix 1).

As a separate analysis of optimization stability, the entire model-development pipeline—automatic discretization, weight optimization, and integer projection—was repeated 1000 times using different random seeds for patient-level data partitioning and weight initialization. The original SOFA-derived score was evaluated on the identical test patients in each replication to enable paired comparison. All 1000 attempted replications satisfied the prespecified convergence criteria; none was rejected or repeated. Across replications, the recalibrated score achieved a median AUROC of 0.814 (IQR 0.808‐0.821; range 0.784‐0.839) and a median AUPRC of 0.502 (IQR 0.489‐0.514; range 0.438‐0.569), compared with 0.754 (IQR 0.747‐0.761; range 0.717‐0.780) and 0.387 (IQR 0.376‐0.397; range 0.341‐0.443), respectively, for the original score. The median paired improvement was 0.060 AUROC (IQR 0.055‐0.066; range 0.027‐0.088) and 0.116 AUPRC (IQR 0.106‐0.125 range 0.061‐0.153), with both differences positive in all 1000 replications. These replications characterize the stability of the optimization procedure under patient-level resampling and are distinct from the hospital-grouped cross-validation used for the primary performance assessment (Table 5).

The direction of the estimated interval weights was highly stable across the 1000 replications. Vasoactive drug weights also showed limited variability: norepinephrine had a median weight of 6 points (IQR 6‐6), epinephrine 5 (IQR 5‐6), dopamine 2 (IQR 2‐2), and dobutamine 0 (IQR 0‐1); mechanical ventilation contributed a median of 2 points (IQR 2‐2).

Cutpoint stability varied by predictor and ordinal position. Several boundaries were highly consistent, including the first GCS boundary (median 4.5, IQR 4.5‐4.5), the first MAP boundary (40.5 mm Hg, IQR 40.5‐40.5), and the first and fourth creatinine boundaries (1.09 mg/dL, IQR 1.09‐1.10; and 5.07 mg/dL, IQR 5.06‐5.13, respectively). Other boundaries, particularly for PaO_₂/FiO_₂, bilirubin, and platelets, showed greater variability or were selected in only a subset of replications.

Variation in individual cutpoints had a more limited effect on the resulting point functions. Across the evaluated predictor ranges, the mean interquartile width of the assigned points was 0.27 points for platelets, 0.55 for MAP, 0.65 for PaO_₂/FiO_₂, 0.71 for GCS, and 0.95 for creatinine. Bilirubin showed greater variability, with a mean interquartile width of 2.29 points.

At the patient level, total scores remained highly concordant across replications. The median within-patient IQR was 2 points on a median score range of 0 to 62, and 63.1% (9367/14,848) of patients had a score IQR of ≤2 points. The median pairwise Spearman correlation between scores generated by different replications was 0.970 (IQR 0.964‐0.975), indicating high stability of patient ranking despite variability in some individual cutpoints. Full distributions of cutpoints and weights are provided in Tables S7 and S8 in Multimedia Appendix 1.

The sensitivity analysis included 61,469 encounters from 186 hospitals, including 46,621 encounters with 1 imputed determination. We used the same 10-fold nested cross-validation procedure as in the primary analysis, with folds grouped by hospital. As the number of patients differed, the folds were redefined for this analysis. The recalibrated score achieved a higher AUROC than the SOFA-derived reference score in all 10 outer folds. Median AUROC was 0.853 (IQR 0.847‐0.857) for the recalibrated score and 0.818 (IQR 0.811‐0.821) for the original SOFA-derived score, corresponding to a median paired difference of 0.037 (IQR 0.034‐0.043; range 0.028‐0.046). Median AUPRC was 0.395 (IQR 0.372‐0.425) versus 0.300 (IQR 0.282‐0.320), with a median paired difference of 0.096 (IQR 0.085‐0.110). As the expanded cohort differed substantially from the primary complete-case cohort in case mix and event prevalence (8% vs 17%), this analysis was interpreted as a robustness assessment of the paired comparison rather than as a direct comparison of absolute performance between cohorts.

To facilitate clinical adoption and prospective validation, we developed an open-source, web-based calculator that implements the recalibrated SOFA-derived score. The calculator is freely available online [25] and performs real-time computation of both the recalibrated score and mortality probability estimates. The interface accepts the same clinical variables used in the SOFA-derived reference score and applies the data-driven thresholds and weights derived in this study. All calculations are performed locally in the web browser, ensuring data privacy and enabling offline use. The calculator provides color-coded risk stratification and is designed for integration into existing clinical workflows, although this integration has not been evaluated prospectively. Source code for the calculator is included in the project repository to enable local deployment and adaptation to institution-specific requirements.

This study demonstrates that data-driven recalibration of a SOFA-derived score using optimization and ML components can substantially improve mortality prediction in contemporary ICU populations while preserving an additive, transparent, and interpretable structure. The recalibrated SOFA-derived score achieved consistent improvements in discrimination across all LOS strata, with AUROC gains of 6 percentage points and AUPRC improvements of 10 percentage points compared with the SOFA-derived reference score. Although global calibration was adequate after unweighted logistic recalibration, calibration was not consistently superior to that of the reference score. By learning data-adaptive thresholds and interval-specific contributions for each SOFA component, the proposed approach provides a mechanism for updating expert-defined scoring rules as clinical populations evolve [26,27] while preserving an additive structure grounded in SOFA-related clinical domains. As the recalibrated score uses its own 0 to 71 scale, treats several predictors separately, and does not reconstruct the original 0 to 24 SOFA score, it is not directly comparable with the canonical SOFA. The proposed framework provides an interpretable approach to updating SOFA-based risk scores.

Previous studies applying ML to ICU mortality prediction have reported AUROC values of 0.80 to 0.90, although often at the cost of interpretability [28]. Black box models face significant barriers to clinical adoption due to their opacity [18]. Interpretable approaches such as AutoScore [19] and FasterRisk [20] have demonstrated that transparent, point-based models can achieve competitive performance, but these methods construct entirely new scores rather than recalibrating existing instruments.

This work differs from de novo score construction in that it recalibrates a prespecified SOFA-derived structure rather than learning both the clinical framework and scoring architecture from scratch. Clinicians have decades of experience interpreting SOFA scores, and this familiarity represents institutional knowledge that would be lost with an entirely novel scoring system. The recent SOFA 2.0 [3] combined expert consensus with empirical analysis of more than one million ICU admissions, representing a major advance. However, its development required substantial expert involvement and resource-intensive consensus processes. Our automated, outcome-driven framework offers a complementary approach enabling more frequent updates as populations evolve.

From a methodological perspective, this work highlights the value of combining supervised ML with explicit structural constraints. Rather than replacing SOFA with a complex black box model, the proposed approach retains its familiar point-based logic and clinical interpretability [18], while allowing both thresholds and weights to be optimized from data [19]. The key innovation lies in the optimization of discretization boundaries and interval weights under constraints that preserve the additive, transparent structure of classical clinical scores. The integer constraint ensures manual computability; the scale-anchoring penalty prevents weight collapse; and automatic interval selection prevents overfitting while allowing variable-specific granularity. To our knowledge, this represents one of the first applications of a constrained optimization framework with ML components to recalibrate a prespecified SOFA-derived organ dysfunction score while preserving additivity and transparency.

Clinically, this recalibrated SOFA-derived model may offer a practical pathway for modernizing severity assessment in the ICU. As it relies on routinely collected SOFA-related variables and retains an additive structure, it could be integrated into existing clinical workflows or EHR systems with minimal friction.

The potential applications extend across multiple domains: (1) prognostication and family communication, where accurate risk estimates support shared decision-making [2]; (2) clinical trial enrollment and stratification, where improved discrimination enhances trial efficiency [29]; (3) quality benchmarking, where risk adjustment supports institutional comparisons [2]; and (4) sepsis identification, given SOFA’s central role in Sepsis-3 definitions [6].

Improved risk stratification achieved through transparent recalibration supports evidence-based decision-making without requiring clinicians to rely on opaque prediction models [18]. This alignment with the growing emphasis on explainable AI in medicine is crucial: regulatory bodies and frontline clinicians increasingly demand transparency in algorithmic decision support tools [18,30].

From an implementation perspective, the model’s computational requirements are minimal—only threshold evaluations and integer summations—enabling real-time deployment within existing EHR infrastructure without specialized hardware [19]. Unlike deep learning models requiring graphics processing unit acceleration, the recalibrated SOFA-derived score can be implemented as a simple lookup table, ensuring compatibility with legacy systems.

This additive structure may facilitate regulatory approval and clinical adoption under existing clinical decision support governance frameworks [18]. Future implementations should consider automated pipelines for periodic recalibration using local institutional data [27], enabling adaptive risk scores that maintain calibration as patient populations evolve [26]. The COVID-19 pandemic illustrated how rapidly populations can change; an automated recalibration framework could enable rapid adaptation to such shifts [31].

SOFA serves purposes beyond mortality prediction, including standardized description of organ dysfunction, sequential assessment, and its role in established clinical definitions such as Sepsis-3. This study addresses mortality risk estimation rather than these broader functions. Accordingly, the recalibrated SOFA-derived score is not intended to replace the canonical SOFA for organ dysfunction assessment. A potential future use would be to report the recalibrated risk score alongside conventional SOFA assessment when the required SOFA-related inputs are available electronically, without requiring a separate set of physiological measurements. This complementary use has not been evaluated in this study and would require external and prospective validation before clinical implementation.

It is important to distinguish periodic model updating from system-level adaptability. This study demonstrates the former: the score parameters can be reestimated using contemporary data while preserving the underlying additive framework. System-level adaptability, whereby the surrounding digital system can modify its behavior, content, or configuration in response to changing program requirements and user needs, represents a distinct property described in adaptable mobile health architectures [32] and was not evaluated in this study. Whether the recalibrated score can be effectively incorporated into clinical workflows, and whether the accompanying calculator supports real-world clinical use, requires prospective implementation and workflow evaluation.

Patterns of measurement availability can themselves carry predictive information in routinely collected ICU data [33,34]. We did not incorporate missingness indicators into the recalibrated score because the objective of this study was to update a prespecified SOFA-derived structure using physiological and treatment-related variables. Adding indicators of whether specific measurements were obtained would introduce a distinct measurement-process signal that may reflect local testing and documentation practices in addition to patient condition, potentially limiting transportability across hospitals [34,35]. Instead, we assessed the robustness of the findings to an alternative missing data assumption through the prespecified sensitivity analysis in the expanded cohort, in which the paired performance advantage was preserved. Explicitly modeling measurement processes together with physiology remains a distinct and potentially valuable direction for future prediction research.

Complete-case selection is an important limitation of this study. Encounters in which all required measurements were available differed from those in which 1 or more were unavailable. The high frequency of missing PaO_₂/FiO_₂ and bilirubin measurements, together with these between-group differences, is consistent with clinically selective measurement; however, the database does not establish why an individual measurement was or was not obtained. Complete-case selection may therefore have enriched the primary analytic cohort for more severely ill and intensively monitored patients and may affect both absolute and comparative performance estimates. Accordingly, the primary performance estimates should not be extrapolated directly to unselected ICU admissions. Nevertheless, the direction of the paired improvement over the SOFA-derived reference score was preserved in a prespecified sensitivity analysis that expanded the cohort to 61,469 encounters under an alternative normal-range imputation assumption for a single missing determination.

Real-world implementations of ML-based risk scores demonstrate that successful clinical integration depends on both technical performance and alignment with existing workflows [36]—tools perceived as disrupting established workflows face resistance regardless of predictive accuracy. By retaining an additive SOFA-derived structure, our approach may reduce this friction.

Future work will focus on external validation in independent cohorts, including international datasets such as MIMIC-IV, AmsterdamUMCdb, and HiRID, as well as extending the framework to longitudinal SOFA trajectories. Additional directions include periodic recalibration to address concept drift [26,27] and prospective evaluation of clinical utility in decision support and trial enrollment contexts. The framework will also be investigated for other additive scores, including APACHE, SAPS, the quick Sequential Organ Failure Assessment, National Early Warning Score, Model for End-Stage Liver Disease, and CURB-65 (confusion, urea, respiratory rate, blood pressure, and age ≥65 years).

This work should be interpreted in light of several limitations. First, the study is retrospective and based on a single database. Although eICU-CRD encompasses 208 US hospitals, external validation in independent international cohorts is required to confirm generalizability [24]. Second, the analysis focuses on static, day 1 SOFA assessment and does not capture temporal dynamics or trajectories of organ dysfunction. Trends in SOFA scores over time carry prognostic information beyond single-point assessments that our framework does not yet incorporate. Third, mortality was the sole outcome; impact on length of stay, functional outcomes, or treatment decisions remains unexplored. Fourth, we could not evaluate SOFA 2.0 as required variables were not consistently captured. Fifth, learned thresholds may reflect dataset-specific practice patterns [26], underscoring the need for validation in other settings. Sixth, because vasoactive-agent dose units could not be harmonized across hospitals and urine output was not incorporated, the reference comparator was a SOFA-derived operationalization rather than an exact reconstruction of the original SOFA. Its cardiovascular component used a prespecified ordinal mapping of MAP and binary vasoactive-agent administration, whereas its renal component used creatinine alone. Accordingly, performance comparisons in this study should be interpreted against this SOFA-derived reference score rather than the canonical original SOFA. Finally, this study evaluates discriminative performance but does not assess clinical utility prospectively.

In summary, we have demonstrated that interpretable ML-based recalibration offers a practical and robust strategy for modernizing established clinical scores such as SOFA. By improving discrimination, achieving adequate global calibration after unweighted logistic recalibration, and maintaining transparency through an additive structure grounded in SOFA-related clinical domains, this approach provides a viable path toward responsible integration of ML into routine critical care assessment.

The framework demonstrates that the dichotomy between predictive performance and interpretability is false: with appropriate methodological design, meaningful improvements in risk prediction can be achieved without sacrificing transparency. This evolutionary approach—leveraging ML to optimize parameters of established instruments while preserving trusted structures—offers a promising route to realizing AI benefits in critical care while maintaining the accountability that patient safety demands.