The formal study of DTRs emerged from the causal inference literature in biostatistics, building on Robins’s [28] pioneering work on g-computation and structural nested models. Murphy [3] formalized the DTR optimization problem within a potential outcomes framework and introduced A-learning as a semiparametric approach for estimating optimal treatment rules at each stage. Subsequently, Murphy [6] demonstrated the connection between DTR estimation and Q-learning from the RL literature, establishing backward induction with regression-based Q-functions as a practical estimation strategy. This Q-learning approach was further developed by Zhao et al. [29], who studied its theoretical properties, and by Moodie et al. [30], who extended it to accommodate continuous-valued treatments.

The Sequential Multiple Assignment Randomized Trial (SMART) design (Lei et al. [31], Almirall et al. [32]) was developed as the experimental gold standard for comparing embedded DTRs, providing unbiased estimation of regime-specific means under sequential randomization. However, SMART designs are expensive and ethically constrained, motivating the development of methods that can learn DTRs from observational data under the assumption of sequential ignorability (no unmeasured confounders).

OWL [9] and its extensions recast DTR estimation as a weighted classification problem, connecting to support vector machine methods and enabling the use of flexible nonparametric classifiers. Zhang et al. [33] proposed robust OWL methods to handle misspecification, while Zhao et al. [34] extended OWL to multi-stage settings. More recently, tree-based methods [35], [36] have offered interpretable DTR estimation, though they may sacrifice optimality for interpretability.

The Chakraborty and Moodie [5] monograph provides a comprehensive treatment of statistical methods for DTR estimation, covering both parametric and semiparametric approaches, while Tsiatis et al. [37] offer an updated treatment incorporating recent methodological developments. Despite this rich methodological landscape, classical DTR methods have limited capacity to handle the very high-dimensional state spaces and complex nonlinear dynamics characteristic of modern clinical data, motivating the adoption of deep learning-based approaches.

The application of RL to clinical decision-making has accelerated dramatically since the mid-2010s, driven by advances in deep RL and the increasing availability of large-scale EHR datasets. The seminal work of Raghu et al. [10] applied fitted Q-iteration with deep neural networks to optimize vasopressor and intravenous fluid dosing for sepsis patients in the ICU, using retrospective data from the MIMIC-III database. This study demonstrated that RL-derived policies could potentially reduce estimated mortality compared to clinician behavior, though the evaluation relied entirely on off-policy estimation methods. Komorowski et al. [38] extended this line of work with a more clinically nuanced sepsis RL model, incorporating physician expertise into the reward design and performing extensive OPE.

RL has been applied across numerous clinical domains, including: mechanical ventilation weaning [39], where RL policies identified potentially superior weaning strategies; chemotherapy dosing optimization [11], where model-based RL approaches planned drug administration schedules to balance tumor reduction and toxicity; heparin dosing [40], using actor-critic methods for anticoagulation management; and glycemic control in diabetes [13], [41], where RL agents learned insulin dosing policies from continuous glucose monitoring data. The critical challenge of learning policies from fixed, previously collected clinical datasets has driven significant interest in offline (batch) RL methods. Batch-constrained Q-learning (BCQ) [16] addresses distributional shift by constraining the learned policy to remain close to the behavior policy observed in the data. Conservative Q-learning (CQL) [17] takes an alternative approach, learning a lower bound on Q-values to penalize out-of-distribution actions. Behavior-Regularized Actor-Critic (BRAC) [42] regularizes the actor update to stay close to the behavior policy. These methods have shown improved stability and reliability in offline settings, though they can be overly conservative, potentially failing to discover beneficial treatment strategies that deviate from historical clinician practice.

Safety considerations in RL-based clinical decision support have received increasing attention. Constrained MDP formulations [18], [43] incorporate explicit safety constraints into the optimization objective, while action masking approaches [44] prevent the selection of clinically contraindicated actions. Gottesman et al. [14] provided essential guidelines for the responsible development and evaluation of RL-based treatment policies, emphasizing the limitations of OPE and the dangers of overinterpreting estimated policy improvements.

Metaheuristic optimization algorithms have a long history of application in healthcare and biomedical engineering, though primarily for single-stage optimization problems rather than sequential decision-making. Genetic Algorithms (GAs) have been applied to treatment scheduling problems, including radiotherapy fractionation optimization [22], surgical scheduling [45], and chemotherapy regimen selection [21]. Particle Swarm Optimization (PSO) has found applications in medical image segmentation [46], intensity-modulated radiation therapy planning [47], and pharmacokinetic parameter estimation [48].

DE, introduced by Storn and Price [19], has been applied to drug dosing optimization [23], PID controller tuning for glucose regulation systems [49], and model parameter estimation in physiological systems [50]. DE is particularly attractive for policy optimization because it handles continuous parameter spaces efficiently, requires minimal algorithm-specific tuning, and maintains population diversity through its mutation-crossover-selection mechanism [20].

Multi-objective evolutionary optimization has also been applied to clinical problems where multiple competing clinical objectives must be balanced, such as simultaneous minimization of treatment toxicity and maximization of efficacy in cancer treatment [51]. However, applying metaheuristic methods to sequential clinical decision-making under the DTR framework remains largely unexplored.

The integration of evolutionary search with gradient-based learning has emerged as a productive research direction in the broader machine learning community. Neuroevolution approaches [24] evolve both the architecture and weights of neural networks, with methods such as NEAT [52] and HyperNEAT [53] demonstrating success on RL benchmarks. Salimans et al. [26] showed that Evolution Strategies (ES) can serve as a scalable alternative to traditional RL algorithms for policy optimization, achieving competitive performance on Atari and MuJoCo tasks while offering superior parallelizability.

PBT [25] combines evolutionary hyperparameter selection with gradient-based weight optimization, allowing simultaneous optimization of the training process itself. Khadka and Tumer [27] proposed Evolutionary Reinforcement Learning (ERL), which maintains a population of actors evolved using genetic operators alongside a gradient-based RL agent, with periodic knowledge transfer between the two components. Pourchot and Sigaud [54] introduced CEM-RL, combining the cross-entropy method with twin delayed deep deterministic policy gradient (TD3) for continuous control. These hybrid approaches consistently demonstrate that combining global population-based search with local gradient-based refinement yields superior performance compared to either approach alone.

Despite this promising body of work, no prior study has systematically developed a hybrid metaheuristic–RL framework specifically designed for clinical DTR optimization. The present work fills this gap by introducing METRO-RL, which adapts the hybrid paradigm to the unique requirements of healthcare applications, including clinical safety constraints, patient heterogeneity, offline evaluation, and the need for interpretable and trustworthy policy recommendations.

We formulate the personalized treatment policy optimization problem as a finite-horizon Markov Decision Process (MDP) defined by the tuple M = (S, A, P, R, γ, T), where S denotes the state space representing patient clinical features, A denotes the discrete action space of available treatment options, P: S × A × S → [0, 1] represents the state transition dynamics capturing disease progression under treatment, R: S × A × S → R is the reward function encoding clinical outcomes, γ ∈ [0, 1] is the discount factor, and T is the treatment horizon (number of decision stages).

A treatment policy π: S → A maps the current patient state to a treatment action. The objective is to find the optimal policy π* that maximizes the expected cumulative discounted reward:

π ∗ = arg maxπ Eπ [ Σt = 0T γt R(st, at, st + 1)], (1)

The policy is parameterized by a neural network with weights θ, yielding πθ(s) = arg maxa Qθ(s, a), where Qθ is the action-value function approximated by the network.

The METRO-RL framework employs a two-stage architecture that exploits the complementary strengths of metaheuristic global search and deep RL local refinement. Stage 1 uses a modified DE algorithm to search the policy parameter space globally, identifying promising regions that are subsequently refined by a modified DQN in Stage 2. This two-stage approach is motivated by the observation that gradient-based RL methods are highly sensitive to initialization [55], and that providing a strong initial parameter configuration from global search substantially improves the quality and stability of the final policy.

Stage 1 (Global search via modified DE). We maintain a population of N = 100 candidate policies, each encoded as the weight vector θi of a Small Multilayer Perceptron (MLP) with 2 hidden layers of 64 units each and ReLU activations. The network maps the patient state s ∈ S to action-value estimates Qθ(s, a) for all actions a ∈ A. The DE operates through iterative mutation, crossover, and selection to evolve this population toward higher-fitness policy parameters.

A central innovation of METRO-RL is a patient similarity kernel that guides knowledge transfer across patient subgroups in both optimization stages. For patients pi and pj, the composite similarity kernel is defined as:

K(pi, pj) = wd · Kdemo(pi, pj) + wc · Kclin(pi, pj) + wm · Kmorb(pi, pj), (7)

where wd + wc + wm = 1 are learned kernel weights (default: wd = 0.2, wc = 0.5, wm = 0.3).

Demographic kernel Kdemo computes similarity over static patient features (age, sex, BMI, baseline eGFR) using a Gaussian Radial Basis Function (RBF):

Kdemo(pi, pj) = exp(− ||xidemo − xjdemo||² / (2σd²)). (8)

Clinical trajectory kernel Kclin uses Dynamic Time Warping (DTW) to measure the alignment cost between multivariate clinical time series (vital signs, lab trajectories):

Kclin(pi, pj) = exp(− DTW(Ti, Tj) / σc), (9)

where Ti and Tj are the clinical time series of the two patients, and σc is a bandwidth parameter estimated from the data.

Comorbidity kernel Kmorb computes the Jaccard similarity coefficient over binary comorbidity profiles encoded from ICD-10 diagnosis codes:

Kmorb(pi, pj) = |ICDi ∩ ICDj| / |ICDi ∪ ICDj|, (10)

The composite kernel is used for two purposes. In Stage 1, spectral clustering on the kernel matrix partitions patients into K = 10 subgroups, and successful policy parameters from high-performing individuals in one cluster are probabilistically transferred to similar clusters during the DE mutation step (lines 17–19 in Algorithm 1). In Stage 2, the kernel guides the construction of patient-specific replay buffers. When training the DQN on a particular patient subgroup, transitions from patients with high kernel similarity are prioritized in the experience replay buffer, effectively creating a personalized training dataset that enriches the learning signal for underrepresented patient types.

The design of clinically meaningful and computationally tractable reward functions is critical for the success of any RL-based treatment policy. We design domain-specific reward functions for each clinical application that balance treatment efficacy, safety, and patient burden.

Since prospective clinical deployment is infeasible without prior regulatory approval and clinical trial validation, we employ a comprehensive battery of OPE methods to estimate the value of candidate policies using the historical clinical dataset. Our OPE battery includes four complementary estimators:

We formulate clinical safety within the CMDP framework. Let Ck(s, a) denote the k-th constraint cost function, with corresponding threshold dk. The constrained optimization problem is:

maxπ Eπ[Σ γt Rt] subject to Eπ[Σ γt Ck(st, at)] ≤ dk, ∀ k. (14)

We address constraints at two levels. Hard constraints represent absolute clinical contraindications and are enforced via action masking, which sets Q(s, a) = −∞ for contraindicated actions in both the DE fitness evaluation and DQN action selection. The action-masking mechanism ensures zero probability of selecting a dangerous action at any point during optimization or deployment. Soft constraints represent relative clinical warnings and are handled via Lagrangian relaxation, augmenting the reward with penalty terms:

Raug(s, a, s ) = R(s, a, s ) − Σk μk · Ck(s, a), (15)

where μk ≥ 0 are Lagrange multipliers updated via dual gradient ascent during Stage 2 training to adaptively tighten or relax soft constraints based on observed constraint violations.

The T2D dataset comprises a synthetic cohort of 10,000 virtual patients generated using a disease progression model inspired by the UVA/Padova glucose-insulin simulator [61], [62]. While the original UVA/Padova simulator focuses on Type 1 Diabetes with minute-level glucose dynamics, we extend the approach to model T2D progression at a monthly granularity, incorporating long-term disease evolution features including progressive beta-cell function decline, insulin resistance dynamics, weight changes, and renal function trajectories. Patient parameters are sampled from distributions calibrated to the United Kingdom Prospective Diabetes Study (UKPDS) cohort demographics [63]: age uniformly distributed between 40 and 75 years, BMI between 25 and 40 kg/m², baseline HbA1c between 7.0% and 10.5%, and baseline eGFR between 45 and 120 mL/min/1.73m². Comorbidities (cardiovascular disease, neuropathy, nephropathy, retinopathy) are assigned probabilistically based on age, diabetes duration, and baseline HbA1c, using prevalence rates from the Diabetes Control and Complications Trial and UKPDS studies.

The treatment horizon is 12 months with monthly decision points. At each decision point, the clinician (or policy) selects from 20 discrete treatment options: 5 metformin dose levels (0, 500, 1000, 1500, 2000 mg/day) combined with 4 insulin dose levels (0, 10, 20, 40 units/day). The simulator propagates the patient state forward by one month, incorporating stochastic noise to model individual variability in drug response, adherence fluctuations, and lifestyle factors. The dataset is split into 7,000 patients for training, 1,500 for validation, and 1,500 for testing.

The sepsis cohort is extracted from the MIMIC-IV database [64], a freely available critical care database comprising de-identified health-related data associated with patients admitted to the Intensive Care Units (ICUs) of the Beth Israel Deaconess Medical Center. We identify sepsis episodes using the Sepsis-3 criteria [65]: suspected infection (defined by concurrent antibiotic administration and body fluid cultures) combined with an acute change in SOFA score of ≥2 points. After applying inclusion criteria (age ≥18, ICU stay ≥24 hours, sufficient data completeness ≥80% of features recorded) and exclusion criteria (readmissions, transfers from other ICUs, patients with do-not-resuscitate orders at admission), we obtain approximately 18,500 sepsis episodes.

The 48-hour treatment window following sepsis onset is discretized into 12 four-hour intervals. The 46- dimensional state vector is computed at each interval using the most recent available measurements, with forward-filling for missing values and population-median imputation for features never recorded. States are normalized to zero mean and unit variance using training set statistics. Treatment actions (vasopressor dose and IV fluid volume within each 4-hour interval) are discretized into 5 bins each using quintiles of the clinician dose distribution, yielding 25 discrete actions. The primary outcome is 90-day all-cause mortality, determined from hospital records and the Social Security Death Index. The dataset is split by admission year: 2008–2017 for training (approximately 14,800 episodes), 2018–2019 for validation (approximately 1,850 episodes), and 2020–2022 for testing (approximately 1,850 episodes).

We compare METRO-RL against eight baseline methods spanning pure RL, offline RL, and pure evolutionary approaches:

I. Standard DQN (Mnih et al. [66]): DQN with experience replay and target network, initialized with random weights.

II. Dueling DQN (Wang et al. [56]): DQN with value-advantage stream decomposition, without CQL penalty or DE initialization.

III. BCQ (Fujimoto et al. [16]): BCQ with generative model for action proposal, designed for offline RL. IV. CQL (Kumar et al. [17]): CQL with pessimistic value estimation for out-of-distribution actions.

V. A2C (Mnih et al. [67]): advantage actor-critic with offline training via importance-weighted policy gradient. VI. OpenAI ES (Salimans et al. [26]): evolution strategy with antithetic sampling and fitness shaping, same network architecture as METRO-RL.

VII. CMA-ES (Hansen [68]): covariance matrix adaptation evolution strategy for direct policy search, adapted for neural network policy parameterization.

VIII. Clinician policy: the observed behavior policy estimated from the training data via behavior cloning (multinomial logistic regression).

Additionally, a random policy baseline is included for reference. All neural network-based methods use the same architecture (2 hidden layers, 64 units, ReLU) for fair comparison. All methods (except the clinician and random baselines) are trained on the same training data and evaluated on the same held-out test set.

We evaluate all methods using both OPE metrics and domain-specific clinical outcome metrics.

OPE metrics: WIS estimated value, FQE estimated value, and DR estimated value, each computed on the test set with 95% confidence intervals via bootstrap resampling (1,000 replicates).

T2D clinical metrics: 1) mean HbA1c reduction from baseline over the 12-month horizon, 2) hypoglycemia rate (percentage of monthly time steps with FPG < 70 mg/dL), and 3) composite glycemic control score, defined as: Score = 0.5 · norm(HbA1c_reduction) + 0.3· (1 − norm(hypo_rate)) + 0.2· norm(adherence), where norm denotes min-max normalization across methods.

Sepsis (Clinical metrics): 1) estimated 90-day mortality rate under the learned policy (via WIS), 2) estimated SOFA score reduction (mean change from admission to 48h), 3) estimated vasopressor-free hours within the 48-hour window, and 4) estimated ICU length of stay.

Safety metrics: 1) rate of contraindicated actions (percentage of state-action pairs violating hard clinical constraints), and 2) treatment volatility index, defined as the mean absolute change in treatment intensity between consecutive decision points, normalized by the clinician policy’s volatility.

Table 1 summarizes the hyperparameters for METRO-RL and all baselines. All experiments are repeated over 30 independent runs with different random seeds to assess statistical reliability.

All experiments are conducted on a workstation equipped with an NVIDIA A100 GPU (40 GB), AMD EPYC 7763 processor (64 cores), and 256 GB RAM. The framework is implemented in Python 3.10 using PyTorch 2.1 for neural network components and a custom DE implementation optimized with NumPy vectorized operations.

Table 2 presents the comprehensive results for T2D management across all compared methods. METRO-RL achieves the highest composite glycemic control score (0.782 ± 0.024), representing a 16.2% relative improvement over standard DQN (0.673 ± 0.031) and an 18.8% relative improvement over the clinician policy (0.658 ± 0.035). The mean HbA1c reduction achieved by METRO-RL (1.82% ± 0.14%) substantially exceeds that of the best standalone RL baseline, CQL (1.56% ± 0.18%), while maintaining a lower hypoglycemia rate (3.2% vs. 4.1%).

Notably, METRO-RL achieves both superior efficacy (higher HbA1c reduction) and superior safety (lower hypoglycemia rate) simultaneously, a combination that proves elusive for standalone RL methods, which typically improve efficacy at the cost of increased adverse events. The pure evolutionary methods (OpenAI ES, CMA-ES) perform worse than the offline RL methods, suggesting that gradient-based refinement is essential for achieving high-quality policies in this domain, but that evolutionary initialization (as in METRO-RL) provides a crucial advantage over random initialization.

Table 3 presents the sepsis treatment results on the MIMIC-IV test set. METRO-RL achieves the lowest estimated 90-day mortality (18.7% ± 2.1%), representing an 11.4% relative reduction compared to the clinician policy (24.3% ± 1.9%) and a 7.0% relative reduction compared to CQL (20.1% ± 2.5%). The estimated SOFA score reduction is also most favorable for METRO-RL (−3.8 ± 0.6 points), indicating superior organ function recovery. METRO-RL policies also yield higher estimated vasopressor-free hours (38.2 ± 2.4 out of 48 hours), suggesting that the learned policy achieves better outcomes with less aggressive vasopressor administration.

The performance gap between METRO-RL and the clinician policy is particularly notable given that sepsis management is a domain where experienced intensivists already employ sophisticated, context-dependent treatment strategies. The 5.6 percentage point absolute reduction in estimated mortality suggests that METRO-RL identifies treatment patterns that are systematically superior to historical practice, particularly in the timing and dosing of vasopressor initiation and fluid management during the critical early hours of sepsis treatment.

Table 4 presents the OPE results using three complementary estimators. The consistency of METRO-RL’s superiority across all three OPE methods provides strong evidence that the observed performance gains are not artifacts of a particular estimation methodology.

METRO-RL not only achieves the highest estimated value across all OPE methods but also exhibits the lowest variance, with standard deviations 20–40% smaller than those of standalone DQN methods. This reduced variance reflects the stabilizing effect of DE-based initialization: by starting from a globally optimized policy parameter region, the subsequent DQN refinement is less sensitive to random seed and mini-batch sampling effects. The consistency across WIS, FQE, and DR estimates is encouraging, as disagreement among estimators would indicate potential unreliability in the evaluation. The Spearman rank correlation between WIS and FQE rankings across methods is 0.98, and between FQE and DR is 0.97, indicating strong agreement among the estimators.

Table 5 presents the safety analysis comparing contraindicated action rates and treatment volatility across all methods. METRO-RL achieves the lowest rate of contraindicated actions (0.3% on T2D, 0.4% on sepsis), demonstrating that the integrated safety framework combining action masking, clinical constraint-masked crossover in DE, and the Lagrangian penalty effectively prevents dangerous treatment recommendations. The residual 0.3–0.4% contraindication rate arises from edge cases at the boundary of the masking criteria (e.g., patients transitioning into hypoglycemia within a decision interval).

The treatment volatility index reveals another important advantage of METRO-RL: its policies recommend smoother treatment transitions than both standalone RL methods and evolutionary approaches. The volatility index of 0.82 (T2D) and 0.91 (sepsis) indicates that METRO-RL produces treatment trajectories that are 8– 18% smoother than clinician behavior. This reduced volatility is clinically desirable, as frequent, large changes in medication dosing are associated with adverse events, reduced patient adherence, and increased nursing workload. In contrast, standard DQN exhibits a volatility index of 1.34–1.45, and evolutionary methods show even higher volatility (1.78–1.91), reflecting the lack of temporal consistency in policies that are not explicitly optimized for smooth treatment transitions.

Table 6 presents the ablation study results, quantifying the contribution of each METRO-RL component to overall performance. Results are reported for the sepsis domain (MIMIC-IV test set), with similar patterns observed for T2D (omitted for brevity).

The ablation reveals several important findings. First, both stages are essential: removing the DE stage (using random initialization for DQN) increases estimated mortality by 3.1 percentage points, while removing the DQN stage (using DE-only policies) increases mortality by 3.8 percentage points. The ablation results confirm the complementary nature of global search and local refinement. Second, the patient similarity kernel contributes a 1.6 percentage point mortality reduction, validating the benefit of knowledge transfer across similar patient subgroups. Third, removing safety constraints slightly improves the estimated mortality (by 0.4 percentage points), reflecting the inherent trade-off between safety and aggressive optimization; however, the contraindicated action rate increases dramatically from 0.4% to 4.8%, making this configuration clinically unacceptable. Fourth, the CQL penalty and state-aware mutation each contribute approximately 1 percentage point of mortality reduction, demonstrating their individual value within the framework.

To assess whether METRO-RL’s advantages are uniform across patient populations or concentrated in specific subgroups, we stratify the sepsis test set by baseline disease severity into three groups based on admission SOFA score: mild (SOFA 2–5, n = 624), moderate (SOFA 6–9, n = 781), and severe (SOFA ≥10, n = 445).

METRO-RL’s advantage is most pronounced in the severe subgroup, where it achieves an estimated 34.8% mortality compared to 44.7% for the clinician policy, a 9.9 percentage point absolute reduction. This finding is clinically significant because severely ill patients represent the population with the greatest unmet therapeutic need and the highest potential for benefit from optimized treatment policies. The patient-state-aware mutation mechanism in Stage 1, which increases exploration intensity for severe states, likely contributes to this disproportionate benefit by encouraging the discovery of unconventional but effective treatment strategies for the most critically ill patients. For mild patients, the advantage is more modest (2.6 percentage point reduction), consistent with the expectation that standard clinician practice is already reasonably effective for lower-acuity patients.

Table 8 reports the wall-clock training time for each method on the sepsis dataset, measured on the hardware described in Section 4.5.

METRO-RL’s total training time of approximately 2 hours and 45 minutes is roughly 5 times longer than standalone DQN methods but substantially faster than CMA-ES (4 hours 5 minutes) and comparable to OpenAI ES (3 hours 10 minutes). The computational overhead is concentrated in Stage 1 (DE), where 100 candidate policies must be evaluated across 200 generations. Importantly, the DE stage is embarrassingly parallel: each candidate policy evaluation is independent, and the current implementation uses 16 parallel workers for fitness evaluation, reducing the effective time per generation. With full utilization of the 64-core processor, Stage 1 could be further accelerated to under 1 hour. Given that METRO-RL training is a one-time offline computation, the 2.75-hour training time is negligible compared to the clinical timescales of the target applications (months for T2D, days for sepsis).

The fundamental advantage of METRO-RL lies in the complementary nature of its two stages. DE’s population-based search provides effective global exploration of the policy parameter space, maintaining diversity and avoiding premature convergence to local optima, a well-known failure mode of gradient-based RL methods [55]. The ablation study (Table 5) confirms that removing DE initialization and relying on random initialization for DQN increases estimated sepsis mortality by 3.1 percentage points, demonstrating that the gradient-based optimizer alone frequently converges to inferior local optima. Conversely, DE alone (without DQN refinement) produces policies that are coarser and less precise than those refined by gradient-based optimization, as evidenced by the 3.8 percentage point mortality increase when the DQN stage is removed. The DE stage efficiently identifies the “basin of attraction” containing the global optimum, while the DQN stage performs the fine-grained descent within that basin.

The patient similarity kernel serves as a structured prior that regularizes both the DE search and the DQN training by encoding the clinical intuition that similar patients should receive similar treatments. Without this kernel, the framework must independently optimize policies for each patient cluster, wasting computational resources on redundant exploration. The kernel’s DTW-based clinical trajectory component is particularly valuable because it captures temporal patterns in disease progression that static demographic features cannot, enabling the identification of patient subgroups that follow similar clinical courses despite potentially different baseline characteristics. The 1.6 percentage point mortality improvement attributable to the kernel (Table 5) validates this design choice.

A critical prerequisite for any AI-based treatment recommendation system is the assurance of patient safety. METRO-RL’s multi-layered safety framework combining hard constraint action masking with soft constraint Lagrangian penalties achieves the lowest rate of contraindicated actions among all evaluated methods (0.3– 0.4%), approaching the level of human clinician practice (0.5–0.6%). The low residual rate of the clinician policy reflects the reality that even experienced physicians occasionally deviate from guidelines in edge cases. We emphasize that the safety framework is a necessary but not sufficient condition for clinical deployment. Before any real-world use, METRO-RL policies would require: 1) extensive prospective validation in randomized clinical trials, 2) regulatory approval under applicable frameworks (e.g., FDA’s Software as a Medical Device guidance), 3) integration into clinical decision support systems that present recommendations to clinicians for review rather than autonomous execution, and 4) continuous monitoring for distributional shift between the training data and the deployment population.

This study has several important limitations that should be considered when interpreting the results. First, all evaluations rely on offline policy estimation using retrospective data and simulated environments. Despite using a battery of four OPE methods, the true clinical impact of METRO-RL policies can only be determined through prospective clinical trials. OPE methods are known to exhibit bias under violations of the positivity assumption (when the behavior policy assigns zero probability to actions selected by the evaluation policy) and under model misspecification [69]. Second, the T2D dataset is synthetically generated from a simplified disease progression model. While the model is calibrated to UKPDS parameters, it cannot capture the full complexity of real-world T2D management, including patient-specific pharmacokinetics, lifestyle variability, and comorbidity interactions. Third, the discrete action space (20 actions for T2D, 25 for sepsis) is a simplification of the continuous dosing decisions that clinicians make in practice. Finer discretization or continuous action spaces could improve clinical granularity but would also increase the computational and statistical complexity of the optimization. Fourth, the framework does not currently model medication adherence, patient preferences, treatment costs, or social determinants of health, all of which substantially influence real-world treatment outcomes. Fifth, the MIMIC-IV sepsis cohort represents a single tertiary academic medical center, and the generalizability of learned policies to other healthcare systems, patient populations, and clinical practices remains to be established.

AI-driven treatment policy optimization raises important ethical concerns. First, there is a risk that policies optimized on historical data may perpetuate or amplify existing disparities in care if the training data reflects systemic biases in treatment allocation across demographic groups [70]. Fairness-aware policy optimization, incorporating group-level equity constraints, is an important direction for future work. Second, the “black box” nature of neural network policies limits clinician understanding of the rationale behind specific treatment recommendations. While our framework does not directly address interpretability, the patient similarity kernel provides partial transparency by linking recommendations to the treatment histories of similar patients. Third, deploying AI-based treatment recommendations requires careful attention to liability, informed consent, and the appropriate division of responsibility between the AI system and the supervising clinician [71].

Since the inception of this work, several offline RL methods have been proposed that could potentially be integrated into the METRO-RL framework. Decision transformer (Chen et al. [72]) recasts RL as a sequence modeling problem using transformer architectures, while implicit Q-learning (IQL) (Kostrikov et al. [73]) avoids querying out-of-distribution actions entirely by learning Q-values through expectile regression. Diffusion-based offline RL methods (Wang et al. [74]) model complex policy distributions using denoising diffusion models. These methods could replace or complement the DQN component in Stage 2, and exploring such integrations represents a promising direction for future work.

Several extensions of METRO-RL merit investigation. First, extending to continuous action spaces using actor-critic architectures in Stage 2 would enable finer-grained dosing recommendations. Second, incorporating real-time adaptation mechanisms would allow policies to update as new patient data becomes available during the treatment course, moving toward fully online personalization. Third, extending the framework to multi-morbidity settings, where patients simultaneously manage multiple chronic conditions with potentially interacting treatments, would enhance clinical applicability.

Fourth, integration with clinical workflow systems, including natural language explanations of treatment recommendations and visualization of expected outcomes under alternative treatment strategies, would improve clinician trust and adoption. Fifth, prospective clinical trials, initially as non-interventional studies comparing METRO-RL recommendations to actual treatment decisions, would provide the essential evidence base for clinical translation.