2026-5-18
Prospective multi-pathogen disease forecasting using autonomous LLM-guided tree search Sarah Martinson1,2 , Michael P. Brenner1,2 , Martyna Plomecka3 , Brian P. Williams1 , Nicholas G. Reich1,4,† and Zahra Shamsi1,† 1 Google Research, 2 School of Engineering and Applied Sciences, Harvard University, 3 Google Deepmind, 4 University of
arXiv:2605.16238v1 [cs.AI] 15 May 2026
Massachusetts
Probabilistic forecasting of infectious diseases is crucial for public health but relies on labor-intensive manual model curation by expert modeling teams. This bespoke development bottlenecks scalability to granular geographic resolutions or emerging pathogens. Here, we present an autonomous system using Large Language Model (LLM)-guided tree search 1 to iteratively generate, evaluate, and optimize executable forecasting software. In a fully prospective, real-time evaluation during the 2025–2026 US respiratory season, the system autonomously discovered methodologically diverse models for influenza, COVID-19, and respiratory syncytial virus (RSV). Aggregating these machine-generated models yielded an ensemble that consistently matched or outperformed the gold-standard, human-curated Centers for Disease Control and Prevention (CDC) hub ensembles out-of-sample. The system successfully navigated data-scarce "cold start" scenarios for RSV. Moreover, controlled retrospective ablations revealed that optimizing log-scale distance metrics prevents reward hacking, while an automated judge-in-the-loop ensures structural fidelity to complex scientific theories. By autonomously translating epidemiological theory into accurate, transparent code, this framework overcomes the modeling labor bottleneck, enabling rapid deployment of expert-level disease forecasting at unprecedented scales.
Introduction Forecasting disease outbreaks is of critical importance for public health monitoring and response planning. In the past 15 years, the sophistication of probabilistic forecasting has increased dramatically in the United States (US). With some centralization of epidemiological forecasting efforts beginning prior to 2020 2–5 , the COVID-19 pandemic accelerated research and integration of forecasts into public health communication and decision-making 6,7 . Much of the innovation in these fields has been driven by collaborative modeling challenges run by governmental agencies or academic groups. These open challenges, many of which follow an emerging set of data standards for the field 8 , accept submissions from any research team, focus on transparency and careful evaluation, and typically integrate with public health decision-makers 6 . This allows different research groups with their own modeling expertise and opinions to contribute to a collective effort to improve disease forecasting and public health monitoring. The United States Centers for Disease Control and Prevention (CDC) runs such collaborative forecasting hubs for major seasonal respiratory diseases, including influenza (FluSight) COVID-19 (COVIDHub), and, more recently, respiratory syncytial virus (RSV; RSVHub). Each week, updated surveillance data is released and participating teams submit standardized probabilistic forecasts. These hubs share a common operational structure, with weekly submissions and jurisdictional targets, with the main outcome incident hospital admissions. Prior seasons have shown that models contributed by academic, governmental and industry research teams are routinely able to outperform simple reference models such as seasonal average and persistence models, a result that has been shown consistently across infectious disease forecasting challenges in the US, United Kingdom and France 7,9–13 .
‡ To whom correspondence should be addressed: [email protected],[email protected]
Prospective multi-pathogen disease forecasting using autonomous LLM-guided tree search
Probabilistic ensemble forecasts that combine these individual contributed models consistently rank near the top of leaderboards 9,11,12,14–16 , and thus serve as the state-of-the-art forecasts for the public health community. There is a huge diversity of different models submitted by individual teams: while some teams implement explicit mechanistic models of disease transmission (e.g. compartmental state-space formulations representing different stages of infection), others use methods from classical statistical frameworks (e.g., auto-regressive or spatio-temporal regression models) or machine learning (ML; e.g., random forests or neural networks). Additionally, many models hybridize these approaches, taking inspiration from a mechanistic understanding of disease transmission (without explicitly building a compartmental model), but with implementations built with standard statistical or ML techniques. Careful post-hoc evaluation has not shown consistent patterns of one kind of model outperforming others 9–11,16 . In this paper, we examine whether Large Language Models (LLMs) in agentic harnesses can create epidemiological prediction models at the same level of skill as those created by teams in the forecast hubs. If true, this could have widespread implications for epidemiological forecasting: the current methodology is labor intensive, requiring many groups to build and maintain operational individual models. Under this system, expanding CDC-style hubs to countries around the world would be highly labor intensive, with the main bottleneck being humans who can code appropriate and accurate models. Instead, automation could make it possible to do this at a vastly larger scale by allowing public health officials with epidemiological expertise but not expert-level programming skills to orchestrate the creation of a diverse array of models for a specific predictive challenge. We build on our previous work 1 , which showed that an LLM-based tree search algorithm called Empirical Research Assistance (ERA) demonstrated competitive performance on COVID-19 forecasting tasks compared to models written by experts when evaluated retrospectively. ERA uses an agentic harness based on Monte Carlo Tree search, inspired by AlphaGo 17 . By framing the creation of empirical software as a “scorable task,” the system iteratively generated, evaluated, and refined Python code to minimize historical forecasting error. When seeded with high-level descriptions of established epidemiological methods extracted from the literature, the search autonomously discovered and hybridized modeling strategies, ultimately generating 14 models that outperformed the official COVIDHub ensemble in retrospective analyses. However, rigorous evaluation in time-series forecasting settings requires prospective, out-of-sample forecasting. Retrospective studies inherently risk data leakage and hindsight bias regarding a pathogen’s overall trajectory, and they often fail to capture the chaotic realities of live forecasting, such as data reporting delays, mid-season shifts in disease dynamics, holiday anomalies, and subsequent data revisions (backfill). For that reason, in this paper we report a fully prospective analysis of the 2025–2026 season across all three CDC Forecast Hubs—FluSight, COVIDHub, and RSVHub. For each hub, ERA generated weekly forecasts that were submitted in real time with auditable time stamps, ensuring no information from future observations could leak into any model. In total, we constructed 142 unique model candidate prompts and executed ERA tree searches for each, producing over 207,500 individual candidate models, and selected 54 of these for an internal prospective hub. Of these, 19 in total were further designated as ensemble components, yielding one submitted ensemble forecast per pathogen (Supplementary Fig. 17). This paper demonstrates that ERA is able to generate, adapt and replicate epidemiological forecasting models across a range of methodologies that, when ensembled, perform competitively with current state-of-the-art methods and with the gold-standard CDC hub ensembles. Additionally, this automated development strategy demonstrates the advantages of a wider and deeper iterative exploration of modeling approaches than current standard practice of manual, bespoke model development allow. It also highlights the value of expertise in guiding model selection and development.
2
Prospective multi-pathogen disease forecasting using autonomous LLM-guided tree search
Results We first outline our main results, showing performance of different models on each of the CDC’s forecasting competitions. We then present results from additional retrospective ablation experiments run on influenza, including: (a) changing the evaluation metric used for scoring; (b) including an explicit LLM-as-judge fidelity gate to encourage instruction-following; and (c) changing the core LLM used to generate the model code. Part I: Real-Time Prospective Model Performance Task, evaluation, and model selection. CDC runs three real-time collaborative forecasting challenges for respiratory diseases in the US: FluSight for influenza, COVIDHub for COVID-19, and RSVHub for RSV. Each week, updated surveillance data are released and participating teams submit standardized probabilistic forecasts of incident hospital admissions across 52 jurisdictions at horizons of zero to three weeks ahead, with each forecast represented as a predictive distribution over 23 quantiles (Fig. 1a–c). Forecasts are evaluated against observations released one to four weeks after submission, making the setting a leak-proof prospective test: no model can access the outcomes it is predicting, and all submissions carry a publicly auditable time stamp on GitHub. We treated this forecasting task as an automated model discovery problem, deploying Empirical Research Assistance (ERA) to generate, evaluate, and iteratively refine code for candidate models that are then used to produce forecasts (Fig. 2). Given a natural language task specification, access to historical surveillance data (including historical values of the target variable and ILINet data, which predates the current target series by more than two decades and thereby offers a longer record of seasonal trends), an evaluation harness with a defined objective function, and a computational sandbox with access to relevant scientific libraries, ERA uses an LLM-guided tree search to autonomously explore the space of possible modeling approaches. At each node of the search tree, the LLM modifies the previous model’s code; the modified model is executed and scored on a validation period; and the score is used to guide subsequent branching decisions. To guard against overfitting to the validation period, all candidate nodes were additionally evaluated on a separate retrospective test set (Fig. 1d), and final models were selected on the basis of combined performance across both periods. The models selected through this process were deployed prospectively starting in November 2025, generating real-time forecasts submitted to the CDC hubs and published weekly with time stamps to the google-research GitHub repository for transparent evaluation over the season (https://github.com/google-research/google-research/tree/master/epi_forecasts). To promote architectural diversity, we supplied ERA with a structured set of prompts spanning three strategic families of modeling approaches: • Single-model adaptations of established methods from the CDC hubs and published literature, whose descriptions are largely based on the corresponding model metadata file submitted to the relevant forecast hub. • “Recombinations”, or double-model adaptations that combine pairs of distinct methodological frameworks. This takes two high performing models and asks the LLM to combine their best features to create a model with a higher score 1 . • “Novel” methods, including unconstrained searches where the LLM received no specific instruction for how to build the model and was free to discover effective approaches, and complex methods suggested for this task by Gemini Deep Research.
3
Prospective multi-pathogen disease forecasting using autonomous LLM-guided tree search
Influenza benefited from the richest existing modeling ecosystem resulting in 107 unique prompts and modeling strategies: its prompt set spanned active FluSight and COVIDHub submissions, historically high-performing methods identified by experts from the published literature, novel architectures generated via Gemini Deep Research, and a structured grid of mechanistic–statistical recombinations. For COVID-19, we used prompts, primarily targeted active COVIDHub submissions, supplemented by models carried over from the prior retrospective study 1 . The RSVhub was new this year so there was a much smaller number of options. We used 16 prompts targeting the small number of available RSVHub methods alongside novel Deep Research architectures; the limited pool of viable base models precluded systematic recombination. From the selected pool, a smaller high-performance subset was designated inclusion in the official ensemble submission to each CDC hub (Supplementary Fig. 17). Models were scored during search using the Weighted Interval Score (WIS; Eq. (1)), a proper scoring rule well-matched to the 23-quantile submission format. 18 All leaderboard comparisons reported in this paper use WIS applied to log-transformed targets and quantiles (referred to throughout as “log WIS”) , a choice that is recommended for epidemiological forecast evaluation 19 and is consistent with the official FluSight evaluation framework. We examine the implications of different metric choices in the controlled retrospective analysis in Part II. At the beginning of the prospective season, we selected 34 influenza models, 19 COVID-19 models, and 2 RSV models from the full pool of generated candidates on the basis of retrospective validation and test performance, with an additional emphasis on maintaining methodological diversity across the selected set. We submitted the predictions from these models weekly with time stamps to the google-research GitHub repository for transparent evaluation throughout the season (https: //github.com/google-research/google-research/tree/master/epi_forecasts). Ensembles of ERA-generated models achieve competitive accuracy across three pathogens. The Google-SAI ensembles—models submitted by Google Research’s Science AI team achieved top-tier performance across all three pathogens (Fig. 3b). The current analysis, completed before the official end of the respiratory season, considers hub submissions up to and including 2 May 2026. We started analysis from the first valid submission of the Google-SAI ensemble to each CDC hub. The task space for all analyses covers 52 jurisdictions (50 states, plus Washington DC and Puerto Rico) and four horizons. For influenza, we evaluate models on a task space including 24 reference dates. Eligible models were those submitting predictions for ≥ 80% of tasks (≥ 3, 744 tasks). The Google-SAI-FluEns ranked first among 43 eligible submissions (out of 57 total models submitting at least one prediction to FluSight) by season-average mean log WIS and ahead of all hub ensembles including the official FluSight-ensemble, which is created by taking the median across all eligible team submissions. For COVID-19, we consider 21 reference dates with the same task space coverage as for influenza, resulting in a minimum required task space of 3, 244 tasks. The Google-SAI-Ensemble ranked first across 12 eligible COVIDHub submissions (out of 17 total), outperforming the CovidHub-ensemble by a clear margin. RSV models are evaluated over 18 reference dates, with the same eligibility criteria, i.e., models must have scorable predictions for ≥ 2, 745 tasks. The Google-SAI-RSVEns ranked first among the four eligible RSVHub submissions (out of seven total), outperforming both the RSVHub-ensemble and the RSVHub-baseline. Jurisdiction-level forecast plots for all three pathogens (Supplementary Figs. 7, 10, 14) illustrate the behavior of these ensemble models across the analyzed reference dates. Nearly all locations saw a clear ‘seasonal wave’ for influenza during the forecasting period. For COVID and RSV, the magnitude of a seasonal rise and fall varied by location. For influenza, the ensemble median roughly tracked the observed seasonal trajectory, although it tended to under-predict as the season initially progressed and over-predict after the season peak. Predictive intervals were wider during the rapid epidemic ascent 4
Prospective multi-pathogen disease forecasting using autonomous LLM-guided tree search
in December 2025 and tightened through the peak and decline phases from January 2026 onward. COVID-19 forecasts were generated over a more gradual seasonal peak in most locations, with the ensemble median following the declining trend reliably across jurisdictions. RSV exhibited the most heterogeneous state-level dynamics, with greater variation in the timing and shape of epidemic curves, and correspondingly wider predictive intervals throughout the season. The standardized rank distribution of FluSight-submitting models (Supplementary Fig. 8) illustrates the variability in individual model performance across tasks. The distribution of standardized ranks show how often a model ranks at the top (rank of 1) or bottom (rank of 0) for predictions made for a specific task. The Google_SAI-FluEns model showed the most consistent performance of any eligible individual model (excluding the FluSight-ensemble), with the highest first quartile rank of any model other than the FluSight-ensemble: this means it had a lower rank less frequently than all other individual models. Several other models were more frequently higher ranked, but also had more frequent lower ranks. As previously shown 11 , the FluSight-ensemble had a narrower distribution of standardized ranks than other models, indicating that it was rarely the best or worst model, which is consistent with it being a model that is a combination of all other models. Ensembles are rarely among the worst performers on any given task–more consistent across conditions than dominant on any one–a profile characteristic of ensembles that aggregate methodologically diverse component models. The Google-SAI ensembles for COVID-19 (Google_SAI-Ensemble, Supplementary Fig. 11) and RSV (Google_SAI-RSVEns, Supplementary Fig. 15) display similar profiles, ranking second among COVIDHub submissions and first among RSVHub submissions on having the smallest first quartile of standardized ranks while also having a high third quartile of ranks. Broken down by forecast horizon, this pattern of competitive consistency is preserved across all prediction windows from 0 to 3 weeks ahead and across all three pathogens (Supplementary Tables 5, 7, and 9). For influenza, on relative and absolute WIS and log WIS, the Google-SAI-FluEns tracks within the top cluster of submitted models at every horizon, well below the CDC baseline and ensemble models, and maintains its relative position without degradation at longer lead times. Discovery of diverse and effective forecasting architectures A superset of Google_SAI-FluEns component models were submitted as individual forecasts to the Google Research GitHub repository throughout the prospective season. Figure 4 shows the relative log WIS of every model in this repository against the respective CDC hub ensemble (values below one indicate better-than-CDC-ensemble performance), alongside several of the top-performing independently submitted models from each hub for context. Across all three hubs, the Google SAI ensemble forecast achieved the best relative log WIS score across all internal component models and all top-performing hub-submitting models (Fig. 4). For influenza, two ERA-generated component models (not including the Google_SAI-FluEns) from our prospective time-stamped repository outperformed the FluSight-ensemble model. Of these two models, one is a single-model adaptation of an existing method (G-Cornell_JHU-hierarchSIR), and the other is a double-model recombination, a hybrid of LANL-DBM and LANL-Inferno. Both of these models were included as component models in Google_SAI-FluEns. For COVID-19, no individual Google model out performed the COVIDHub-Ensemble. For RSV, both individual Google models out-performed the RSVHub-Ensemble. The leading recombination model for influenza merged two methods: the LANL Discrepancy-Based Model (LANL-DBM), which pairs a compartmental Susceptible-Infected-Recovered (SIR) structure with a statistical discrepancy term and LANL-Inferno, which is a Bayesian probabilistic forecasting framework by the same authors 20,21 . Notably, neither model has submitted to CDC hubs in recent
5
Prospective multi-pathogen disease forecasting using autonomous LLM-guided tree search
seasons. Analysis of the recombination model’s source code reveals that although ERA’s generated code does not faithfully implement the mechanistic SIR compartmental structure at the core of the original LANL-DBM framework that is used to generate a baseline ‘expected’ shape of a given season—instead, this is replaced by a cyclical regression model—it preserves several of the core statistical ideas of the target papers. The ERA-generated model mirrors the original formulation by modeling in logit space and adopting a similar three-part additive structure to decompose baseline trends, global discrepancies, and local residuals. Additionally, it generates probabilistic forecasts through simulated sampling and borrows seasonal information from the longer historical ILI record to compensate for the short target hospitalization time series. For COVID-19, the model pool combined newly generated ERA models with five models carried over and fine-tuned from our prior retrospective study 1 : two novel architectures discovered through unconstrained Gemini Deep Research searches (G-DeepResearch_CounterfactualSimulation and G-DeepResearch_RegimeSwitchingDetection), and three double-adapted hybrids combining pairs of established COVIDHub methods (G-CMU_TimeSeries-UMass_gbqr; G-CMU_climate-baselineUMass_ar6_pooled; G-CEPH_Rtrend_covid-CMU_climate_baseline). In the prospective analysis, none outperformed the COVIDHub ensemble, and one of each of the Deep Research and double-adapted models performed comparably to the ensemble. Their continued strong performance on the current prospective season provides direct evidence that ERA-generated models can generalize beyond the data on which they were originally trained and validated. As reported in Table 14, the Deep Research–originated models largely failed for RSV forecasting, with the majority performing worse than the CDC’s flat-line baseline—a naïve persistence forecast that projects the most recent observed value forward. The single-model adaptations likewise underperformed. Our best-performing RSV models ultimately came from unconstrained searches, suggesting that for extremely data-sparse pathogens the search benefits more from architectural freedom than from transferred methodological priors. Despite these constraints, both novel models included in the RSV pool outperformed the RSVHub-ensemble in this pathogen’s inaugural hub season. A distinctive feature of our approach is the explicit re-writing of existing CDC hub models, which enables a direct comparison between the prospective performance of ERA-generated adaptations and the original methods submitting to the same hub in real time. In general, there was not a clear correspondence between the score of the original model and the ERA-generated version of the model (Fig. 5). Every adaptation for influenza was judged to be a “Match”, “Partial Match” or “No Match” of the original model (see Methods). Adaptations judged as Matches or Partial Matches sometimes performed better and sometimes worse than the original model. In the comparison with models generated for influenza, all models that did not match the original implementation improved on the mean log WIS.
6
Prospective multi-pathogen disease forecasting using autonomous LLM-guided tree search
Part II: Controlled Retrospective Analysis Comparison of models generated using different LLMs We observed empirically that ERA, when operating without explicit methodological instructions, exhibits a systematic bias toward modeling approaches that are heavily represented in LLM training corpora—most notably gradient-boosted trees, random forests, and similar standard machine learning regressors. While such architectures can perform competitively in isolation, the historical strength of CDC hub ensembles derives precisely from their methodological diversity, which encompasses mechanistic compartmental models, classical time-series methods, and domain-specific statistical frameworks 9,11 . Reproducing that diversity with an automated system therefore requires not only instruction guidance but also a mechanism to enforce compliance, particularly for model classes that are structurally more difficult to implement and less well-represented in the coding knowledge available to the underlying LLM. To characterize both the limits and the enablers of instruction-following in ERA systematically, we designed a controlled retrospective ablation study. We selected four target methods from established influenza forecasting hubs, deliberately spanning a range of implementation difficulties: (1) a pure machine learning approach UMass-gbqr, which uses gradient boosting with feature engineering, (2) a package-dependent spatial time-series model UGA_flucast-INFLAenza requiring the pyinla library, (3) a highly structured mechanistic model Cornell_JHU-hierarchSIR, requiring a C++ SIR simulation bound to Python via pybind11, and (4) a highly specialized probabilistic framework NU-PGF_FLUH. We evaluated each configuration using three independent ERA-adaptations of each method to ensure robustness, with all experiments conducted under a standardized computational budget of 2,500 nodes or 2,500 hours of sandbox runtime, whichever limit was reached first. Three tiers of the Gemini model family were compared: Gemini 2.5 Flash, Gemini 3 Flash, and Gemini 3 Pro. The results show a clear gradient in instruction-following difficulty driven by how well each target method’s core computation patterns are represented in publicly available scientific software(Table 1). The UMass-gbqr method, which relies on pandas and lightgbm libraries ubiquitous in publicly available data science code, was implemented successfully across nearly all tiers and replicas. In contrast, the NU-PGF_FLUH method proved consistently challenging to implement across all model tiers. This method is highly technical, requiring formulating disease spread as a stochastic branching process and analytically characterizing its offspring distribution via Probability Generating Functions (PGFs) to approximate a Susceptible-Latent-Infectious-Recovered (SLIR) compartmental model. This is a highly specialized framework from mathematical epidemiology where little publicly available implementation code exists. Table 10 shows that none of the models were able to achieve full methodological fidelity for NU-PGF_FLUH, with all replicas assessed as partial matches by an external judge (see Methods). Strikingly, however, the inability to follow the specific method instructions did not prevent stronger models from achieving better predictive accuracy: Gemini 3 Pro attained a substantially lower mean WIS (∼140.09) compared to Gemini 2.5 Flash (∼236.06), suggesting that models generated by latergeneration LLMs discover more effective heuristic approximations even when the precise statistical framework remains out of reach. Similarly, the cross-language engineering demands of the Cornell_JHU-hierarchSIR method further exposed sharp capability differences between model tiers. This method requires a hierarchical SIR model implemented in C++ and integrated with Python via pybind11—a form of polyglot software engineering that demands both epidemiological domain knowledge and advanced systems programming proficiency. Gemini 2.5 Flash and Gemini 3 Flash largely failed to navigate this integration, producing mostly Partial Match outcomes across all replicas. Gemini 3 Pro, by contrast, achieved full fidelity success across all three replicas, discovering its best-performing solution notably
7
Prospective multi-pathogen disease forecasting using autonomous LLM-guided tree search
early within the search trajectory—reflecting the higher-tier model’s superior capacity to reason about cross-language API bindings and system-level constraints, which translated directly into a more efficient traversal of the solution space (Table 10). Comparison of approaches to enforce instruction following To characterize the fidelity of instruction following, we used two methodologically distinct judging mechanisms. The first is a generic in-loop LLM judge, which reads the natural-language method instruction and evaluates at each search node in ERA whether the generated code sufficiently follows the described approach. The resulting compliance score gates whether a candidate is passed to full validation or penalized with a worst-case WIS of 1,000 for not successfully following the given required instruction. The second uses an external post-hoc judge, deployed after the search to assess the final best-node models. Judging mechanisms are described in further detail in the Methods section. Unlike the in-loop judge, this external judge uses method-specific evaluation rubrics, developed with care and calibrated against assessments by expert epidemiologists, to assign graded fidelity outcomes (Match, Partial Match, No Match). These two judges serve different roles: the in-loop judge shapes the search, while the external judge provides the post-hoc ground-truth fidelity evaluation of unconstrained search results. Because ERA is a single-objective optimizer, simultaneously rewarding instruction fidelity and predictive accuracy requires combining these two objectives into a single signal. The judge-in-theloop implements one such strategy—a binary gate that zeroes out non-compliant branches—but its effectiveness depends critically on both the inherent difficulty of implementing the requested method and the degree of alignment between the generic in-loop judge and the calibrated external judge. Examining the four target methods from the previous section reveals four empirically distinct outcomes of this interaction (Table 11). For UMass-gbqr, where the method is straightforward to implement, the in-loop judge provided no meaningful benefit: fidelity and performance were already high without it. For UGA_flucast-INFLAenza, the in-loop judge successfully steered the search toward valid spatial time-series implementations, improving both fidelity (from one to three matches) and mean WIS (from ∼316.16 to ∼196.51), indicating strong alignment between the in-loop and external judges on what constitutes compliance for this method. For Cornell_JHU-hierarchSIR, the judge enforced some structural compliance—producing one fidelity match—but may have over-constrained the search space of the less capable Gemini 2.5 Flash, resulting in a marked degradation of predictive performance (mean WIS increasing from ∼179.09 to ∼411.97). For NU-PGF_FLUH, even with judgein-the-loop steering, ERA was unable to produce a method-faithful implementation—an outcome consistent with the inherent complexity of the method (described in the previous section). Taken together, these findings demonstrate that automated judging is most beneficial in an intermediate regime: the target method must be architecturally achievable by the underlying LLM, and the in-loop judge must be well-calibrated against expert standards. When the in-loop judge’s pass threshold (a score of 7 out of 10 in our design) is too permissive for a given method, it can accept implementations that the expert-calibrated external judge would classify as non-compliant. In these cases, the gating mechanism provides a false sense of fidelity without meaningfully steering the search. Tuning the judge prompt and the mechanism by which the judgment is incorporated into the search score is a critical design decision that needs further study. Internal optimization metrics shape the evolution of forecasting logic The choice of hill-climbing metric in an automated tree search fundamentally determines the kinds of models the search gravitates toward, how stable the discovery process is across independent runs, and how robustly the resulting models generalize out of sample. The evaluation harness for 8
Prospective multi-pathogen disease forecasting using autonomous LLM-guided tree search
ERA returns a single scalar score that serves as the sole reward signal guiding the LLM agent’s exploration; consequently, the mathematical properties of the chosen metric are directly encoded into the architecture and behavior of every model the search produces. To test how this affects model performance, we evaluated three candidate metrics drawn from the two principal families of proper scoring rules for probabilistic forecasting 22 . The Continuous Ranked Probability Score (CRPS) measures the integrated squared difference between the predictive cumulative distribution function (CDF) and the empirical CDF of the observation (the Weighted Interval Score (WIS) used by the CDC is a quantile-based approximation of CRPS 18 ). The Logarithmic Score (Log Score) is a classical proper scoring rule that evaluates the log-probability assigned to the observed outcome. Finally, the Log-scale CRPS is a variant in which both predictions and observations are log-transformed prior to scoring, combining the distributional breadth of CRPS with greater sensitivity to relative errors at low counts. Each metric was used to direct five independent replicas of 2,500-node searches, with all three metrics evaluated on the resulting best-node models (Fig. 6). The search trajectories produced by each metric differed markedly in both stability and their generalization behavior. Figure 6 tracks the cumulative best score across discovered nodes—the score of the best model found so far at each point in the search—averaged across five independent replicas (bold lines), with individual replica trajectories shown as faint lines. CRPS-driven searches exhibited steady, near-monotonic improvement across all replicas, converging smoothly toward competitive final scores. Log-scale CRPS searches were similarly stable, and additionally showed relatively close correspondence between validation and test performance throughout the search, reflecting consistent generalization across replicas. Log Score-driven searches were strikingly more volatile: individual replicas diverged substantially early in the search, with dramatic transient score spikes before partial recovery, and the spread across replicas was substantially wider than that observed under the other two metrics (Fig. 6, right column). This instability has a clear mathematical origin. Whereas CRPS and its approximations operate as distance-based metrics—penalizing forecast errors continuously as a function of magnitude—the Log Score measures the log-probability that the model’s predictive distribution assigns to the observed outcome. Because these models produce sample-based distributions, a candidate that assigns zero probability mass to the region containing the ground truth incurs an effectively infinite Log Score penalty. This produces an extremely noisy reward signal: many candidate architectures receive near-identical, near-maximal penalties regardless of how close to the ground truth their predictions were, making it difficult for the LLM agent to distinguish promising directions for further exploration from dead ends. The high between-replica variance visible in Figure 6 (right column) is consistent with this mechanism, as small random differences in early architectural choices lead to dramatically different search trajectories when the reward signal cannot reliably discriminate among candidate solutions. Distance-based metrics such as CRPS and Log-scale CRPS structurally avoid this pathology, providing smoother and more consistently navigable reward landscapes. A second important divergence concerns out-of-sample generalization. The persistent gap between validation (blue) and test (red) lines throughout Figure 6 reflects an intrinsic difference in epidemiological difficulty between the two seasons rather than a modeling failure; the diagnostic signal for overfitting is instead a rising test curve concurrent with a still-falling validation curve. Under CRPS optimization, this pattern is evident across multiple evaluation rows: in the Log CRPS panel (d), test performance worsens from ∼0.40 back toward ∼0.43 as validation continues declining to ∼0.34, and mild overfitting is also visible in panel a, where the test curve stabilizes and slightly upticks in the later stages of the search. Log-scale CRPS optimization (panels b, e, h) produced relatively tight co-movement between validation and test scores across all three evaluation metrics, and particularly for CRPS (panel b) and Log Score (panel h). The Log Score column exhibits a qualitatively distinct
9
Prospective multi-pathogen disease forecasting using autonomous LLM-guided tree search
pathology: rather than systematic overfitting, it is characterized by high between-replica variance, with individual runs diverging sharply, particularly in the first 500 nodes. However, after this, the mean cumulative score trajectory appears to stabilize, and does so at mean score values lower than or comparable to those of other optimization metrics. Individual test score trajectories also stabilize at lowest test score values across optimization metrics, with at least one trajectory below the minimum mean value (dashed red horizontal). Notably, these competitive aggregate scores mask qualitatively poor forecast behavior, as documented below: Log Score-optimized models produce distributions that are systematically too wide and biased downward. The choice of optimization metric in general shapes the qualitative character of the resulting forecasts in ways that extend beyond aggregate accuracy measures. As summarized in Table 2, models discovered under Log Score optimization exhibited substantially wider 50% prediction intervals (mean width 188.8, compared with 136.8 and 149.0 for CRPS and Log-scale CRPS respectively), a pronounced tendency to underestimate observed values (66.1% of predictions falling below the observed count, versus 54.5% and 55.7%), and a mean bias of −93.2 weekly hospital admissions— more than twice the underprediction of the other two conditions (−38.0 and −39.0, respectively). This pattern is consistent with a ‘timid’ forecasting profile: the infinite-penalty structure of the Log Score rewards models that hedge with wide distributions, systematically sacrificing responsiveness to epidemic surges in favor of distributional conservatism. CRPS and Log-scale CRPS produced broadly comparable forecast quality on these metrics, with similar MAE (∼96–98) and prediction interval widths. Taken together, these results identify Log-scale CRPS as the most reliable optimization target — it combines CRPS’s smooth reward landscape with closer validation/test agreement than CRPS itself, and avoids the wide, downward-biased forecasts induced by Log Score’s infinite-penalty structure
Discussion This work demonstrates that an LLM-guided tree search system can generate epidemiological forecasting models that compete with the state of the art across three respiratory pathogens in a fully prospective, leak-proof evaluation setting. The key finding is not that any single ERA-generated model outperforms all alternatives, but rather that the system can rapidly produce a large, methodologically diverse pool of competitive models whose ensemble matches or exceeds the gold-standard CDC hub ensembles—themselves aggregations of forecasts from dozens of expert teams 9,11,14 . Crucially, every model in this pool is a complete, executable Python program whose forecasting logic can be read, audited, and modified by any domain scientist. All model source code and weekly forecast submissions are publicly available under the Apache 2.0 open-source license1 —a property that distinguishes ERA from black-box AutoML 23 systems and aligns with the transparency requirements increasingly emphasized in public health forecasting. The practical implications of this automation extend far beyond the specific pathogens studied here. Currently, collaborative forecasting hubs rely on the sustained, labor-intensive participation of expert teams running and maintaining bespoke models week after week. This manual paradigm scales poorly to highly granular geographic resolutions (e.g., US counties), new geographic regions, or novel pathogens 2,8 . By automating model ideation, implementation, and evaluation, ERA democratizes access to expert-level forecasting, providing a concrete blueprint for deploying ensemble-grade capacity at vastly larger scales. Additionally, ERA could substantially lower the barrier to entry for hub participation, enabling rapid deployment of competitive forecasting capacity in settings that currently lack the specialized modeling workforce required. Moreover, this potentially allows expanding the forecast “hub” concept 8 to jurisdictions and granularities where it is harder to attract expert attention. 1 https://github.com/google-research/google-research/tree/master/epi_forecasts
10
Prospective multi-pathogen disease forecasting using autonomous LLM-guided tree search
The RSV results provide a direct test of a “cold start” scenario. With minimal historical data, no mature modeling literature to draw from, and a hub in its first season, ERA’s best-performing models came from unconstrained architectural search in which the system discovered effective approaches without domain-specific guidance. While the resulting models were necessarily simpler than those generated for influenza (both were gradient-boosted pipelines fed with cross-pathogen auxiliary inputs), their top-ranked performance among the other RSVHub submissions, ensemble and baseline model demonstrates that automated search can provide useful forecasting capacity even in data-scarce regimes. One limitation of the ERA-generated forecasts presented and evaluated in this work is that, while they achieve high accuracy scores relative to existing models, they suffer from similar systematic biases observed across multiple seasons of these forecasting efforts. Forecasts of seasonal influenza and COVID-19 have struggled to capture dynamics of rapid increase and decrease during seasonal epidemics 11,15,16 : often models underestimate the rate of increase as the epidemic is rising and the rate of decrease as it falls. Forecasts generated by ERA are no exception to this. However, this reflects not so much a limitation of the LLM model generation pipeline as it does a structural challenge with developing a set of scorable tasks for epidemic forecasting that incentivizes forecasts which can be used for improving situational awareness during seasonal outbreaks or used for decision-support. Developing evaluation rubrics that are aligned with public health decision-making is an active area of research 24–26 , and systems like ERA highlight the opportunities available for training models that have been optimized for specific contexts. A second limitation concerns the ensemble aggregation strategy itself. In this study, the submitted ensemble forecast was computed as a simple, equally weighted median across component models—an approach that, while robust and well-established in forecasting practice 14 , makes no attempt to learn from the relative strengths and weaknesses of its components. Yet the forecasting task is inherently high-dimensional: individual models may excel in different jurisdictions, at different forecast horizons, or under different epidemic conditions such as rapid growth versus post-peak decline. The optimal combination strategy is therefore not a single set of static weights but a contextdependent mapping—itself a challenging optimization problem 27,28 . Framing ensemble construction as a scorable task and applying ERA’s tree search to find effective adaptive weighting schemes or selection logic is a natural direction for future work. Additionally, we note that the current results are based on preliminary data as of 2026-05-02 and will be updated as final data for the season becomes available. However, the observed values in the spring months tend to be lower than in the fall and winter months and are not likely to substantially change to overall message of the present work, even if the specific scores or rankings change in a few places. A further limitation of the prospective evaluation is that not all models are assessed on identical sets of forecasting tasks. Models join hubs at different points in the season, occasionally skip submission weeks, or cover different subsets of jurisdictions, meaning that each pairwise comparison is computed over a different intersection of shared tasks. While the pairwise relative scoring framework we adopt (following established practice 18 ) mitigates this by computing score ratios only over mutually completed tasks and then averaging across all opponents, it does not fully eliminate the issue: a model that begins submitting after a particularly difficult or easy phase of the epidemic is effectively evaluated on a different forecasting problem than one that was active throughout. We partially address this by imposing an 80% task-coverage eligibility threshold, but acknowledge that strict like-for-like comparison across all models remains infeasible in any real-time collaborative forecasting setting. The controlled retrospective experiments in Part II reveal both the capabilities and the structural limitations of the current system. The bias toward gradient-boosted implementations—observed 11
Prospective multi-pathogen disease forecasting using autonomous LLM-guided tree search
even when the search was explicitly prompted with mechanistic or statistical instructions—reflects a fundamental asymmetry in the training distribution of the Google Gemini LLMs: standard software patterns are more reliably generated than highly specialized scientific frameworks. We showed that two complementary interventions can partially mitigate this bias: using more capable foundation models (Gemini 3 Pro achieved full fidelity on some of the more complex methods where Gemini Flash models could not), and deploying an LLM judge-in-the-loop to gate non-compliant implementations. However, for the most specialized methods (NU-PGF_FLUH) neither intervention was sufficient to achieve full methodological fidelity—an outcome reflecting the inherent complexity of the target method. An important question remains open: what is the role of expertise and expert judgment in guiding model development in a real-world challenge such as epidemiological forecasting? The Google team that drove this project included several scientists, analysts, and engineers (without explicit infectious disease or epidemiological training) and one expert consultant from academia who has participated in epidemiological forecasting challenges for around 10 years. While the tree search experiments were orchestrated, compiled and analyzed by the Google engineers, the expert provided input on which existing models to prioritize exploring, and provided some judgments on how faithful the ERA-generated re-implementations of models were. It is challenging to quantify the contribution of domain-specific expertise to these forecasting results. However, in one case, the hybrid of LANLDBM and LANL-Inferno, the expert identified two models that had been high-performing models several years ago but were no longer being submitted. ERA implemented a version of this model that captured important elements of the original model, and this adaptation was highly performant in the prospective evaluation. Therefore, we suggest there is some anecdotal evidence that even in the presence of a highly automated LLM-driven model-generation architecture, there remains a place for human expertise to highlight opportunities for model exploration. This work demonstrates the large potential for LLMs in agentic harnesses to be used in real-world, real-time, public good settings. By democratizing access to expert-level models over a wide array of methodological approaches, systems like ERA can provide humans in a variety of fields with model code that has been optimized for use specific scientific settings, lowering the barrier to entry and development in resource-limited cases.
Acknowledgements We are grateful to our colleagues in Google Research and Google DeepMind for the environment in which to do this work. We are grateful for our coauthors in the original ERA paper 1 , and in particular would like to thank Matthew Abraham, Erica Brand, Marc Coram, Lizzie Dorfman for help and important discussions, as well as John Platt, Yossi Matias and James Maniyka for their support and encouragement.
Code Availability ERA-generated model source codes for models that were submitted prospectively and their weekly submissions are publicly available at https://github.com/google-research/google-research/tree/ master/epi_forecasts.
12
Prospective multi-pathogen disease forecasting using autonomous LLM-guided tree search
Data Availability All surveillance data used in this study are publicly available. Hospital admission data were obtained from the CDC’s National Healthcare Safety Network (NHSN) Weekly Hospital Respiratory Data reporting system 29 . Outpatient syndromic surveillance data were obtained from the CDC’s ILINet system. Forecast submissions from all participating teams are publicly available through the CDC FluSight, COVIDHub, and RSVHub repositories.
Competing Interests S.M. carried out this work as part of a student researchership at Google Research. M.P.B. holds appointments at both Harvard University and Google Research. M.P. is an employee of Google DeepMind. B.P.W. and Z.S. are employees of Google Research. N.G.R. serves as a faculty member at the University of Massachusetts Amherst and was paid as a scientific consultant by Google to provide expertise for and advice on the experiments presented in this paper. His effort on this project was supported by Google, and aligned with UMass policies on external consulting. His engagement with Google was disclosed to and approved by the UMass Research Compliance office.
13
Prospective multi-pathogen disease forecasting using autonomous LLM-guided tree search
Methods 1. The Epidemiological Forecasting Task The forecasting tasks addressed in this study focus on three major seasonal respiratory pathogens in the United States: influenza, COVID-19, and RSV. Unlike retrospective "hindcasting," this work was implemented as a real-time prospective study, requiring models to generate predictions before the ground-truth observations were recorded. This approach subjects the forecasting process to real-time surveillance conditions, including reporting delays and data revisions typical of live national healthcare systems. Epidemiological datasets are frequently revised and back-filled over time, meaning that preliminary data available during a live outbreak is often less complete than the finalized historical record. While our prospective experiments were evaluated in real-time using these early data snapshots, the retrospective analysis used the updated versions of the data. This distinction is not a threat to the study’s validity because retrospective models are compared solely against each other in internal benchmarks rather than being measured against historical real-time hub submissions. The primary data source for inputs and targets was the CDC’s Weekly Hospital Respiratory Data (HRD) from the National Healthcare Safety Network (NHSN) 29 . As shown in Figure 1a, the study follows a rolling-origin submission scheme where new data are released weekly, triggering a corresponding forecast submission. During the study period, this resulted in a total of 22 weekly submissions for influenza, 19 for COVID-19, and 18 for RSV, of which 17 were included in the analysis (the first submission was excluded due to a technical operational error in the submission pipeline) (Fig 1c). In our pipeline, preprocessing consisted of extracting jurisdiction-level weekly admission counts, cleaning non-numeric entries, and linearly interpolating any missing weeks to produce a complete time series grid across all 52 jurisdictions. The specific input features provided to each pathogen’s models differed based on data availability and relevance: influenza models received weekly incident influenza hospitalizations supplemented with a longer history of outpatient syndromic surveillance data from the CDC’s ILINet system; COVID-19 models received weekly incident COVID-19 hospitalizations; and RSV models received weekly hospitalizations for all three pathogens along with ILINet data, reflecting the cross-pathogen auxiliary input strategy described in the Results. Beyond this standardized preprocessing, each generated model was free to implement its own additional feature engineering within the ERA sandbox. The scale of the forecasting task is high-dimensional, requiring full geographic coverage of the United States. Models provide predictions for 52 jurisdictions, including all 50 states, the District of Columbia, and Puerto Rico. For each jurisdiction, forecasts must cover a 4-week time horizon (Fig. 1c). Additionally, to characterize the uncertainty of future disease trajectories, forecasts are produced in a probabilistic format rather than as single-point estimates. As illustrated in Figure 1b, for each target week, the models generate a predictive distribution represented by 23 specific quantiles (ranging from 0.01 to 0.99). This allows the system to represent the full range of potential outcomes, from the expected median trajectory to the low-probability tails of the distribution that represent extreme surge events. Together, these requirements result in the generation of 4,784 unique prediction points per infection type (calculated as 52 jurisdictions × 4 weeks × 23 quantiles) every week. The temporal data for each pathogen were divided into three non-overlapping splits to ensure rigorous evaluation: Validation, Retrospective Test, and Prospective Held-out (Fig. 1d). The Validation block was used during search for hill-climbing (Section 4.2); the Retrospective Test block was used post-search for node selection; and the Prospective Held-out period represents the live evaluation
14
Prospective multi-pathogen disease forecasting using autonomous LLM-guided tree search
phase, where model performance was recorded in real-time against data that did not exist during the system’s development phase. The Prospective Held-out period represents the live evaluation phase, where model performance was recorded in real-time against data that did not exist during the system’s development phase. 2. Evaluation Metrics The primary metric for evaluating accuracy across all phases is the Weighted Interval Score (WIS). The WIS is a proper scoring rule that penalizes both lack of precision (wide intervals) and poor calibration (ground truth falling outside the predicted quantiles). Small WIS values imply predictions are both sharp and well-calibrated. Formally, let 𝐹 denote a predictive distribution of the quantity of interest 𝑌 with median 𝑚 and central prediction intervals at levels 𝛼1 , . . . , 𝛼 𝐾 , where 𝑙𝛼𝑘 and 𝑢𝛼𝑘 are the 𝛼𝑘 /2 and 1 − 𝛼𝑘 /2 quantiles of 𝐹 . The WIS for an observed value 𝑦 is defined as 18 : ! 𝐾 ∑︁ 1 WIS( 𝐹, 𝑦 ) = 𝑤0 | 𝑦 − 𝑚 | + 𝑤𝑘 · IS𝛼𝑘 ( 𝐹, 𝑦 ) , (1) 𝐾 + 0.5 𝑘=1 where 𝑤0 = 0.5, 𝑤𝑘 = 𝛼𝑘 /2, and the interval score for level 𝛼𝑘 is IS𝛼𝑘 ( 𝐹, 𝑦 ) = (𝑢𝛼𝑘 − 𝑙 𝛼𝑘 ) +
2 𝛼𝑘
( 𝑙 𝛼𝑘 − 𝑦 ) 1( 𝑦 < 𝑙 𝛼𝑘 ) +
2 𝛼𝑘
( 𝑦 − 𝑢𝛼𝑘 ) 1( 𝑦 > 𝑢𝛼𝑘 ) .
In this study, forecasts are represented by 23 quantiles corresponding to 𝐾 = 11 central prediction intervals plus the median, yielding 𝛼𝑘 ∈ {0.02, 0.05, 0.10, . . . , 0.90, 0.98}. The log WIS variant used for official CDC leaderboard rankings applies the same WIS formula (1) to log-transformed predictions and observations, i.e., inputs are transformed following ˜𝑦 = log(1 + 𝑦 ), where log is the natural logarithm and the offset of one accommodates possible zero counts. Concretely, let 𝐹˜ denote the predictive distribution of log(1 + 𝑌 ) induced by 𝐹 , or, equivalently, the distribution with quantiles log(1 + 𝑙 𝛼𝑘 ) , log(1 + 𝑢𝛼𝑘 ), and median log(1 + 𝑚). Then, log WIS = WIS( 𝐹˜, ˜𝑦 )). Because not all models submit predictions for the same set of tasks—some join later in the season, omit certain reference dates, or cover fewer jurisdictions or prediction horizons—direct comparison of mean scores can be misleading. We therefore report pairwise relative WIS 18 as the primary model comparison metric in aggregated settings. For each pair of models ( 𝑖, 𝑗), the ratio of their mean scores is computed over the intersection of tasks both models submitted; the relative score for model 𝑖 is then the geometric mean of these pairwise ratios across all opponents 𝑗 ≠ 𝑖. We rescale scores so the relevant CDC hub ensemble reports a score of one, and scores below this indicate better-than-ensemble performance. We compute this metric for WIS and log WIS, and report both in Supplementary Tables 4, 6, and 8. Pairwise relative log WIS serves as the primary ranking metric throughout. To complement these aggregate scores with a measure of consistency, we report standardized ranks 11 . For each individual forecasting task (a unique combination of jurisdiction, horizon, and reference date), all submitting models are ranked by log WIS and the ranks linearly rescaled to [0, 1], where 1 denotes the best-performing model on that task and 0 the worst. The resulting distribution of standardized ranks across tasks summarizes how frequently a model performs near the top or bottom of the field, capturing variability in performance that aggregate mean scores alone would obscure. 3. The ERA Framework ERA is an autonomous system that automates the development of empirical software for scorable 15
Prospective multi-pathogen disease forecasting using autonomous LLM-guided tree search
scientific tasks 1 . It searches over the space of code implementations using a tree search algorithm: each node in the tree corresponds to a concrete implementation–including specific combinations of feature engineering, data preprocessing, and model architecture–and is scored against a predefined objective function. Node selection for tree expansion balances exploitation of the highest-scoring candidates with exploration of less-visited ones via a PUCT (Predictor + Upper Confidence bound applied to Trees) algorithm. The search is guided by a Large Language Model (LLM)—specifically the Gemini family of models—acting as an agent that generates, tests, and refines Python code. Starting from a natural language problem statement, the agent proposes initial architectures and iteratively improves them based on execution feedback and performance scores. This branching structure allows the system to explore multiple modeling paradigms simultaneously, learning from the successes and failures of previous nodes to discover increasingly effective solutions. To manage computational resources and ensure system stability, ERA search trajectories are bounded by dual constraints: a maximum number of explored nodes and a cumulative execution time limit within the sandbox environment. To maintain scientific rigor and prevent data leakage, ERA operates within a secure, immutable evaluation harness. This harness contains the datasets and the evaluation logic but is entirely hidden from the LLM agent. The agent can submit code to the harness and receive performance metrics and error logs, but it cannot modify the evaluation protocol or access ground-truth data from the held-out periods. This ‘blind’ evaluation ensures that discovered models are robust and prevents the system from ‘hacking’ evaluation criteria. 4. ERA Configuration for Epidemiological Forecasting This section describes three task-specific design choices made to apply ERA to multi-pathogen forecasting. The first is a structured prompting strategy across multiple search modes, and the second is a two-stage model selection scheme that separates search-time scoring from post-search node selection. Both are part of the pipeline used to generate models for prospective submission to the CDC hubs. The third is an LLM-as-judge gate on methodological fidelity, used in a retrospective analysis of ERA’s ability to follow given instructions. 4.1. Prompt Sources and Search Modes Each independent ERA search ran for up to 2, 500 generated nodes or until a cumulative sandbox execution time limit was reached, whichever came first. Within these bounds, the LLM agent generated candidate code that was executed inside a sealed evaluation harness (Section 3); only the resulting score and execution logs were returned to the agent, with all data and evaluation logic immutable from the agent’s perspective. For each pathogen and forecast origin, searches were initialized under one of three prompt regimes. Single Model Adaptations, provide the agent with a description of one existing methodology—drawn from a forecasting hub (e.g., CDC FluSight) or a published paper—and instruct it to implement and refine that method. Double Model Adaptations provide descriptions of two distinct existing methods and instruct the agent to combine their logic into a single hybrid combining strengths of both. Novel Model Search provides only the task description, leaving the search unconstrained over architectures. The full prompts used to seed Single and Double Model Adaptation runs are listed in Supplementary Section D, and the prompt outline for the Novel Model Search runs in Supplementary Section C.
16
Prospective multi-pathogen disease forecasting using autonomous LLM-guided tree search
4.2. Two-stage model selection Hill-climbing metric and moving-window aggregation. During search, candidates were scored using WIS on the validation block of the relevant pathogen-season split (Fig. 1d, blue). To avoid rewarding candidates that fit narrowly to a single forecast origin, each candidate’s validation WIS was aggregated across multiple rolling origins within the validation block. This biases the search toward models that perform consistently across epidemic phases (growth, peak, and decline) rather than to any single time slice. The aggregated validation WIS drove node selection within the PUCT tree search (Section 3). Post-search node selection. A single search produces a tree of candidates whose aggregated validation WIS scores are known by construction. To choose one representative model per search while protecting against overfitting to the validation block, we evaluated every node post-hoc on a held-out retrospective test block (Fig. 1d, red)—a period never seen by the agent during the search—and selected the node minimizing Selection Score = Validation WIS + 2 × Retrospective Test WIS. The 2× weighting on the retrospective-test term places greater emphasis on out-of-search generalization than on in-search fit. 5. Real-time Prospective Implementation 5.1. Large-Scale Candidate Model Generation To construct a high-capacity library of forecasting candidates, we executed a wide array of independent ERA search experiments prior to the start of the prospective season. These experiments spanned several strategic categories: adaptations of established expert models from literature and forecasting hubs, Double Model Adaptations merging distinct methodologies, and unconstrained searches driven by LLM exploratory hypotheses and "Deep Research" agents (Supplementary Fig. 17). This comprehensive approach ensured that the system began the prospective phase with a broad and methodologically diverse library of models for each pathogen (Tables 12, 13, 14). 5.2. Pathogen-Specific Implementation Logic The application of ERA varied across the three target pathogens due to differences in data availability and hub maturity. Influenza benefited from a mature forecasting ecosystem (CDC FluSight), allowing for a large pool of candidate models. For COVID-19, in addition to new ERA-generated models, five models were carried over from previous retrospective research and fine-tuned on the larger validation dataset available for this prospective season to ensure calibration across evolving variants. In contrast, RSV presented a "cold start" challenge due to the lack of mature literature and established forecasting hubs. For RSV, ERA relied on unconstrained architectural searches and the inclusion of cross-pathogen indicators (COVID-19 and Influenza trends) to identify predictive signals, resulting in a smaller but specialized pool of models specifically tailored for this pathogen. RSV presented a qualitatively different challenge. Historical surveillance data for this pathogen are extremely limited, the published forecasting literature is sparse, and the short available hospitalization record constrained the size of meaningful validation and test splits. To compensate, we provided ERA models with concurrent COVID-19 and influenza hospitalization time series as auxiliary inputs to harness for generalizable cross-pathogen patterns. We pursued three complementary strategies to overcome data scarcity. First, we used Gemini Deep Research to survey the epidemiology, econometrics, and deep learning literature for methods suited to forecasting with sparse target data and richer 17
Prospective multi-pathogen disease forecasting using autonomous LLM-guided tree search
auxiliary signals. This survey produced seven candidate architectures, each addressing data sparsity through a distinct mechanism spanning Bayesian transfer learning, mechanistic-neural hybrids, multi-task deep learning, and domain adaptation. Each strategy was translated into a step-by-step implementation prompt and provided to ERA for automated search. Second, we attempted singlemodel adaptations of methods that had begun submitting to the RSV hub earlier in its inaugural season. Third, we ran unconstrained ERA searches without specific methodological instructions. 5.3. Tiered Selection and Ensemble Aggregation Model selection followed a rigorous tiered hierarchy. First, the optimal "best node" from each independent search tree was identified using the compound selection score (𝑉 𝑎𝑙𝑖𝑑𝑎𝑡𝑖𝑜𝑛 𝑊 𝐼𝑆 + 2 × 𝑇𝑒𝑠𝑡 𝑊 𝐼𝑆). From this pool of winners, a subset was curated for inclusion in an internal hub—hosted publicly on GitHub—based on a combination of their performance scores and methodological diversity. Finally, a smaller, high-performance subset was selected to form the official "Google-SAI" submission ensemble for each pathogen. The final forecasts were generated using a median ensemble approach, where the prediction for each jurisdiction and horizon was calculated as a simple average of the 23 quantiles across all participating component models. 5.4. Operational Protocol and Real-time Adaptation The prospective evaluation phase consisted of approximately 20 weeks of live submissions to official and internal forecasting hubs. Throughout this period, the system ingested real-time surveillance data—including weekly hospital admissions and syndromic ILINet indicators—as they were released. This live deployment required models to generate forecasts using the preliminary data versions available at the time of submission, providing a true test of utility in an active public health response. For Influenza forecasts, a strategy shift occurred around January 24, 2026; upon learning that the CDC utilizes the Log WIS metric for official evaluations, we updated our selection process to prioritize models based on Log WIS performance to better align with the official evaluation framework. We note that the first RSVHub submission was excluded from evaluation due to a technical operational error; all RSV performance results therefore reflect the remaining 17 submission weeks. 6. Retrospective Experimental Design 6.1. Evaluation of Search Agent Fidelity and Instruction Following To evaluate how the choice of the underlying AI "brain" affects model discovery, we conducted a systematic ablation study comparing Gemini 2.5 Flash, Gemini 3.0 Flash, and Gemini 3.0 Pro. Success was measured across two dimensions: the predictive accuracy (WIS) of the resulting models and their "Judgment" score (Match, Partial Match, or No Match), representing the agent’s ability to satisfy structural constraints, such as the specific implementation of hierarchical SIR architectures. We further validated the Automated Judge by executing parallel search experiments—with and without the judge-in-the-loop—to quantify the success rate in implementing non-straightforward compartmental logic. This phase also served to calibrate the LLM-Judge against expert human audits. To ensure a fair comparison between different experimental configurations in the retrospective study, we enforced a standardized computational budget for each search trajectory. Each experiment was terminated upon reaching either a total of 2,500 explored nodes or a cumulative sandbox runtime of 2,500 hours, whichever criterion was met first. This design ensures that experiments producing computationally intensive models—which require longer training or execution times—are evaluated against an equivalent resource expenditure as experiments producing simpler, faster-executing
18
Prospective multi-pathogen disease forecasting using autonomous LLM-guided tree search
architectures. This standardized bounding prevents performance gains from being confounded by disparate computational allocations. The optimization of ERA search trajectories requires a strategic navigation of the inherent tradeoffs between node depth and cumulative sandbox execution time. While the framework typically bounds searches by both parameters, the effective limit is dictated by the architectural complexity of the candidate models rather than their specific scientific paradigm. High-capacity architectures—such as the feature-intensive machine learning models or mechanistic simulations—tend to be timebound, frequently exhausting the sandbox budget while exploring significantly fewer nodes than computationally lighter counterparts. Consequently, users must calibrate search budgets to account for the fact that per-node computational costs for training and cross-validation directly constrain the total depth of the discovery path. 6.2. Analysis of Optimization Metrics and Search Stability We investigated how the choice of the internal “hill-climbing” metric dictates the stability and evolution of the forecasting logic. Using a controlled retrospective setup, we compared the effects of optimizing for three candidate metrics drawn from the two principal families of proper scoring rules for probabilistic forecasting 22 . The Continuous Ranked Probability Score (CRPS) measures the integrated squared distance between the predictive cumulative distribution 𝐹 and the step function at the observed value 𝑦 : ∫ ∞ 2 CRPS( 𝐹, 𝑦 ) = 𝐹 ( 𝑥 ) − 1( 𝑥 ≥ 𝑦 ) 𝑑𝑥. (2) −∞
CRPS generalizes mean absolute error to full distributions and, like WIS (which is itself an intervalbased approximation to CRPS), is a distance-based metric that penalizes predictions continuously as a function of their deviation from the observation. The Log-scale CRPS variant applies the same formula after log-transforming both predictions and observations: Log CRPS = CRPS( 𝐹log , log(1 + 𝑦 ), placing greater weight on relative errors at low counts. The Logarithmic Score (Log Score) evaluates the negative log-density assigned to the observed outcome: Log Score( 𝐹, 𝑦 ) = − log 𝑓 ( 𝑦 ) ,
(3)
where 𝑓 is the predictive density. Unlike CRPS, the Log Score assigns an infinite penalty when the model places zero probability mass at the observation, producing a sharply discontinuous reward landscape. As illustrated in Figure 6, this analysis allowed us to track the “Cumulative Best” progress over 2,500 nodes, revealing how specific metrics steer the search process toward either conservative, stable models or more aggressive, high-risk trajectories. To enable the calculation of Log Scores and other non-standard metrics during these optimization experiments, we modified the predictive output format. Rather than generating the standard 23 quantiles, models in this subset were required to output 1,000 predictive samples per target. This sample-based format allowed for the utilization of the scoringrules library, which provides robust implementations for Log Score, CRPS, and Log-CRPS calculations. By maintaining this consistent sample-based framework, we ensured that the comparison between objective functions was mathematically rigorous. To guarantee data integrity, all retrospective experiments were conducted using an immutable evaluation harness that isolated the LLM agent from the ground-truth data.
19
Prospective multi-pathogen disease forecasting using autonomous LLM-guided tree search
Main Text Figures and Tables
Figure 1 | Overview of the epidemiological forecasting task and evaluation framework. a, Rolling submission scheme. Historical surveillance data (black line) are available up to the forecast origin (vertical dashed lines). At each weekly submission point, forecasts are generated for four subsequent target weeks (horizontal arrows). b, Probabilistic forecast for an example submission week. Taking the end of the observed data as the origin, the model produces a full predictive distribution via 23 quantiles for each of the four target horizons. The dotted line indicates the median trajectory. c, Task dimensions and scale for each pathogen: 52 US jurisdictions, 4-week horizon, 23-quantile distributions, with weekly updates. d, Sample time series of weekly hospital admissions for influenza (top), COVID-19 (middle), and RSV (bottom) for California, with non-overlapping temporal splits for validation (blue), retrospective testing (red), and prospective evaluation (green).
20
Prospective multi-pathogen disease forecasting using autonomous LLM-guided tree search
Figure 2 | The ERA System Architecture. The ERA workflow automates the translation of natural language scientific hypotheses into optimized executable code. Inputs: ERA ingests a natural language problem statement, datasets, and an evaluation harness. Processing: An LLM-driven tree search explores a high-dimensional solution space, generating and refining Python implementations of epidemiological models. Output: ERA identifies the best-scoring solution—a mathematically optimized model that minimizes the defined loss function while adhering to architectural constraints.
21
Prospective multi-pathogen disease forecasting using autonomous LLM-guided tree search
Figure 3 | Prospective Performance of ERA-Generated Ensembles on CDC Leaderboards. A Weekly observed hospital admissions (black) and Google-SAI ensemble forecast submissions for California, with median and 50% prediction intervals (shaded), for influenza, COVID-19, and RSV. B Pairwise relative log WIS ranking of the Google-SAI ensemble (dark blue) among eligible CDC hub submissions (lower is better). For each pair of models ( 𝑖, 𝑗), the ratio of their mean log WIS scores is computed over tasks both models submitted; the relative score for model 𝑖 is the geometric mean of these pairwise ratios across all other models, rescaled so the CDC hub ensemble equals one (dashed horizontal line). Values < 1 indicate better performance than the CDC ensemble. Eligible models were those submitting scorable predictions for at least 80% of covered reference dates, horizons, and jurisdictions. The Google-SAI ensembles achieved top-tier placement among eligible models across all three pathogens.
22
Prospective multi-pathogen disease forecasting using autonomous LLM-guided tree search
Figure 4 | Component Models Performance. The ERA pipeline generated a broad and diverse pool of models during the prospective season. Each bar represents the season-average pairwise relative log WIS of Google internal component models and CDC-submitted models. For each pair of models ( 𝑖, 𝑗), the ratio of their mean log WIS scores is computed over tasks both models submitted; the relative score for model 𝑖 is then the geometric mean of these pairwise ratios across all other models, rescaled so the CDC hub ensemble equals one (dashed horizontal line). Values < 1 indicate better performance than the CDC ensemble. CDC-submitted models are restricted to the top 20 best-performing eligible models (≥ 80% task coverage), plus the CDC hub ensembles, hub baseline, and the Google SAI ensemble submission. The ERA-generated component models (light blue) are hatched by methods origin: Google Internal Single Model Adaptations from hubs and relevant publications (solid blue), Google Internal Double Model Adaptations (diagonal-hatched blue), Google Internal Novel Models (cross-hatched blue). Additionally, CDC-submitted models (solid red) and CDC hub ensemble and baseline (diagonal-hatched red) are included. Asterisks denote components of the Google SAI ensemble submission (dark blue). Note that relative log WIS scores may be different than other results (e.g., in Figure 3 because they are computed on a different set of models, including internal Google models. top Influenza: 34 models generated, with two individual ERA-generated models and the Google SAI ensemble outperforming the FluSight-ensemble. bottom left COVID-19: 17 models generated with the Google_SAI-Ensemble outperforming the COVIDHub-ensemble. bottom right RSV: two novel models generated, both outperforming the RSVHub-ensemble in the hub’s inaugural season.
23
Prospective multi-pathogen disease forecasting using autonomous LLM-guided tree search
Figure 5 | Methodological fidelity and prospective performance of ERA-generated adaptations for influenza forecasting. Paired comparison of prospective WIS (lower is better) on a logarithmic scale between the target method from the CDC Hub (deep charcoal circles) and its corresponding ERA-generated adaptation (diamonds). The ERA adaptations are color-coded based on their structural scientific fidelity as assessed by an external post-hoc judge: green diamonds indicate a full Match to the target methodology; yellow diamonds indicate a Partial Match with substantive algorithmic deviations; and red diamonds indicate a No Match where the search abandoned the requested framework. Sample sizes (𝑛) are integrated directly into the model labels on the y-axis and denote the intersection of forecasting tasks—defined by unique combinations of jurisdiction, horizon, and reference date—mutually completed by both the target method and the ERA adaptation, serving as the shared baseline over which the pairwise mean WIS was calculated.
24
Prospective multi-pathogen disease forecasting using autonomous LLM-guided tree search
Table 1 | Summary of LLM search agent performance and instruction-following fidelity, for a selection of models generated retrospectively. The Judge column indicates whether the automated judge-in-theloop was active during the search (Yes) or whether the search relied solely on the forecasting objective (No). Outcomes are reported as counts of Match (M), Partial Match (PM), and No Match (NM) across runs. Values are reported as Mean ± SD where applicable. Asterisks (*) indicate where data was incomplete or experiments are still running. Note that WIS values in this table are reported on the natural scale and are not directly comparable to log WIS scores reported in the main prospective results. Method
Gemini Version Judge Outcomes (M/PM/NM) WIS (Mean ± SD)
UMass-gbqr
2.5 Flash 2.5 Flash 3 Flash 3 Pro
No Yes No No
3/0/0 3/0/0 2/1/0 2/1/0
291.38 ± 208.43 167.62 ± 17.72 166.95 ± 8.91 155.08 ± 5.47
UGA_flucast-INFLAenza
2.5 Flash 2.5 Flash 3 Flash 3 Pro
No Yes No No
1/2/0 3/0/0 1/2/0 0/3/0
316.16 ± 110.19 196.51 ± 37.95 176.55 ± 4.53 137.78 ± 18.58
2.5 Flash 2.5 Flash Cornell_JHU-hierarchSIR 3 Flash 3 Pro
No Yes No No
0/3/0 1/1/1 1/2/0 3/0/0
179.09 ± 11.53 411.97 ± 207.36 158.60 ± 15.33 140.62 ± 11.93
2.5 Flash 2.5 Flash 3 Flash 3 Pro
No Yes No No
0/3/0 0/3/0 0/3/0 0/3/0
236.06 ± 83.69 205.11 ± 82.40 153.97 ± 5.22 140.09 ± 12.76
NU-PGF_FLUH
Table 2 | Forecast quality statistics for models discovered under each optimization metric, evaluated on the retrospective influenza test set. Values are means across five independent search replicas. Mean Bias: average difference between predicted median and observed admissions (negative = underprediction). Prop. Under/Over: proportion of predictions where the median falls below or above the observed value. Mean 50% CI Width: average width of the 50% prediction interval (lower = sharper). MAE: mean absolute error of the predicted median. Optimization Metric Mean Bias Prop. Under Prop. Over Mean 50% CI Width
MAE
CRPS Log-scale CRPS Log Score
96.3 97.6 135.7
−37.96 −39.01 −93.16
0.545 0.557 0.661
0.429 0.415 0.324
136.8 149.0 188.8
25
Prospective multi-pathogen disease forecasting using autonomous LLM-guided tree search
Figure 6 | Influence of Hill-Climbing Metrics on Search Stability and Generalization. Each column represents one set of ERA experiments, with the objective function set to be optimizing for one of three metrics: CRPS, log-scale CRPS and Log Score. Five experiments were conducted for each metric (narrow lines) and, for each experiment, the resulting forecasts were scored using each of the three metrics on both an in-sample validation set (blue lines) and out-of-sample test set (red lines). Each line tracks the “cumulative best" score of the Tree Search across 2,500 nodes (x-axis). The mean value of five experiments are highlighted (wider lines). The dashed horizontal lines represent the best mean score achieved across all optimization experiments for that specific metric. The mean validation score lines are monotonically decreasing when the optimized score and the evaluation metric agree because we track the cumulative best-scoring node. However, the test score of that same node need not decrease monotonically, and a rising test curve concurrent with a still-falling validation curve is a diagnostic signal of overfitting to the validation period. The persistent vertical gap between validation and test lines reflects a difference in epidemiological difficulty between the two seasons rather than a modeling failure.
26
Prospective multi-pathogen disease forecasting using autonomous LLM-guided tree search
Supplementary Figures & Tables A. Prospective Model Performance Analysis This section contains additional analysis figures for the models submitting prospectively to the three CDC hubs and for component models submitted to the Google Research internal hub. All analyses is restricted to 52 jurisdictions (50 states, Washington D.C., and Puerto Rico; national-level forecasts are excluded) and to forecast horizons 0–3 weeks ahead (horizon -1 excluded). Table 3 | Validation and Test Periods for Each Pathogen. Pathogen
Validation Period
Test Period
Influenza
2023-10-13 to 2024-05-15 & 2024-10-15 to 2025-05-15 2023-06-01 to 2024-06-01 2024-10-01 to 2025-07-01
2022-10-15 to 2023-05-15
COVID-19 RSV
2024-06-01 to 2025-09-01 2025-07-07 to 2025-11-07
A.1. Influenza Here, we present individual analysis for influenza forecasts, both for our component models and for those models submitting to the FluSight Forecast Hub. Models submitting influenza forecasts are analyzed over reference dates from 22 November 2025 to 2 May 2026 (24 reference dates, covering a total of 4, 680 forecasting tasks). Analyses are restricted to models submitting for ≥ 80% of the required task space, i.e. ≥ 3, 744 (unless otherwise specified). Under this criterion, 43 hub-submitting models are eligible for analysis out of a total of 57 models which submitted at least one task prediction to the FluSight Forecast Hub this season. Table 4 | Summary performance of ERA-generated component models and eligible models submitting to the FluSight Forecast Hub. Models are ranked by pairwise relative log WIS (ascending) across submitted tasks, where a value of 1.0 indicates performance equal to the hub ensemble; values below 1.0 (green) indicate better performance and values above 1.0 (red) indicate worse performance than the ensemble. For each pair of models, score ratios are calculated over the intersection of tasks both models submitted, then geometrically averaged across all opponents, such that the metric is not affected by differences in the set of tasks each model chose to forecast 18 . The Mean log WIS and Mean WIS columns are shaded from green (lowest, best) to red (highest, worst) across models. Due to pairwise relative WIS being computed on prediction task overlap between pairs of models, which varies across all pairs, and the geometric mean being employed in the calculation, one model can have a lower mean log WIS than another model but a higher pairwise relative log WIS.
Model
n tasks
Pairwise Pairwise Mean Rel. log Rel. log WIS WIS WIS
Mean WIS
Type
Google_SAI-FluEns
4680
0.917
0.3002
0.941
71.37
CDC
*G-LANL_DBM x LANL_Inferno
4680
0.931
0.2979
0.797
58.23
Hybrid
CMU-TimeSeries
4524
0.944
0.3090
0.872
70.08
CDC
FluSight-trained_mean
4680
0.947
0.3124
0.966
74.39
CDC
UGA_flucast-INFLAenza
4264
0.954
0.3191
0.910
63.30
CDC
FluSight-HJudge_ensemble
4680
0.955
0.3164
0.954
73.48
CDC 27
Prospective multi-pathogen disease forecasting using autonomous LLM-guided tree search
Model
n tasks
Pairwise Pairwise Mean Rel. log Rel. log WIS WIS WIS
Mean WIS
Type
OHT_JHU-nbxd
4680
0.958
0.3028
1.059
79.25
CDC
*G-Cornell_JHU-hierarchSIR
4472
0.982
0.3101
0.930
69.55
Adapted
UMass-flusion
4472
0.984
0.3200
0.959
71.03
CDC
NAU-vulPES
4680
0.990
0.3271
0.950
72.31
CDC
FluSight-trained_med
4680
0.992
0.3292
0.969
74.49
CDC
FluSight-ensemble
4680
1.000
0.3361
1.000
78.20
CDC ensemble
G-PSI-PROF_MOA
4056
1.008
0.3034
0.985
71.98
Adapted
G-UMass-gbqr
4680
1.012
0.3309
1.004
76.67
Adapted
*G-NU-PGF_FLUH v3
2808
1.014
0.2756
0.959
30.14
Adapted
UVAFluX-FS_OptimWISE
3848
1.015
0.2996
1.070
73.88
CDC
*G-LANL_DBM x UMass_Flusion
4056
1.019
0.3058
1.003
72.26
Hybrid
FluSight-lop_norm
4680
1.022
0.3388
0.978
75.23
CDC
G-LANL_DBM v2
2808
1.039
0.2818
1.324
41.32
Adapted
*G-CU_SIRS x CMU_climate_baseline
4472
1.041
0.3395
1.201
94.89
Hybrid
MIGHTE-Joint
3848
1.054
0.3098
1.098
74.70
CDC
G-NU-PGF_FLUH
4472
1.055
0.3357
1.029
76.83
Adapted
MIGHTE-Nsemble
4680
1.069
0.3519
1.156
86.62
CDC
UGA_flucast-Scenariocast
4420
1.073
0.3424
1.012
78.33
CDC
*G-UGA_flucast-INFLAenza
4680
1.093
0.3493
1.282
93.91
Adapted
NEU_ISI-AdaptiveEnsemble
4294
1.104
0.3528
1.091
82.56
CDC
G-LANL_DBM v3
2808
1.118
0.3038
1.008
31.67
Adapted
*G-UMass_Flusion
3640
1.121
0.3143
0.978
58.98
Adapted
G-UGA_flucast_Copycat
4680
1.152
0.3795
1.161
89.57
Adapted
G-CMU_climate_baseline x UGA_flucast_Copycat
2912
1.168
0.3932
1.227
130.98 Hybrid
G-UMass_KCDE
3640
1.169
0.3291
1.327
78.55
Adapted
G-Cornell_JHU-hierarchSIR_2
3224
1.174
0.3306
1.465
70.13
Adapted
NAU-FourCAT
4472
1.175
0.3825
1.000
71.86
CDC
G-CU_SIRS x UVA_Gaussian_ processes
4680
1.175
0.3755
1.013
75.99
Hybrid
G-multi-layer-SE
4472
1.179
0.3905
1.339
103.09 Novel
CU-ARNB_Net
4056
1.184
0.3617
1.022
76.87
CDC
G-NU-PGF_FLUH v2
3016
1.190
0.3238
1.103
40.45
Adapted
CU-ensemble
4680
1.202
0.3876
1.163
88.81
CDC
PSI-PROF_MOA
4680
1.208
0.3969
1.384
106.05 CDC
G-CMU_climate_baseline x UGuelph_ CompositeCurve
4472
1.214
0.4154
1.142
94.25
Hybrid
G-CMU_timeseries
3640
1.224
0.3486
1.396
87.78
Adapted 28
Prospective multi-pathogen disease forecasting using autonomous LLM-guided tree search
Model
n tasks
Pairwise Pairwise Mean Rel. log Rel. log WIS WIS WIS
Mean WIS
Type
G-CMU_climate_baseline
3640
1.273
0.3576
1.198
70.19
Adapted
G-CFA_Pyrenew_ Pyrenew_H_Flu
4680
1.307
0.4533
1.396
110.81 Adapted
G-Cornell_JHU-hierarchSIR v3
2808
1.344
0.3638
1.156
36.18
Adapted
G-UMass-ar6_pooled
4472
1.349
0.4564
1.347
99.42
Adapted
G-Time_Series_to _Vision _Transfer _Learning
4472
1.358
0.4274
1.888
140.72 Novel
G-PSI-PROF_MOA v2
3432
1.368
0.3826
1.194
63.05
Adapted
G-Cornell_JHU-hierarchSIR v4
2808
1.409
0.3803
1.682
52.42
Adapted
G-PSI_PROF
4680
1.420
0.5072
1.313
104.16 Adapted
FluSight-baseline
4680
1.566
0.5299
1.488
114.54 CDC
G-UGuelph_ CompositeCurve
4680
1.593
0.5524
1.502
116.54 Adapted
G-UVA_Gaussian_ processes
3640
1.736
0.4864
1.819
109.63 Adapted
G-LANL_DBM
3224
1.787
0.4950
1.685
77.50
Adapted
G-NU-PGF_FLUH v4
2808
2.382
0.6329
2.004
62.02
Adapted
29
Prospective multi-pathogen disease forecasting using autonomous LLM-guided tree search
Figure 7 | National weekly observed influenza hospitalizations (black) and Google_SAI-FluEns forecast submissions by jurisdiction over the 2025–26 season with median (dark blue points) and corresponding 50% and 95% prediction intervals (blue shaded regions). Only every second set of submitted horizon forecasts is plotted for improved readability.
30
Prospective multi-pathogen disease forecasting using autonomous LLM-guided tree search
Figure 8 | Standardized log WIS rank distributions of top 15 ranked eligible models submitting scorable predictions for ≥ 80% of the task space to the CDC’s FluSight Forecast Hub for the 2025-26 season. A standardized rank of one indicates that the model had the best log WIS for that particular task (location, target, and horizon for that reference date), and a value of zero indicates it had the worst log WIS of submitting models. Density plots show interpolated distributions of standardized ranks achieved by each model for every forecast. Quartiles are colored from purple (bottom quarter, i.e., worst model ranks), through blue and light green to bright green (top quarter, i.e., best model ranks). Medians are represented by vertical white lines. Models are ordered by lowest quartile, with models that rarely had a low rank near the top. Google_SAI-FluEns is outlined in dark blue.
31
Prospective multi-pathogen disease forecasting using autonomous LLM-guided tree search
Figure 9 | FluSight Forecast Hub mean model performance (season average mean WIS; lower is better) by jurisdiction.
32
Prospective multi-pathogen disease forecasting using autonomous LLM-guided tree search
Table 5 | Mean log WIS by model and prediction horizon for models submitting to the FluSight Forecast Hub. Models are ranked by the sum of mean log WIS at horizons 0 and 1 (ascending). Within each horizon column, cells are shaded independently from green (lowest log WIS, best performance) to red (highest log WIS, worst performance); colors are not comparable across columns. The n tasks column reports total number of scored forecast tasks submitted across all horizons. Two CFA_Pyrenew models that do not submit forecasts for all four horizons are included; empty cells indicate horizons for which no forecast was submitted. Model
n tasks
0
1
2
3
CMU-TimeSeries
4,524
0.1950
0.2784
0.3558
0.4192
Google_SAI-FluEns
4,680
0.2057
0.2758
0.3380
0.3953
FluSight-trained_mean
4,680
0.1979
0.2853
0.3581
0.4248
UGA_flucast-INFLAenza
4,264
0.1981
0.2866
0.3686
0.4433
FluSight-HJudge_ensemble
4,680
0.2022
0.2895
0.3649
0.4254
OHT_JHU-nbxd
4,680
0.2097
0.2878
0.3444
0.3819
UMass-flusion
4,472
0.2129
0.2869
0.3590
0.4385
MIGHTE-Joint
3,848
0.2066
0.2966
0.3555
0.3976
FluSight-trained_med
4,680
0.2077
0.3020
0.3790
0.4455
UVAFluX-FS_OptimWISE
3,848
0.2218
0.2900
0.3384
0.3609
NAU-vulPES
4,680
0.2108
0.3036
0.3760
0.4345
FluSight-ensemble
4,680
0.2089
0.3057
0.3900
0.4585
FluSight-lop_norm
4,680
0.2166
0.3095
0.3919
0.4551
CFA_Pyrenew-Pyrenew_HE_Flu
2,083
0.2109
0.3167
UGA_flucast-Scenariocast
4,420
0.2215
0.3229
0.3941
0.4428
MIGHTE-Nsemble
4,680
0.2328
0.3244
0.3949
0.4731
UGA_CEID-auto_AVG_LB
4,472
0.2267
0.3502
0.4978
0.6391
NEU_ISI-AdaptiveEnsemble
4,294
0.2498
0.3337
0.4009
0.4417
PSI-PROF_MOA
4,680
0.2340
0.3604
0.4632
0.5534
CEPH-Rtrend_fluH
4,680
0.2497
0.3471
0.4610
0.5765
UGA_flucast-Copycat
4,680
0.2384
0.3627
0.4659
0.5521
CU-ARNB_Net
4,056
0.2495
0.3562
0.4212
0.4358
PSI-PROF
4,680
0.2446
0.3660
0.4649
0.5512
NAU-FourCAT
4,472
0.2546
0.3574
0.4427
0.4942
Gatech-ensemble_stat
4,680
0.2364
0.3764
0.4977
0.6161
Gatech-ensemble_prob
4,560
0.2411
0.3750
0.4900
0.5823
CU-ensemble
4,680
0.2636
0.3619
0.4389
0.5037
MOBS-EpyStrain_Flu
4,138
0.2495
0.3792
0.5153
0.6330
UVAFluX-Ensemble
4,628
0.2624
0.3757
0.4811
0.5450
NEU_ISI-FluBcast
4,106
0.2501
0.4020
0.5225
0.6188
MOBS-GLEAM_RL_FLUH
4,574
0.2634
0.3938
0.5093
0.6098
Cornell_JHU-hierarchSIR
4,680
0.2886
0.3732
0.4470
0.5099
33
Prospective multi-pathogen disease forecasting using autonomous LLM-guided tree search
Model
n tasks
0
1
2
3
CFA_Pyrenew-Pyrenew_H_Flu
2,355
0.2661
0.4221
NAU-epymorph
4,680
0.2927
0.3977
0.4812
0.5364
UMass-AR2
4,628
0.2925
0.4277
0.5445
0.6550
UNC_IDD-InfluPaint
4,590
0.3085
0.4303
0.5283
0.5966
UI_CompEpi-EpiGen
3,927
0.3027
0.4407
0.5638
0.6337
FluSight-baseline
4,680
0.3075
0.4715
0.6232
0.7500
VTSanghani-PRIME
4,116
0.3248
0.4731
0.5964
0.7137
NIH-Flu_ARIMA
3,958
0.3436
0.4599
0.5387
0.6099
UM-DeepOutbreak
4,420
0.3415
0.5493
0.6543
0.7069
LosAlamos-DoSiDo
4,680
0.3489
0.5466
0.7105
0.8541
LosAlamos-ThinMint
4,680
0.4078
0.5447
0.6661
0.7691
LosAlamos_NAU-CModel_Flu
4,676
0.4527
0.6275
0.7997
0.9486
UGA_CEID-Walk
4,420
0.4586
0.6293
0.7753
0.8964
A.2. COVID-19 This section includes analyses for COVID-19 forecast models, analyzed between 13 December 2025 and 2 May 2026 (21 reference dates, covering a total of 4,056 possible prediction tasks). Analyses are restricted to models submitting for ≥ 80% of the required task space, i.e., ≥ 3, 244 tasks (unless otherwise specified). Under this criterion, 12 hub-submitting models are eligible for analysis out of a total of 17 models which submitted at least one task prediction to COVIDHub this season. Table 6 | Summary performance of ERA-generated component models and eligible models submitting to COVIDHub. Models are ranked by pairwise relative log WIS (ascending) across submitted tasks, where a value of 1.0 indicates performance equal to the hub ensemble; values below 1.0 (green) indicate better performance and values above 1.0 (red) indicate worse performance than the ensemble. For each pair of models, score ratios are calculated over the intersection of tasks both models submitted, then geometrically averaged across all opponents, such that the metric is not affected by differences in the set of tasks each model chose to forecast 18 . The Mean log WIS and Mean WIS columns are shaded from green (lowest, best) to red (highest, worst) across models.
Model
n tasks
Pairwise Pairwise Mean Rel. log Rel. log WIS WIS WIS
Mean WIS
Type
Google_SAI-Ensemble
4056
0.978
0.2225
0.876
14.78
CDC
CovidHub-ensemble
4056
1.000
0.2296
1.000
17.04
CDC ensemble
*G-CMU_TimeSeriesUMass_gbqr
4056
1.005
0.2289
0.946
15.92
Hybrid
*G-CADPH-CovidCAT_Ensemble
2392
1.014
0.2554
0.906
12.67
Adapted
UGA_flucast-INFLAenza
4056
1.016
0.2344
0.966
16.76
CDC
*G-DeepResearch _Counterfactual Simulation
4056
1.020
0.2313
0.892
14.80
Novel
34
Prospective multi-pathogen disease forecasting using autonomous LLM-guided tree search
Model
n tasks
Pairwise Pairwise Mean Rel. log Rel. log WIS WIS WIS
Mean WIS
Type
*G-UM-DeepOutbreak
3224
1.049
0.2421
0.898
12.61
Adapted
*G-UMass-gbqr
3224
1.050
0.2426
0.956
13.37
Adapted
CMU-TimeSeries
4056
1.062
0.2421
1.101
18.31
CDC
*G-DeepResearch _Regime Switching Detection
4056
1.090
0.2465
0.964
16.20
Novel
G-JHU_CSSE-CSSE_Ensemble
3016
1.107
0.2601
1.071
14.18
Adapted
*G-CMU_climate_baselineUMass_ar6_pooled
4056
1.125
0.2568
1.093
18.52
Hybrid
*G-MOBS-GLEAM_COVID
3224
1.134
0.2611
0.990
13.82
Adapted
NEU_ISI-AdaptiveEnsemble
3515
1.151
0.2485
1.182
20.37
CDC
CEPH-Rtrend_covid
4056
1.160
0.2671
1.159
20.04
CDC
G-CMU-TimeSeries
3016
1.180
0.2766
1.153
15.22
Adapted
G-UGA_flucast-INFLAenza
3016
1.180
0.2766
1.131
14.93
Adapted
UMass-gbqr
4004
1.187
0.2738
1.203
20.42
CDC
G-UMass-ar6_pooled
3016
1.198
0.2820
1.083
14.44
Adapted
G-Metaculus-cp
3016
1.207
0.2829
1.165
15.39
Adapted
*G-CEPH_Rtrend_covidCMU_climate_baseline
4056
1.210
0.2732
1.127
18.93
Hybrid
CFA-EpiAutoGP
3549
1.241
0.2845
1.330
21.53
CDC
G-CFA-EpiAutoGP
3016
1.254
0.2940
1.176
15.54
Adapted
UMass-ar6_pooled
4004
1.255
0.2868
1.332
21.75
CDC
*G-NEU_ISI-AdaptiveEnsemble
3224
1.294
0.2975
1.042
14.52
Adapted
OHT_JHU-nbxd
4056
1.379
0.3162
1.481
25.57
CDC
CovidHub-baseline
4056
1.390
0.3231
1.240
21.37
CDC
G-CEPH-Rtrend_covid
3016
1.397
0.3262
1.323
17.42
Adapted
UM-DeepOutbreak
3796
1.621
0.3695
1.511
25.00
CDC
35
Prospective multi-pathogen disease forecasting using autonomous LLM-guided tree search
Figure 10 | National weekly observed COVID-19 hospitalizations (black) and Google_SAI-Ensemble forecast submissions by jurisdiction over the 2025–26 season with median (dark blue points) and corresponding 50% and 95% prediction intervals (blue shaded regions). Only every second set of submitted horizon forecasts is plotted for improved readability.
36
Prospective multi-pathogen disease forecasting using autonomous LLM-guided tree search
Figure 11 | Standardized log WIS rank distributions of top 12 ranked eligible models submitting scorable predictions for ≥ 80% of the task space to the CDC’s COVIDHub for the 2025-26 season. A standardized rank of one indicates that the model had the best log WIS for that particular task (location, target, and horizon for that reference date), and a value of zero indicates it had the worst log WIS of submitting models. Density plots show interpolated distributions of standardized ranks achieved by each model for every forecast. Quartiles are colored from purple (bottom quarter, i.e., worst model ranks), through blue and light green to bright green (top quarter, i.e., best model ranks). Medians are represented by vertical white lines. Models are ordered by lowest quartile, with models that rarely had a low rank near the top. Google_SAI-Ensemble is outlined in dark blue.
37
Prospective multi-pathogen disease forecasting using autonomous LLM-guided tree search
Figure 12 | COVIDHub mean model performance (season average mean WIS; lower is better) by jurisdiction.
38
Prospective multi-pathogen disease forecasting using autonomous LLM-guided tree search
Figure 13 | Prospective performance of ERA-generated adaptations for COVID-19 forecasting. Paired comparison of prospective Mean Weighted Interval Score (WIS; lower is better) on a logarithmic scale between the target method from the CDC COVIDHub (deep charcoal circles) and its corresponding ERA-generated adaptation (old gold diamonds). Models are ranked sequentially from lowest to highest target mean WIS. Horizontal ash-grey lines span the performance delta between each paired target and adaptation, while alternating light-grey background bands anchor the rows. Sample sizes (𝑛) are integrated directly into the model labels on the y-axis and denote the intersection of forecasting tasks mutually completed by both the target method and the ERA adaptation, serving as the shared baseline over which the pairwise mean WIS was calculated. A.3. RSV Here we provide analyses for RSV forecast models, analyzed between 3 January and 2 May 2026 (18 reference dates, covering a total of 3, 432 maximum forecast tasks per submitting model). Analyses are restricted to models submitting for ≥ 80% of the required task space, i.e., ≥ 2, 745 tasks (unless otherwise specified). Under this criterion, four hub-submitting models are eligible for analysis out of a total of seven models which submitted at least one task prediction to the RSVHub this season. Table 8 | Summary performance of ERA-generated component models and eligible models submitting to the RSVHub. Models are ranked by pairwise relative log WIS (ascending) across submitted tasks, where a value of 1.0 indicates performance equal to the hub ensemble; values below 1.0 (green) indicate better performance and values above 1.0 (red) indicate worse performance than the ensemble. For each pair of models, score ratios are calculated over the intersection of tasks both models submitted, then geometrically averaged across all opponents, such that the metric is not affected by differences in the set of tasks each model chose to forecast 18 . The Mean log WIS and Mean WIS columns are shaded from green (lowest, best) to red (highest, worst) across models.
Model
n tasks
Pairwise Rel. log WIS
Mean log WIS
Pairwise Rel. WIS
Mean WIS
Type
Google_SAI-RSVEns
3432
0.887
0.2354
0.875
17.95
CDC
*G-general_instruction_3
3432
0.928
0.2443
0.880
18.03
Novel
*G-general_instruction_1
3432
0.961
0.2515
0.993
19.95
Novel
RSVHub-ensemble
3432
1.000
0.2600
1.000
20.06
CDC ensemble 39
Prospective multi-pathogen disease forecasting using autonomous LLM-guided tree search
Model
n tasks
Pairwise Rel. log WIS
Mean log WIS
Pairwise Rel. WIS
Mean WIS
Type
RSVHub-baseline
3432
1.209
0.3045
1.087
21.50
CDC
CEPH-Rtrend_rsv
3432
1.210
0.3048
1.242
24.04
CDC
40
Prospective multi-pathogen disease forecasting using autonomous LLM-guided tree search
Figure 14 | National weekly observed RSV hospitalizations (black) and Google_SAI-RSVEns forecast submissions by jurisdiction over the 2025–26 season with median (dark blue points) and corresponding 50% and 95% prediction intervals (blue shaded regions). Only every second set of submitted horizon forecasts is plotted for improved readability.
41
Prospective multi-pathogen disease forecasting using autonomous LLM-guided tree search
Table 7 | Mean log WIS by model and prediction horizon for models submitting to the CDC’s COVIDHub. Models are ranked by the sum of mean log WIS at horizons 0 and 1 (ascending). Within each horizon column, cells are shaded independently from green (lowest log WIS, best performance) to red (highest log WIS, worst performance); colors are not comparable across columns. The n tasks column reports total number of scored forecast tasks submitted across all horizons. Two CFA_Pyrenew models that do not submit forecasts for all four horizons are included; empty cells indicate horizons for which no forecast was submitted. Model
n tasks
0
1
CFA_Pyrenew-Pyrenew_HE_COVID Google_SAI-Ensemble CovidHub-ensemble UGA_flucast-INFLAenza CMU-TimeSeries NEU_ISI-AdaptiveEnsemble CEPH-Rtrend_covid UMass-gbqr CFA-EpiAutoGP UMass-ar6_pooled CFA_Pyrenew-Pyrenew_H_COVID CovidHub-baseline OHT_JHU-nbxd UM-DeepOutbreak
1,974 4,056 4,056 4,056 4,056 3,515 4,056 4,004 3,549 4,004 2,061 4,056 4,056 3,796
0.1687 0.1825 0.1778 0.1874 0.1919 0.2002 0.2069 0.2008 0.2063 0.2145 0.2247 0.2405 0.2604 0.2598
0.2121 0.2048 0.2113 0.2221 0.2280 0.2317 0.2446 0.2507 0.2588 0.2610 0.3055 0.2997 0.3025 0.3436
2
3
0.2390 0.2542 0.2562 0.2605 0.2659 0.2968 0.3089 0.3175 0.3186
0.2715 0.2845 0.2798 0.2971 0.3056 0.3308 0.3521 0.3703 0.3710
0.3607 0.3450 0.4139
0.4057 0.3662 0.4739
B. Prospective Model Inventory The following tables list the models generated or tracked during the prospective phase of the study, extracted from the operational records.
42
Prospective multi-pathogen disease forecasting using autonomous LLM-guided tree search
Figure 15 | Standardized log WIS rank distributions of ranked eligible models submitting scorable predictions for ≥ 80% of the task space to the CDC’s RSVHub for the 2025-26 season. A standardized rank of one indicates that the model had the best log WIS for that particular task (location, target, and horizon for that reference date), and a value of zero indicates it had the worst log WIS of submitting models. Density plots show interpolated distributions of standardized ranks achieved by each model for every forecast. Quartiles are colored from purple (bottom quarter, i.e., worst model ranks), through blue and light green to bright green (top quarter, i.e., best model ranks). Medians are represented by vertical white lines. Models are ordered by lowest quartile, with models that rarely had a low rank near the top. Google_SAI-RSVEns outlined in dark blue. Note that, with only four models considered, an invisible Q4 (Google_SAI-RSVEns) implies that in more than 25% of cases, the model ranked first, while an invisible Q1 (RSVHub-Baseline) implies that the model had zero rank in more than 25% of considered cases.
Table 9 | Mean log WIS by model and prediction horizon for models submitting to the CDC’s RSVHub. Models are ranked by the sum of mean log WIS at horizons 0 and 1 (ascending). Within each horizon column, cells are shaded independently from green (lowest log WIS, best performance) to red (highest log WIS, worst performance); colors are not comparable across columns. The n tasks column reports total number of scored forecast tasks submitted across all horizons. Two CFA_Pyrenew models that do not submit forecasts for all four horizons are included; empty cells indicate horizons for which no forecast was submitted. Model
n tasks
0
1
2
3
CFA_Pyrenew-Pyrenew_HE_RSV Google_SAI-RSVEns RSVHub-ensemble CEPH-Rtrend_rsv RSVHub-baseline CFA_Pyrenew-Pyrenew_H_RSV
1,627 3,432 3,432 3,432 3,432 1,775
0.1756 0.1839 0.1829 0.2078 0.2244 0.2098
0.2280 0.2238 0.2339 0.2759 0.2750 0.2921
0.2552 0.2976 0.3454 0.3345
0.2893 0.3421 0.4106 0.4020
43
Prospective multi-pathogen disease forecasting using autonomous LLM-guided tree search
Figure 16 | RSVHub mean model performance (season average mean WIS; lower is better) by jurisdiction.
Figure 17 | Flowchart showing the model selection process for the prospective study, from candidate prompt construction for ERA through internal hub model selection and final ensemble inclusion and submission for influenza, COVID and RSV.
44
Prospective multi-pathogen disease forecasting using autonomous LLM-guided tree search
Table 10 | Comparison of different LLMs in following instructions and performance WIS. Method
Model
Nodes
Runtime (hr)
WIS
Judgment
gemini-2.5-flash
927 341 1017 1267 747 1440 344 394 446
2488.21 2495.12 2496.71 2498.22 2494.13 2499.55 2497.31 2499.68 2490.9
176.14 191.81 169.33 148.98 176.27 150.54 138.28 153.54 130.03
Partial Match Partial Match Partial Match Partial Match Partial Match Match Match Match Match
2500 2500 2500 2500 2500 2500 2500 2500 2500
219.94 187.12 254.04 106.35 125.49 120.35 280.73 1930.41 1478.84
332.31 195.42 180.44 150.59 159.98 151.35 126.65 141.59 152.04
Partial Match Partial Match Partial Match Partial Match Partial Match Partial Match Partial Match Partial Match Partial Match
2498 2500 2500 2458 2370 2386 2500 2500 956
120.64 105.51 79.68 89.8 73.82 114.83 156.67 115.5 2493.22
207.47 427.79 313.21 181.22 172.17 176.26 145 151.67 116.67
Partial Match Partial Match Match Partial Match Partial Match Match Partial Match Partial Match Partial Match
912 1121 2500 851 997 790 358 442 353
2492.3 2499.94 62.93 2498.52 2497.47 2498.76 2496.17 2497.78 2495.51
183.02 159.44 531.67 159.9 176.96 164 149.91 160.8 154.53
Match Match Match Match Match Partial Match Partial Match Match Match
gemini-3-flash-preview Cornell_JHU-hierarchSIR gemini-3-pro-preview
gemini-2.5-flash
gemini-3-flash-preview NU-PGF_FLUH gemini-3-pro-preview
gemini-2.5-flash
gemini-3-flash-preview UGA_flucast-INFLAenza gemini-3-pro-preview
gemini-2.5-flash
gemini-3-flash-preview UMass-gbqr gemini-3-pro-preview
45
Prospective multi-pathogen disease forecasting using autonomous LLM-guided tree search
Table 11 | Comparison of Gemini 2.5 Flash performance with and without judge in the loop. Method
No Judge Nodes Runtime Score Judgment
927 Cornell_JHU-hierarchSIR 341 1017
Judge in Loop Nodes Runtime Score Judgment
2488.21 176.14 Partial Match 2031 2495.12 191.81 Partial Match 357 2496.71 169.33 Partial Match 701
2511.67 172.54 No Match 2489.63 531.69 Partial Match 2499.10 531.69 Match
NU-PGF_FLUH
2500 2500 2500
219.94 187.12 254.04
332.31 Partial Match 1742 195.42 Partial Match 872 180.44 Partial Match 276
2499.13 166.94 Partial Match 2497.58 148.72 Partial Match 2494.35 299.68 Partial Match
UGA_flucast-INFLAenza
2498 2500 2500
120.64 105.51 79.68
207.47 Partial Match 2500 427.79 Partial Match 2500 313.21 Match 2500
1020.40 185.45 Match 1172.44 238.76 Match 989.20 165.32 Match
UMass-gbqr
912 1121 2500
2492.3 183.02 Match 2499.94 159.44 Match 62.93 531.67 Match
427 467 632
2491.09 147.17 Match 2493.92 177.14 Match 2494.53 178.55 Match
Table 12 | Inventory of Influenza Models in the Prospective Study. Cells in Val and Test columns are colored on a gradient from light green (best) to light red (poor), with scores capped at 80 for color calculation. Model Name
Category
LANL_DBM_and_UMass_Flusion Double Adapt. LANL_DBM_and_LANL_Inferno Double Adapt. Cornell_JHU-hierarchSIR Single Adapt.: FluSight UGA_flucast-INFLAenza Single Adapt.: FluSight UMass_Flusion Single Adapt.: FluSight CU_SIRS_and_CMU_climate_baseline Double Adapt. NU-PGF_FLUH (node 2027) Single Adapt.: FluSight CMU_climate_baseline_and_UGA_flucast Double Adapt. PSI-PROF_MOA Single Adapt.: FluSight NU-PGF_FLUH (node 838) Single Adapt.: FluSight CU_SIRS_and_UVA_Gaussian_processes Double Adapt. multi-layer-SE Novel LANL_DBM (node 717) Single Adapt.: Paper UMass-gbqr Single Adapt.: COVIDHub CU_SIRS Single Adapt.: Paper Cornell_JHU-hierarchSIR (node 2211) Single Adapt.: FluSight UGA_flucast_Copycat Single Adapt.: FluSight NU-PGF_FLUH Single Adapt.: FluSight CMU_timeseries Single Adapt.: FluSight Time_Series_to_Vision Novel NU-PGF_FLUH (node 321) Single Adapt.: FluSight UMass-ar6_pooled Single Adapt.: COVIDHub CMU_climate_baseline_and_UGuelph Double Adapt. CFA_Pyrenew_Pyrenew_H_Flu Single Adapt.: FluSight Cornell_JHU-hierarchSIR (Pro) Single Adapt.: FluSight LANL_DBM (Pro) Single Adapt.: Paper PSI_PROF Single Adapt.: FluSight UGuelph_CompositeCurve Single Adapt.: FluSight UVA_Gaussian_processes Single Adapt.: Paper
Hub? Ens? Y Y Y Y Y Y Y Y Y Y Y Y Y Y Y Y Y Y Y Y Y Y Y Y Y Y Y Y Y
Y Y Y Y Y Y Y N N N N N N N N N N N N N N N N N N N N N N
Val
Test
59.3 63.1 72.5 51.4 60.3 63.9 61.1 62.0 69.9 64.8 59.9 60.7 73.4 67.8 56.6 73.2 69.1 78.4 65.9 63.5 76.8 70.6 73.8 73.9 79.0 65.3 84.8 81.3 80.3
42.1 42.2 41.0 51.8 52.6 50.9 45.8 50.9 47.2 44.7 45.4 52.6 48.5 51.5 58.6 50.9 53.9 49.5 55.8 59.3 53.5 57.3 55.8 56.0 53.4 62.2 56.3 59.7 69.6
46
Prospective multi-pathogen disease forecasting using autonomous LLM-guided tree search
Table 12 continued from previous page Model Name
Category
UMass_KCDE CMU_climate_baseline PSI-PROF_MOA (v2) Cornell_JHU-hierarchSIR (gemini3) LANL_DBM (gemini 3) NU-PGF_FLUH (with packages) Cornell_JHU-hierarchSIR (pkg) LANL_DBM (with packages) NU-PGF_FLUH (gemini 3) LANL_DBM_and_UGuelph CFA_Pyrenew_and_UMass_Flusion CU_SIRS_and_UGuelph CMU_climate_and_LANL_DBM CU_SIRS_and_UMass-gbqr CFA_Pyrenew_and_CMU_climate Two_Stage_Classifier CMU_climate_and_CU_SIRS CU_SIRS_and_UMass_KCDE CU_SIRS_and_UGA_flucast PSI_PROF_and_UVA_Gaussian CU_SIRS_and_UGA_Copycat CU_SIRS_and_UMass_Flusion CFA_Pyrenew_and_UMass_KCDE PSI_PROF_and_CMU_climate CMU_climate_and_CFA_Pyrenew CMU_climate_and_UMass-ar6 LANL_DBM_and_CMU_timeseries LANL_DBM_and_UGA_Copycat CMU_climate_and_UVA_Gaussian CFA_Pyrenew_and_UGuelph PSI_PROF_and_UMass-gbqr CFA_Pyrenew_and_UGA_Copycat DGP-SAD CFA_Pyrenew_and_LANL_Inferno PSI_PROF_and_UGA_Copycat PSI_PROF_and_UGuelph CU_SIRS_and_UMass-ar6 CMU_climate_and_UMass_KCDE LANL_DBM_and_UMass_KCDE PSI_PROF_and_UMass_KCDE AETF LANL_DBM_and_UVA_Gaussian CMU_climate_and_PSI_PROF LANL_DBM PSI_PROF_and_LANL_Inferno LANL_DBM_and_UMass-gbqr epi-bridge CMU_climate_and_UMass-gbqr CMU_climate_and_LANL_Inferno CU_SIRS_and_CMU_timeseries LANL_DBM_and_CMU_climate neuro-mech-QSE NU-PGF_FLUH
Single Adapt.: Paper Single Adapt.: COVIDHub Single Adapt.: FluSight Single Adapt.: FluSight Single Adapt.: Paper Single Adapt.: FluSight Single Adapt.: FluSight Single Adapt.: Paper Single Adapt.: FluSight Double Adapt. Double Adapt. Double Adapt. Double Adapt. Double Adapt. Double Adapt. Novel Double Adapt. Double Adapt. Double Adapt. Double Adapt. Double Adapt. Double Adapt. Double Adapt. Double Adapt. Double Adapt. Double Adapt. Double Adapt. Double Adapt. Double Adapt. Double Adapt. Double Adapt. Double Adapt. Novel Double Adapt. Double Adapt. Double Adapt. Double Adapt. Double Adapt. Double Adapt. Double Adapt. Novel Double Adapt. Double Adapt. Single Adapt.: Paper Double Adapt. Double Adapt. Novel Double Adapt. Double Adapt. Double Adapt. Double Adapt. Novel Single Adapt.: FluSight
Hub? Ens? Y Y Y Y Y N N N N N N N N N N N N N N N N N N N N N N N N N N N N N N N N N N N N N N N N N N N N N N N N
N N N N N N N N N N N N N N N N N N N N N N N N N N N N N N N N N N N N N N N N N N N N N N N N N N N N N
Val
Test
80.5 89.4 129.9 72.7 120.0 74.7 86.3 224.3 53.4 60.2 64.5 58.3 61.6 65.7 55.9 69.1 56.6 60.6 54.1 68.4 54.4 63.5 58.3 44.7 63.5 56.8 60.3 65.3 63.4 64.3 67.9 62.4 76.2 77.0 63.5 75.0 53.7 58.7 49.7 58.4 67.7 60.6 61.6 63.2 74.4 62.7 81.2 68.3 69.7 67.1 61.8 73.4 68.5
72.0 79.1 148.6 57.9 75.6 48.0 56.3 128.3 68.9 46.0 46.0 49.8 48.9 48.3 53.9 47.9 54.6 53.4 57.5 50.4 58.3 53.8 56.7 63.7 54.4 58.6 57.0 55.2 56.3 57.3 55.9 58.7 52.0 51.6 58.9 53.4 64.4 61.9 67.9 63.7 60.3 63.9 63.4 63.0 58.1 64.0 54.9 61.4 61.1 62.7 66.4 61.7 65.5
47
Prospective multi-pathogen disease forecasting using autonomous LLM-guided tree search
Table 12 continued from previous page Model Name
Category
PSI_PROF_and_UMass-ar6 CFA_Pyrenew_and_UMass-ar6 PI-GTN CFA_Pyrenew_and_CMU_timeseries EpiSERA LANL_Inferno DERAS LANL_DBM_and_UMass-ar6 CFA_Pyrenew_and_UGA_flucast PSI_PROF_and_CMU_timeseries ASST-Flu PSI_PROF_and_UGA_flucast CFA_Pyrenew_and_UVA_Gaussian GEPLABS CMU_climate_and_UGA_flucast CU_SIRS_and_LANL_Inferno
Double Adapt. Double Adapt. Novel Double Adapt. Novel Single Adapt.: Paper Novel Double Adapt. Double Adapt. Double Adapt. Novel Double Adapt. Double Adapt. Novel Double Adapt. Double Adapt.
Hub? Ens? N N N N N N N N N N N N N N N N
N N N N N N N N N N N N N N N N
Val
Test
55.5 72.7 47.7 76.8 86.9 64.1 86.4 64.6 81.4 68.0 64.5 78.2 78.2 72.3 64.0 83.3 47.3 102.4 102.5 86.0 176.1 91.4 55.4 161.0 61.2 232.4 252.6 137.2 53.4 264.4 55.2 15128.0
48
Prospective multi-pathogen disease forecasting using autonomous LLM-guided tree search
Table 13 | Inventory of COVID-19 Models in the Prospective Study. Cells in Val and Test columns are colored on a gradient from light green (best) to light red (poor), with scores capped at 45 for color calculation. Model Name
Category
UM-DeepOutbreak MOBS-GLEAM_COVID DeepResearch_CounterfactualSimulation NEU_ISI-AdaptiveEnsemble UMass-gbqr OHT_JHU-nbxd CMU_TimeSeries-UMass_gbqr CADPH-CovidCAT_Ensemble DeepResearch_RegimeSwitchingDetection CEPH_Rtrend_covid-CMU_climate_baseline CMU_climate_baseline-UMass_ar6_pooled UMass-ar6_pooled CFA-EpiAutoGP Metaculus-cp UGA_flucast-INFLAenza CMU-TimeSeries JHU_CSSE-CSSE_Ensemble CMU-climate_baseline CEPH-Rtrend_covid
Single Adapt.: COVIDHub Single Adapt.: COVIDHub Last Season Search Single Adapt.: COVIDHub Single Adapt.: COVIDHub Single Adapt.: COVIDHub Last Season Search Single Adapt.: COVIDHub Last Season Search Last Season Search Last Season Search Single Adapt.: COVIDHub Single Adapt.: COVIDHub Single Adapt.: COVIDHub Single Adapt.: COVIDHub Single Adapt.: COVIDHub Single Adapt.: COVIDHub Single Adapt.: COVIDHub Single Adapt.: COVIDHub
Hub?
Ens?
Val
Test
Y Y Y Y Y N Y Y Y Y Y Y Y Y Y Y Y Y Y
Y Y Y Y Y N Y Y Y Y Y N N N N N N N N
34.5 35.1 35.9 39.1 36.3 35.8 37.9 41.3 38.4 37.4 41.1 41.1 43.4 42.0 55.1 55.7 52.6 61.1 50.2
26.2 27.8 27.9 26.7 28.3 28.8 28.0 27.0 29.3 30.2 29.4 30.5 30.3 31.1 30.9 30.8 37.6 44.3 116.7
49
Prospective multi-pathogen disease forecasting using autonomous LLM-guided tree search
Table 14 | Inventory of RSV Models in the Prospective Study. Cells in Val and Test columns are colored on a gradient from light green (best) to light red (poor), with Val capped at 30 and Test capped at 8 for color calculation. Model ID/Name
Category
general_instruction_1 general_instruction_2 general_instruction_3 bsts_nowcasting CFA_Pyrenew-Pyrenew_E_RSV UGA_flucast-INFLAenza CFA_Pyrenew-Pyrenew_HE_RSV adversarial_domain_adapt multitask_shared_encoder bayesian_transfer_priors UM-DeepOutbreak CFA_Pyrenew-Pyrenew_H_RSV hybrid_mechanistic_ude deepar_cold_start CEPH-Rtrend_rsv nhits_hierarchical_interp
No specified method No specified method No specified method Novel Single Adapt.: RSVHub Single Adapt.: RSVHub Single Adapt.: RSVHub Novel Novel Novel Single Adapt.: RSVHub Single Adapt.: RSVHub Novel Novel Single Adapt.: RSVHub Novel
Hub?
Ens?
Val
Test
Y N Y N N N N N N N N N N N N N
Y N Y N N N N N N N N N N N N N
23.2 24.9 21.7 28.5 27.5 27.4 27.6 23.0 25.5 57.6 36.1 94.1 94.1 94.1 48.5 65.9
2.5 3.1 3.2 3.2 3.3 3.5 3.6 4.6 8.0 8.9 9.3 10.3 10.4 10.4 15.4 20.6
Table 15 | EPIFORGE 2020 Reporting Checklist 30 . #
Section
Checklist item
Reported in
1
Title/Abstract
Title; Abstract
2
Introduction
3 4
Methods Methods
5 6
Methods Methods
7
Methods
Describe as forecast or prediction research Define purpose and forecasting targets Fully document the methods Prospective, real-time, or retrospective Origin of input source data Provide source data or document why not Data processing procedures
8 9 10 11 12 13 14 15 16 17 18 19
Methods Methods Methods Methods Methods Methods Results Results Results Discussion Discussion Discussion
Model type and assumptions Model code availability Model validation Accuracy evaluation method Benchmark comparator Forecast horizon Uncertainty of results Nontechnical summary Time-stamped data object Weaknesses Public health implications Generalizability
Introduction, ¶3–5 Methods, Sections 1–5 Methods §1 ¶1; Results Part I & II headers Methods §1 ¶4 (NHSN HRD citation) NHSN publicly available; ; code on GitHub Methods §1 ¶4 (truncation, interpolation) Methods §3.1; SI Tables 5–7 Methods §2; code on GitHub Methods §3.3; Table 1; Figure 1d Methods §1 ¶5 (WIS); Part II §2 Part I §2 (CDC ensembles); Figs. 3–4 Methods §1 ¶2–3; Figure 1c Figs. 3a–c; SI forecast figures Abstract GitHub time-stamped submissions Discussion ¶4 Discussion ¶2–3 Discussion ¶5
50
Prospective multi-pathogen disease forecasting using autonomous LLM-guided tree search
Figure 18 | Fidelity versus predictive accuracy across model tiers and experimental conditions. Scatter plots illustrate the relationship between structural scientific fidelity and forecasting performance (WIS) for four epidemiological methods (rows) across different Large Language Model (LLM) search agents and judging configurations (columns). Each point represents an independent search trajectory replica. The horizontal axis indicates the Fidelity Score (0.0 to 1.0) as determined by the expert-calibrated LLM-Judge, representing how faithfully the generated code followed complex methodological instructions (e.g., mechanistic SIR structures). The vertical axis represents the Weighted Interval Score (WIS), where lower values indicate superior predictive accuracy. Configurations include the baseline Gemini 2.5 Flash without judging, Gemini 2.5 Flash with the automated judge-in-the-loop, and the frontier-class Gemini 3 Flash and Gemini 3 Pro models.
51
Prospective multi-pathogen disease forecasting using autonomous LLM-guided tree search
C. Problem Statement Prompt Each ERA experiment begins with a shared problem statement that defines the forecasting task, input/output specification, evaluation metric, and available datasets. Because the three pathogens share the CDC Forecast Hub infrastructure, the problem statement follows a common template; only the target variable, auxiliary data sources, and data-augmentation guidance differ across infections. Below we reproduce the canonical prompt (shown for influenza) in full, followed by a summary of the infection-specific variations for COVID-19 and RSV. C.1. Canonical Problem Statement (Flu) Problem Statement — FLU (FluSight)
Overview. Modelers at the CDC’s FluSight Forecast Hub are tasked with producing unconditional probabilistic forecasts for weekly influenza hospitalizations across the United States. This is a critical public health task, as these forecasts help inform resource allocation and policy decisions during flu season. The challenge is to create forecasts that characterize uncertainty across all reasonable future scenarios, not just a limited set of conditions. This task aims to develop a superhuman forecasting model for this problem. The input data is a time series of historical flu hospitalizations and related signals. The output must be a set of quantile forecasts for all 50 states, Washington DC, and Puerto Rico, across multiple time horizons. The goal is to develop a model that is more accurate and better calibrated than those produced by leading human experts. Problem Statement & Deliverable. Your primary task is to create a forecasting model that predicts probabilistic forecasts of Total Influenza Hospital Admissions for every US state and jurisdiction. The goal is to create a model that achieves the lowest possible Weighted Interval Score (WIS) over a rolling-window evaluation. Your deliverable is a single Python function, fit_and_predict_fn, that takes in training and test data and returns a pandas DataFrame containing the required quantile predictions. Function Signature & Output Requirements. The forecasting model must be encapsulated within a function named fit_and_predict_fn with the following signature: def fit_and_predict_fn( train_x: pd.DataFrame, train_y: pd.Series, test_x: pd.DataFrame, ) -> pd.DataFrame: return test_y_hat_quantiles The returned DataFrame must have its index match the input test_x, columns named by quantile (e.g. quantile_0.01, quantile_0.5, quantile_0.975), and the predicted quantiles for any given row must be monotonically increasing.
Dataset Description. The following data objects are available: • Primary Training Data: train_x (historical features), train_y (historical target: Total Influenza Admissions). • Historical Augmentation Data: ilinet_hhs, ilinet, ilinet_state — DataFrames containing ∼20 years of historical Influenza-Like Illness (ILI) data, available only for dates before 2022-10-15. • Reference & Example Data: locations (geographic/population data), sample_submission_df (correct output format), example_train_x, example_train_y, example_test_x (small example DataFrames).
Feature Definitions. target_end_date: Saturday of the epiweek; location_name: full state name; location: FIPS code; population: total population; Total Influenza Admissions: the target variable (available from late 2020 onward).
52
Prospective multi-pathogen disease forecasting using autonomous LLM-guided tree search
Augmenting Training Data with Historical ILINet Data. The core challenge is the limited history of the target variable. To overcome this, ∼20 seasons of historical ILINet data are provided. While not the same target, it is highly correlated and captures the essential seasonal dynamics of influenza. Two strategies are suggested: Strategy 1 — Standardize and Combine: (1) Apply a standardization method to both datasets independently to make the “shape” of the seasons comparable. (2) Treat the standardized historical ILINet data as additional, independent flu seasons and append them to the training data. (3) Train the model on this combined “library” of seasons. Strategy 2 — Learn a Transformation: (1) Identify the date range where both the target and the historical ILINet data overlap. (2) Use this overlapping period to learn a statistical transformation (e.g. linear regression, quantile mapping) that maps the ILINet data onto the same scale as Total Influenza Admissions. (3) Apply this transformation to the entire 20-year history to create a “synthetic” history for the target variable. (4) Train the final model on this augmented training set.
Key Considerations. • Time Series Awareness: Handle seasonality, trends, and lags appropriately. • Calibration: Ensure predicted quantiles are well-calibrated. • Logging: Configure models to be quiet during training (verbose=0). Do not suppress critical warnings or tracebacks.
Detailed Instructions. An expert has instructed the model to implement a specific method for this forecasting task. Minor improvements to the method are permitted, but the original core principles must be maintained. Before writing code, the model must explicitly list the 3–4 core principles of the prescribed method in a comment block, then implement accordingly.
C.2. Infection-Specific Variations The COVID-19 and RSV problem statements share the same function signature, output format, evaluation metric (WIS), and key considerations as the influenza prompt above. The differences are summarized below. COVID-19 Variations • Target variable: Total COVID-19 Admissions (replaces Total Influenza Admissions). • Forecast Hub: CDC CovidHub Forecast Hub (replaces FluSight). • No auxiliary historical data: The COVID-19 prompt does not include ILINet augmentation data or data-augmentation strategies. The model relies solely on the primary training data (train_x, train_y) and the reference/example data. • Simpler dataset: No ilinet_* DataFrames are provided. RSV Variations • Target variable: Total RSV Admissions (replaces Total Influenza Admissions). • Core challenge — Extreme Data Sparsity: For most states, the RSV target is only available from late 2024 onward, making this the most data-constrained of the three tasks. The prompt explicitly frames data sparsity as the central challenge and emphasizes “transfer of knowledge” from auxiliary datasets. • Richer auxiliary data: In addition to the ∼20-year ILINet history, the RSV prompt provides: – Total COVID-19 Admissions and Total Influenza Admissions columns in train_x (longer history than the RSV target). – Weekly percentage of ED visits due to RSV from the NSSP dataset (available from 2022 onward, but not for all locations, e.g. Puerto Rico and Missouri are missing).
53
Prospective multi-pathogen disease forecasting using autonomous LLM-guided tree search
• Data-augmentation strategies (replacing the ILI-focused strategies in the flu prompt): – Strategy 1 — Transfer Learning from NSSP: Learn a mapping between the NSSP signal and the RSV target in their overlap period, then apply it to create a synthetic history. – Strategy 2 — Multi-Target Modeling: Use the COVID-19, flu, and NSSP columns as features and learn cross-correlations to inform RSV forecasts.
D. Method Prompts The following tables list the natural-language method descriptions provided to ERA as search instructions for each model submitted to the internal Google Research hub during the prospective season. Models sharing the same base prompt (e.g., different search replicas or LLM variants of the same method) are listed once. Prompts are organized by pathogen and by the source of the methodological description. D.1. Method Descriptions for Influenza (FluSight) Source: Published Literature LANL_DBM “A hierarchical Bayesian model for forecasting weekly influenza-like illness proportions. It combines a mechanistic Susceptible-Infectious-Recovered (SIR) model with two discrepancy terms to improve predictions, operating on the logit scale of the true unobserved proportion of illness. The model assumes the observed proportion follows a Beta distribution. The logit of the true proportion is the sum of the logit of the infectious proportion from the SIR model, a common discrepancy term across seasons, and a season-specific discrepancy term. The common discrepancy is modeled as a reverse random walk. The season-specific discrepancy is modeled as an autoregressive reverse random walk, constrained at the final time point, with parameters shared hierarchically across seasons. Priors for the SIR model parameters are also informed by previous seasons. 1. Model the observed weekly weighted influenza-like illness proportion for season *j* and week *t*, *yj,t*, using a Beta distribution: *yj,t* ~ Beta( 𝜆 *𝜋j,t*, 𝜆 (1 - *𝜋j,t*)), where *𝜋j,t* is the true unobservable proportion and 𝜆 is a pre-set concentration parameter. 2. Model the logit of the true proportion as the sum of three components: logit(*𝜋j,t*) = logit(*Ij,t*) + *𝜇 t* + *𝛿j,t*. 3. *Ij,t* is the infectious proportion from a standard SIR model, solved numerically using the fourth-order Runge-Kutta method. The initial susceptible proportion *Sj,0* is fixed at 0.9. Empirical Bayes priors are used for the SIR parameters *Ij,0*, * 𝛽 j*, and * 𝜌j* = 𝛾 j/* 𝛽 j*, based on fits from other seasons, with appropriate truncations. 4. The common discrepancy term *𝜇 t* is modeled as a reverse random walk: * *𝜇 T* ~ N(0, 𝜎2𝜇 T) * *𝜇 t* | *𝜇 t+1* ~ N(*𝜇 t+1*, 𝜎2𝜇 ) for *t* = T-1, ..., 1. * Gamma priors are placed on the precisions 𝜎-2 𝜇 T and 𝜎-2𝜇 . 5. The season-specific discrepancy term *𝛿j,t* is constrained at *𝛿j,T* = -logit(*Ij,T*) and follows an autoregressive reverse random walk: * *𝛿j,t* | *𝛿j,t+1* ~ N(𝛼j *𝛿j,t+1*, 𝜎2𝛿,j) for *t* = T-1, ..., 1. 6. The autoregressive parameters 𝛼j are modeled hierarchically: logit(𝛼j) ~ TN(logit(0.02), logit(0.98))(logit(0.9), 𝜎2𝛼), with a Gamma prior on 𝜎-2𝛼. 7. The season-specific discrepancy precisions 𝜎-2𝛿,j are modeled hierarchically: 𝜎-2𝛿,j ~ Gamma(a𝛿, b𝛿), with Gamma hyperpriors on a𝛿 and b𝛿. 8. Use Markov Chain Monte Carlo (MCMC) to sample from the posterior distributions of the parameters and generate forecasts. The execution environment enables probabilistic programming frameworks including PyMC and PyStan for hierarchical Bayesian modeling and MCMC sampling of the posterior distributions.”
54
Prospective multi-pathogen disease forecasting using autonomous LLM-guided tree search
UMass_Flusion “An ensemble forecasting approach for weekly influenza hospital admissions, integrating three distinct data sources: the target hospital admission counts, historical hospitalization rates from a sentinel network, and influenza-like illness data combined with test positivity. The method combines predictions from three component models: two tree-based quantile regression models trained jointly on all data sources and locations, differing in their feature sets, and an autoregressive time series model with a holiday effect covariate trained only on the target hospitalization data across all locations. Data are preprocessed through rate conversion, a power transformation, and standardization. The final ensemble forecast is produced by averaging the quantile predictions from the component models. 1. **Data Acquisition and Preparation:** * Obtain weekly state and national level data for: * NHSN influenza hospital admissions (the primary target). * FluSurv-NET hospitalization rates. * ILI+ (ILI data combined with influenza test positivity rates). * Adjust FluSurv-NET and ILI+ historical data to improve consistency over time, if applicable (though the study found these adjustments counterproductive). 2. **Data Standardization:** For each data source, location, and time point: * Convert NHSN counts to rates per 100,000 population. * Apply a fourth-root transformation to stabilize variance. * Center and scale the data by dividing by the 95th percentile and subtracting the mean for each location and data source. 3. **Component Model 1: Tree-Based Quantile Regression (Full Features):** * Train separate tree-based models to predict each required quantile of the *change* in the standardized signal from the last observed value. * Use data from all three sources (NHSN, FluSurv-NET, ILI+) and all locations (state, regional, national). * Include features such as: data source, location, week of season, time until Christmas, forecast horizon, and various measures of the recent local level, trend, and curvature of the standardized signal. * Employ bagging by averaging predictions from models trained on random subsets of seasons. 4. **Component Model 2: Tree-Based Quantile Regression (No Level Features):** * Similar to Component Model 1, but exclude features that directly measure the local level of the signal (e.g., most recent value, rolling means, intercepts of polynomial fits). 5. **Component Model 3: Autoregressive with Covariates:** * Fit an autoregressive model of order 8 to the standardized NHSN data only, jointly across all locations. * Share autoregressive coefficients across locations. * Include a covariate representing proximity to Christmas week. * Estimate separate innovation variance parameters for each location. * Generate multi-step ahead forecasts by iterating one-step-ahead predictions. 6. **Ensemble Formation:** * For each forecast target (location, horizon) and each required quantile level, average the predicted quantile values obtained from the three component models. * Invert the standardization steps (mean/scaling, fourth root, population adjustment) to obtain forecasts on the original scale of hospital admission counts. 7. **Prediction Generation:** * Generate predictions for the current week and the next three weeks.” LANL_Inferno “A Bayesian model for forecasting influenza-like illness (ILI) percentages, treated independently for each geographical unit. The model uses historical data to heuristically estimate parameters defining a typical seasonal ILI trajectory on the logit scale, season-specific deviations, and observation noise. Parameters for the mean vector and covariance matrix of the season-specific deviations (modeled as a multivariate normal distribution) are estimated from past seasons. Forecasts for the ongoing season are generated by conditioning on observed data and sampling future ILI values using Markov Chain Monte Carlo (MCMC). 1. **Estimate True ILI (𝜃s,t):** For each historical season *s* and week *t*, smooth the observed ILI/100 (ys,t) using a 3-week moving average. Calculate and add a holiday effect
55
Prospective multi-pathogen disease forecasting using autonomous LLM-guided tree search
adjustment, estimated as the average difference between ys,t and the moving average across seasons. This yields 𝜃ˆs,t. 2. **Estimate Observation Noise (𝛼):** Estimate the Beta distribution precision parameter 𝛼 by maximizing the likelihood of ys,t given 𝜃ˆs,t across all historical data. 3. **Estimate Typical Logit ILI (𝛾𝑡 ):** Compute 𝛾ˆ𝑡 as the average of logit(𝜃ˆ𝑠,𝑡 ) across 2 ):** Calculate all historical seasons for each week 𝑡 . 4. **Estimate Season-level Variance (𝜎𝑚 ˆ ˆ weekly deviations 𝛿𝑠,𝑡 = logit(𝜃𝑠,𝑡 ) - 𝛾ˆ𝑡 . Compute the mean seasonal deviation 𝜇ˆ 𝑠 for each 2 as the sample variance of these 𝜇 season. Estimate 𝜎ˆ𝑚 ˆ 𝑠 values. 5. **Estimate Covariance 2 Parameters (𝜎𝑆 , 𝜆 , 𝜙):** Estimate the total variance 𝜎ˆ𝑆2 of the centered deviations (𝛿ˆ𝑠,𝑡 - 𝜇ˆ 𝑠 ). Estimate the squared exponential covariance function parameters 𝜆 and 𝜙 by maximizing the multivariate normal likelihood of the observed 𝛿ˆ𝑠 vectors, using the estimated 𝜇ˆ 𝑠 and 𝜎ˆ𝑆2 . 6. **Generate Forecasts:** For the current season, given observations up to the current week, 2, sample future ILI values using MCMC. The model uses the estimated parameters {ˆ 𝛼, 𝛾 ˆ, 𝜎ˆ𝑚 2 ˆ ˆ 𝜎 ˆ𝑆 , 𝜆 , 𝜙} to define the priors and likelihood for the unobserved future weeks.” CU_SIRS “A data assimilation system is used to forecast seasonal influenza outbreaks by integrating realtime influenza-like illness (ILI) estimates with an ensemble of simulations from a dynamic transmission model. The model is a compartmental model where individuals transition between Susceptible, Infectious, and Recovered states, with transmissibility influenced by absolute humidity. An ensemble of model simulations, each with slightly different initial conditions and parameters, is run forward in time. Weekly observational data on ILI prevalence are used to adjust the state variables (e.g., number of susceptible and infected individuals) and key epidemiological parameters of each ensemble member. This adjustment process nudges the ensemble mean towards the observations and reduces the ensemble spread, effectively refining the model’s representation of the current state of the epidemic. Forecasts are generated by running the updated ensemble forward in time, and the spread of these forecast trajectories provides an estimate of prediction uncertainty. 1. Initialize an ensemble of 200 simulations of a susceptible-infectious-recovered-susceptible (SIRS) model. Each ensemble member starts with unique, randomly sampled state variables (number of susceptible and infected individuals) and parameters (average duration of immunity, mean infectious period, maximum and minimum basic reproductive numbers) drawn from prior distributions in September. 2. Force each SIRS model simulation with observed daily absolute humidity (AH) data for New York City. 3. On a weekly basis, assimilate weekly Google Flu Trends (GFT) estimates of ILI using an ensemble adjustment Kalman filter. This updates the ensemble’s state variables (S, I) and parameters to be more consistent with the observations. A small multiplicative inflation (e.g., 1.02) is applied during assimilation. 4. After each weekly assimilation, generate a forecast by integrating each member of the updated (posterior) ensemble forward for 300 days. These forecast integrations use the first 5 days of true future AH conditions followed by daily climatological AH conditions for the remaining 295 days. 5. Analyze the ensemble forecasts to predict metrics such as the timing of the peak week of the influenza outbreak. 6. Assess forecast skill by comparing predictions to actual observed peak timing in retrospective analyses. 7. Use the ensemble spread (e.g., variance in predicted peak timing) as a measure of forecast confidence, with lower spread indicating higher confidence.” UVA_Gaussian_processes “This method forecasts weekly disease incidence using only past weekly incidence data. It transforms the incidence data (e.g., modified square root) and models the transformed series. The model structure captures non-linear relationships between weeks within and across seasons based on four derived features: the week of the season, a periodic component, the incidence level at the start of the season, and a categorical indicator of the season’s overall severity (e.g., mild, moderate, severe). This severity indicator is treated as a latent variable
56
Prospective multi-pathogen disease forecasting using autonomous LLM-guided tree search
for the season being forecast and is updated as weekly data becomes available. Crucially, the variance of the errors is allowed to differ depending on the estimated season severity category, making the model adaptable to different noise levels in low vs. high incidence seasons. Predictions for future weeks and seasonal summaries (like peak timing and magnitude) are generated by sampling from the model’s predictive distribution, considering different hypotheses for the latent season severity. 1. Transform the weekly incidence counts using the function ‘f(x) = sqrt(x + 1) - 1‘. 2. For each week in the historical data, create four predictor variables: a. ‘x1‘: Week number within the season (1 to 52). b. ‘x2‘: A sinusoidal function of ‘x1‘ to model cyclic behavior, for example, ‘sin(2*pi*x1/52)‘. c. ‘x3‘: The transformed incidence value of the last week of the *previous* season. This value is constant for all 52 weeks within the current season. For the first season, use the transformed value of its first week. d. ‘x4‘: A categorical variable representing season severity: -1 for mild, 0 for moderate, and 1 for severe, based on whether the maximum weekly incidence in that season crosses predefined thresholds (e.g., <25 for mild, >100 for severe in San Juan). 3. Fit a non-linear regression model to the transformed incidence data using ‘x1, x2, x3, x4‘ as inputs. The model assumes a multivariate normal distribution for the responses, with a covariance structure dependent on the distances between input vectors. The covariance matrix is ‘𝜏^2 * (C + Λ)‘, where ‘C‘ is a correlation matrix derived from a kernel function (e.g., squared exponential) applied to the inputs, and ‘Λ‘ is a diagonal matrix. 4. Implement heteroskedasticity: The diagonal elements of ‘Λ‘ (nugget terms) depend on the severity category ‘x4‘ of the observation’s season. Three separate nugget parameters (‘𝜂_{-1}, 𝜂_0, 𝜂_{+1}‘) are estimated, one for each severity level. 5. Forecast a new season: a. At the start of the season (week 0), consider three hypotheses for the new season’s severity, corresponding to ‘x4 = -1, 0, 1‘. b. For each hypothesis, initialize a continuous latent variable for severity (e.g., to -1, 0, or 1). c. As each new week of data arrives: i. Update the weights of the three severity hypotheses based on the predictive log-likelihood of the observed data under each nugget regime (‘𝜂_{-1}, 𝜂_0, 𝜂_{+1}‘). ii. For each hypothesis, optimize the continuous latent severity variable to maximize the predictive log-likelihood of the data observed so far in the current season. d. Generate Monte Carlo samples of future weekly incidence trajectories from the model’s predictive distribution. This involves drawing from multivariate normal distributions conditioned on the observed data, under each of the three nugget/severity hypotheses, and combining them based on their current weights. e. Derive distributions for peak week, peak incidence, and total season incidence from the Monte Carlo samples. 6. Inverse transform the predictions and summaries back to the original count scale using ‘f^{-1}(y) = (y + 1)^2 - 1‘ for ‘y >= 0‘, and 0 otherwise.” UMass_KCDE “This semi-parametric method predicts infectious disease incidence by first estimating separate probability distributions for incidence in individual future weeks. These estimations are based on recent past incidence values and the time of year, using a flexible, data-driven weighting of historical observations. These individual weekly predictive distributions are then linked together using a statistical model to capture the dependence between incidence levels across different weeks. This combined model allows for the creation of a joint probability distribution over the entire forecast period, enabling predictions not only for specific weeks but also for characteristics of the entire season, such as the week of peak incidence and the intensity of that peak. 1. **Estimate Weekly Conditional Densities:** For each desired future week (prediction horizon *h*), estimate the conditional probability density of disease incidence. This density is conditioned on observed incidence from a few recent past time points and the current time within the season. This estimation uses a kernel-based approach. * The kernel function is a product of two components: * A periodic kernel to account for seasonality, based on the time of year. * A multivariate kernel (e.g., log-normal for continuous data, or a discretized version for count data) for the incidence values, allowing for a fully parameterized bandwidth
57
Prospective multi-pathogen disease forecasting using autonomous LLM-guided tree search
matrix to capture relationships between predictor lags and the prediction target. * Select bandwidth parameters for the kernels by maximizing a cross-validated log-likelihood score on the training data. 2. **Model Joint Dependence with Copulas:** Combine the marginal conditional densities for each week (from Step 1) into a joint predictive distribution for the entire sequence of future weeks in the season. This is done using an isotropic normal copula function, whose parameters are estimated to capture the temporal dependence between weekly incidence levels, typically after the kernel parameters are fixed. 3. **Derive Target Predictions:** * Incidence at horizon h: Use the marginal density from Step 1 directly. * Timing of Peak Week & Peak Incidence: Sample trajectories from the joint distribution constructed in Step 2. Use Monte Carlo integration on these samples to estimate the probability distribution for the week number with the highest incidence and the distribution of the incidence value in that peak week. Binned incidence is used for peak intensity predictions.”
Source: FluSight Hub Cornell_JHU-hierarchSIR “An SIR model with unknown case ascertainment, basic reproduction number, population immunity and a splined effective reproduction number is used to model seasonal influenza dynamics in a given season. Across-season trends (’hyperparameters’) in the SIR model’s parameters are derived by wrapping it in an across-season Bayesian hierarchical model. Hyperparameters are used as priors when forecasting the current season. Disease model integrated in C++ and bound to Python with pybind11, Bayesian hierarchical posterior probability coded in raw Python and sampled using the ensemble sampler of Goodman and Weare available in ‘emcee‘ (motivation: computationally inefficient but amazingly robust). Workflow automatically pulls NHSN HRD data through a timed GH actions and deploys it on a local runner (Dell Optiplex 3050 Micro running Ubuntu Server). Average time from data pull to forecast ready: 3 hours. The environment enables ‘scipy.interpolate‘ for spline modeling and ‘emcee‘ for Bayesian posterior sampling with an ensemble MCMC sampler.” UGA_flucast-INFLAenza “A spatial time-series model for estimating influenza activity. The model is structured to jointly fit several components: seasonality terms to capture regular patterns in flu incidence, statespecific connectivity to account for spatial dependencies and interactions between different locations, and a random walk component to model temporal trends and random fluctuations. Estimation of the model parameters is performed using a Bayesian statistical approach. 1. Employ a spatial time-series framework to analyze and forecast influenza trends across different states. 2. Incorporate components into the model to represent the characteristic seasonal variations inherent in influenza outbreaks. 3. Include parameters that reflect the degree of connection and potential for transmission between different states. 4. Utilize a random walk process to account for underlying temporal evolution and stochastic elements in the data. 5. Jointly estimate all model parameters, including those for seasonality, inter-state connectivity, and the random walk, using a unified Bayesian inference method designed for hierarchical models. 6. Generate probabilistic forecasts of future influenza activity based on the fitted model.” NU-PGF_FLUH “A branching process model analytically solved with probability generating functions (PGFs), enabling Bayesian inference of the time-varying reproduction number and mechanistic modeling of immunity. The branching process approximates the full stochastic evolution of a discrete-time, single-population, SLIR compartmental model, extended to evaluate weekly
58
Prospective multi-pathogen disease forecasting using autonomous LLM-guided tree search
hospitalizations resulting from infections. Enabled packages that may be useful for this approach include SymPy for symbolic manipulation of probability generating functions, and PyMC for Bayesian inference of the time-varying reproduction number and latent epidemic states.” PSI-PROF_MOA “The PROF routines perform a deterministic fit of our compartmental SIR[H]2 (including vaccinated compartments for Susceptibles and Infectious) model to weekly hospitalization incidence profiles. The model includes a hospitalization compartment which is split into two sub-compartments. This split ensures that the model preserves the correct generation time (Tg) and that the ratio between cumulative recovered and hospitalized individuals is determined by the infection-hospitalization-ratio. The transmission-rate coefficient (Beta) is a time-dependent function of 2 or more arc-tangents. The actual data is augmented with ’future’ data from a Method Of Analogs (MOA). The MOA data receives a reduced fitting-weight compared to actual data. The model fit is inferred by an MCMC procedure. It is followed by stochastic simulations through the forecast time-window using the inferred parameter distributions and a daily cadence which is then aggregated to weekly incidence to produce the forecasts. Where there is little-to-no epi signal in the data, a baseline statistical model is substituted for the mechanistic model.” UGA_flucast_Copycat “Matches seasonal growth rate trends against historic growth rate curves to identify the closest matches. Makes forecasts based on nearest neighbor trajectories. 1. Calculate current seasonal growth rate trends from available data. 2. Compare these current growth rate trends against a database of historic growth rate curves. 3. Identify the historic growth rate curves that are the "closest matches" to the current trends. 4. Generate forecasts by projecting forward based on the trajectories of these identified nearest neighbor historic curves.” CMU_timeseries “A basic quantile auto-regression fit using lagged values of influenza-related hospitalization counts (normalized by population). The data are whitened so that the historical flusurv and ILI data can be used to augment the number of training examples. The model is fit jointly across all 50 US states, the District of Columbia, Puerto Rico, and the Virgin Islands, using a window of data 7 weeks centered on the forecast date, including all past years. Each of the 23 quantiles is learned using a separate quantile regression with nonnegativity and quantile sorting constraints applied post hoc. All data signals are available as indicators through the Delphi Epidata API (https://cmu-delphi.github.io/delphi-epidata). This model is averaged with both a climatological model and a simple linear model, with the weight on the windowed seasonal model above being ~6 times that of either the linear or climatological model. 1. **Data Acquisition:** Obtain influenza-related hospitalization counts (e.g., FluSurv-NET from NHSN) and ILI-like data (e.g., percentage of ED visits from NSSP) for all US states, DC, Puerto Rico, and the Virgin Islands from the Delphi Epidata API. 2. **Data Preprocessing:** * Normalize the influenza-related hospitalization counts by population for each location. * Apply a whitening process to the historical time series data to enable the use of both hospitalization and ILI data sources to increase training data size. 3. **Windowed Quantile Auto-Regression:** * For each location and forecast date, select training data from all available past years, restricted to a 7-week window centered on the current week of the year. * Fit quantile auto-regression models using lagged values of the preprocessed data as predictors. * Train separate models for each of the 23 required quantiles (0.01, 0.025, 0.05, ..., 0.95, 0.975, 0.99). The models are fit jointly across all locations. * Apply post hoc adjustments
59
Prospective multi-pathogen disease forecasting using autonomous LLM-guided tree search
to ensure forecast quantiles are non-negative and monotonically increasing. 4. **Baseline Model Generation:** * Generate quantile forecasts using a climatological model based on historical data patterns. * Generate quantile forecasts using a simple linear regression model. 5. **Ensemble Averaging:** * Combine the quantile forecasts from the three models: the Windowed Quantile Auto-Regression, the Climatological model, and the Linear model. * The final ensemble forecast is a weighted average, where the Windowed Quantile Auto-Regression model output receives approximately six times the weight of the Climatological model’s output, and also six times the weight of the Linear model’s output.” CFA_Pyrenew_Pyrenew_H_Flu “Simple renewal model fit only to the target signal. Used to benchmark related multi-signal models the same framework with additional observables. 1. **Model Latent Infections:** Estimates daily latent infections using a renewal equation. The time-varying reproduction number (Rt) is modeled as a latent random walk on the log scale, evolving weekly. A fixed generation interval distribution is used. 2. **Estimate Latent Hospitalizations:** Calculates daily latent hospitalizations by convolving the daily latent infections with a fixed hospitalization interval distribution (delay from infection to hospital admission) and multiplying by a latent infection-to-hospitalization rate (IHR). 3. **Observation Model:** Links the modeled daily latent hospitalizations to the observed weekly hospital admissions data. Daily latent hospitalizations are aggregated to the weekly level (MMWR epi-weeks). A Negative Binomial distribution is used as the observation model to account for overdispersion in the hospitalization counts. 4. **Parameter Inference:** Fits the model parameters, including those for the Rt random walk, initial infection levels, the IHR, and the Negative Binomial dispersion, to the observed hospitalization data using a Bayesian statistical framework.” PSI_PROF “The PROF routines perform a deterministic fit of our compartmental SIR[H]2 model to daily hospitalization incidence profiles. The model includes a hospitalization compartment which is split into two sub-compartments. This split ensures that the model preserves the correct generation time (Tg) and that the ratio between cumulative recovered and hospitalized individuals is determined by the infection-hospitalization-ratio. The transmission-rate coefficient (beta) is a time-dependent function of 2 or more arc-tangents. The model fit is inferred by an MCMC procedure. It is followed by stochastic simulations through the forecast time-window using the inferred parameter distributions and a daily cadence which is then aggregated to weekly incidence to produce the forecasts. Where there is little-to-no epi signal in the data, a baseline statistical model is substituted for the mechanistic model. This model is similar to the previous model ’PSI-DICE’, but in an updated computational framework (https://github.com/predsci/PROF). 1. Ingest publicly-available daily confirmed hospital admission incidence data for the target location. 2. Fit a compartmental SIR[H]2 model to the hospitalization data. The hospitalization compartment is split into two sub-compartments to maintain correct generation times. 3. Define the transmission rate (beta) within the model as a flexible, time-dependent function composed of two or more arc-tangents. 4. Use a Markov Chain Monte Carlo (MCMC) procedure to infer the joint posterior distribution of the model parameters, including those defining the time-varying transmission rate. 5. If there is a sufficient epidemiological signal in the data, generate ensemble forecast trajectories by running stochastic simulations of the compartmental model forward over the forecast time window, sampling parameters from the inferred posterior distributions. Daily simulations are aggregated to weekly incidence. 6. If there is little-to-no epidemiological signal in the recent data, substitute a baseline statistical model for the mechanistic model to generate forecasts.”
60
Prospective multi-pathogen disease forecasting using autonomous LLM-guided tree search
UGuelph_CompositeCurve “We produce a peak-aligned composite curve epidemic curves across past seasons and regions, normalizing the season for percent of cases each week, and then taking the median percentage across the resulting curves for each number of weeks from the peak, following the approach of Schanzer et al., Influenza and Other Respiratory Viruses 2010. Then, for each state, we compute the normalized case count for the preceding four weeks, and compare against the normalized four week windows of the composite curve, following the approach of Morel et al. PLOS Computational Biology, 2023. sample_method: “Each window along the curve is weighted based on its distance from the most recent window and these weights are used to select 100 windows (allowing repetition) along the curve. These windows are used to produce the 100 sample trajectories by applying the subsequent scaling factors to the observed data.” 1. **Construct Composite Curve:** a. Collect historical weekly disease incidence data for multiple past seasons and regions. b. For each historical season/region, identify the week of peak incidence. c. Normalize each season’s curve by calculating the percentage of the total seasonal incidence that occurred in each week. d. Align all the normalized seasonal curves based on their week of peak incidence (week 0 at the peak). e. For each week relative to the peak (-n to +m weeks), calculate the median percentage across all aligned historical curves. This sequence of median percentages forms the "peak-aligned composite curve," representing a typical season’s shape. 2. **Prepare Current State Data:** a. Obtain the most recent four weeks of observed disease incidence data for the target state. b. This four-week window of data is used for comparison. 3. **Compare and Weight Potential Futures:** a. Create sliding four-week windows along the entire length of the composite curve. b. Compare the pattern (shape and relative changes) of the state’s most recent four-week data window to each four-week window from the composite curve. The comparison method assesses the similarity, potentially accounting for differences in scale and minor temporal shifts, inspired by time-warping approaches. c. Calculate a weight for each window on the composite curve based on its similarity to the state’s current four-week window. Higher similarity results in a higher weight. The weight is inversely related to the "distance" between the windows. 4. **Sample Future Trajectories:** a. Based on the weights calculated in step 3c, select 100 four-week windows from the composite curve. This selection is done probabilistically, where windows with higher weights are more likely to be chosen (sampling with replacement). b. Each selected window represents a potential "current stage" of the epidemic within the typical season defined by the composite curve. 5. **Generate Forecasts:** a. For each of the 100 sampled windows: i. Determine a scaling factor. This factor adjusts the magnitude of the composite curve’s projections to match the level of the state’s currently observed incidence data. It’s derived from the ratio between the state’s recent data and the corresponding values in the selected composite curve window. ii. Project future weekly incidence by taking the sequence of median percentages from the composite curve, starting from the week immediately following the selected four-week window, and multiplying them by the scaling factor. b. This process results in 100 plausible sample trajectories for future disease incidence in the state. c. These 100 trajectories are used to generate the probabilistic forecasts for various targets (e.g., weekly incidence, peak timing, peak intensity).”
61
Prospective multi-pathogen disease forecasting using autonomous LLM-guided tree search
Source: COVIDHub CMU_climate_baseline “Using data from 2022 onwards, the climatological model uses samples from the 7 weeks centered around the target week and reference week to form the quantiles for the target week, as one might use climate information to form a meteorological forecast. To get more variation at some potential issue of generalization, one can form quantiles after aggregating across geographic values as well as years (after converting to a rate based case count). This model uses a simple average of the geo-specific quantiles and the geo-aggregated quantiles. 1. **Data Selection:** Use weekly flu hospitalization data, including only data from the 2022-2023 flu season onwards. 2. **Define Sample Window:** For a given target week, establish a 7-week window centered on that week number (target week number ± 3 weeks). 3. **Calculate Geo-Specific Quantiles:** a. For each geographic location, collect the historical hospitalization counts that fall within the 7-week window identified in Step 2, across all included seasons. b. From this location-specific set of historical values, calculate the required forecast quantiles. 4. **Calculate Geo-Aggregated Quantiles:** a. Convert the historical hospitalization counts to rates (e.g., cases per 100,000 population) for all locations and included seasons. b. Pool all historical rates from all locations that fall within the 7-week window identified in Step 2. c. Calculate the required forecast quantiles from this aggregated set of rates. d. Convert these rate-based quantiles back to counts for each specific geographic location using its population. 5. **Combine Forecasts:** For each location and each quantile level, the final forecast value is the simple average of the value from the Geo-Specific calculation (Step 3b) and the value from the Geo-Aggregated calculation (Step 4d).” UMass-gbqr “The model uses gradient boosting to produce probabilistic forecasts of influenza hospital admissions. It is trained on multiple influenza surveillance data streams, including the target signal of hospital admissions reported to the National Healthcare Safety Network (NHSN), as well as influenza-like illness data combined with virologic testing (ILI+) and data on laboratory-confirmed influenza hospitalizations from a sentinel network (FluSurv-NET). The model is trained jointly on data from multiple locations (state, regional, and national levels). Features provided to the model include those summarizing recent signal activity (such as levels, trends, and curvature), properties of the location (e.g., population, spatial scale), information about the timing of forecast creation (e.g., week of the season, proximity to holidays), and the forecast horizon (0 to 3 weeks ahead). To handle the different scales and variances of the input data, transformations are applied, including a power transform and standardization. The model predicts the change in the transformed signal, which is then converted back to the original scale. Separate gradient boosting models are fit to predict each of the 23 required quantiles for the predictive distribution, and bagging is used to improve stability. 1. **Data Collection and Preparation:** Gather weekly data from three sources: NHSN, FluSurv-NET, and ILI+. Adjustments are made to FluSurv-NET and ILI+ to enhance consistency over time. 2. **Data Transformation:** Standardize the data across signals and locations. This involves: * Converting NHSN counts to rates per 100,000 people. * Applying a fourth-root transformation to stabilize variance. * Centering and scaling the transformed data based on historical statistics for each signal and location. 3. **Feature Engineering:** Create features for the model, including: * Indicators for the data source, location, and spatial scale. * Location population. * Week of the season and difference from Christmas week. * Forecast horizon (0, 1, 2, or 3 weeks). * Measures of the most recent signal value, trend, and curvature derived from the past few weeks of data. * Lagged versions of these signal activity measures. 4. **Model Training:** * Train separate models for each of the 23 required predictive quantiles. * Use gradient boosting, optimizing the quantile loss function. * Train
62
Prospective multi-pathogen disease forecasting using autonomous LLM-guided tree search
jointly on historical data from all three surveillance signals and all available locations. * The prediction target is the difference between the transformed signal value at the forecast horizon and the most recently observed value. * Employ bagging: train multiple models on different random subsets of the available seasons and average their quantile predictions. * Exclude data from pandemic influenza seasons and summer off-seasons during training. 5. **Prediction Generation:** * Generate predictions for the change in the transformed signal. * Add the most recent observed transformed value to the predicted change. * Invert the data transformations (scaling, centering, and power transform) to obtain forecasts on the original scale of hospital admissions.” UMass-ar6_pooled “AR(6) model with exogenous regressors for holiday effects, applied after a fourth root data transform. AR and exogenous regressor coefficients are shared across all locations. A separate variance parameter is estimated for each location. 1. Load weekly hospitalization data for each location. 2. Apply a fourth root transformation to the hospitalization data. This involves custom centering and scaling, then taking the fourth root. 3. Create an additional covariate to account for increased hospitalizations around Christmas. 4. Fit a single time series model to the transformed data across all locations. The model includes: * Autoregressive terms up to lag 6 (AR(6)). * The Christmas holiday covariate. * The coefficients for the autoregressive terms and the holiday covariate are shared (pooled) across all locations. * A separate parameter for the error variance is estimated for each location. 5. Generate multi-step ahead forecast distributions from the fitted model on the transformed scale. 6. Calculate the required prediction quantiles from the forecast distributions. 7. Invert the data transformations to return the quantile forecasts to the original scale of hospitalization counts. This involves raising to the power of 4 and reversing the centering and scaling. 8. Format the quantile forecasts into the standard CDC FluSight Hub submission file.”
Source: Deep Research multi-layer-SE “The "Superhuman" strategy is a Multi-Layered Stacking Ensemble (Level-1) built from a portfolio of diverse, SOTA hybrid models (Level-0) designed to solve the Proxy-Target Data Asymmetry and Spatio-Temporal dynamics challenges. The ensemble uses a trained meta-forecaster to learn the optimal, time-varying fusion weights for superior probabilistic forecasting. 1. **Level-0 Component Design: Proxy-Target Transfer (MT-PatchTST):** * **Architecture:** Implement the Multi-Task PatchTST (MT-PatchTST), a time-series Transformer with **Patching** (segmenting the time series into subseries-level patches, e.g., 7-day patches) and **Channel-Independence** (a single shared Transformer encoder applied independently to each variable). * **Implementation:** * **Pre-training (Self-Supervised):** Train the PatchTST backbone on the full ~20-season ILINet dataset using a **masked-patch "infilling" task** to learn fundamental seasonal dynamics. * **Fine-Tuning (Multi-Task):** Fine-tune the pre-trained backbone on the <5 years of paired ILINet and NHSN data, using two output heads to simultaneously predict ILINet (auxiliary task) and the 23 NHSN quantiles (target task). 2. **Level-0 Component Design: Hybrid State-Space (DSSM-Integrator):** * **Architecture:** Implement the Deep State Space Model (DSSM) Integrator, which models ILINet and NHSN as two different observations ($NN_{g1}$ and $NN_{g2}$) of the same underlying, unobserved latent "flu state" ($l_t$). The latent state’s dynamics ($l_t = f(l_{t1})$) are governed by an RNN ($RNN_f$). * **Implementation (Two-Period Training):** * **Period 1 (2000-2020):** Train the model on ILINet data only, learning the dynamics
63
Prospective multi-pathogen disease forecasting using autonomous LLM-guided tree search
($RNN_f$) and the ILI emission function ($NN_{g1}$). * **Period 2 (2021-Present):** Continue training on both data streams, learning the NHSN emission function ($NN_{g2}$) from the paired data while refining the dynamics ($RNN_f$). 3. **Level-0 Component Design: Spatio-Temporal (HeatGNN-Flu):** * **Architecture:** Implement the HeatGNN-Flu model, a hybrid GNN framework that uses a dynamic graph. * **Implementation:** The GNN is trained to learn embeddings that predict the *parameters* of a local SIR model (e.g., the local transmission rate $\beta_i(t)$) for each jurisdiction. The message-passing function uses these **epidemiology-informed embeddings** to model heterogeneous spread. 4. **Generative Probabilistic Decoders:** * Couple each Level-0 model’s feature vector ($h_t$) with a SOTA generative decoder (e.g., **Conditional Diffusion Model/TimeGrad** or **Mixture Density Network/MDN**). * **TimeGrad Sampling:** For Diffusion Models, generate the final quantiles by running the reverse diffusion process **10,000 times** to create an empirical distribution, and then sort and select the 23 required quantiles from this distribution. * **MDN Parameter Output:** For MDNs, the network’s final layer must output the parameters (mean $\mu$, standard deviation $\sigma$, and weight $w$) for a **Gaussian Mixture Model (GMM)**, from which the 23 quantiles can be analytically derived. 5. **Level-1 Stacking Ensemble Training:** * Generate a Level-1 training dataset by using a time-series cross-validation scheme to produce **out-of-sample forecasts** (the 23 quantiles) from all Level-0 component models. * Train the Level-1 meta-forecaster on this dataset, using the Level-0 forecasts as input features and the true NHSN target values as labels. * The objective is to train a model to *learn* the optimal weighted combination to combine the diverse component forecasts, making the L1 output mathematically superior to a simple average.”
64
Prospective multi-pathogen disease forecasting using autonomous LLM-guided tree search
D.2. Method Descriptions for COVID-19 (COVIDHub) Source: COVIDHub UM-DeepOutbreak “A sequence-to-sequence deep neural network with self-attention mechanisms, trained in a multi-task setting and calibrated using adaptive conformal inference. 1. **Tensor Construction:** Prepare 3D tensors of input sequences (hospitalizations, search trends, covariates). 2. **Seq2Seq Encoding:** Pass sequences through Recurrent Neural Network (RNN/GRU) layers and a Self-Attention module to encode historical context into a latent vector. 3. **Multi-Task Training:** Train the network by minimizing a global loss function that aggregates error across all geographic regions simultaneously, allowing weight sharing. 4. **Decoding:** Use the decoder to autoregressively generate point forecast sequences. 5. **Conformal Calibration:** Apply Adaptive Conformal Inference (ACI) to the output. Calculate residuals on a hold-out set and dynamically adjust the width of the prediction intervals based on recent coverage error to guarantee valid uncertainty.” MOBS-GLEAM_COVID “A stochastic, age-structured metapopulation model that simulates disease transmission within geographic patches and spatial spread via real-world air travel and commuting networks. 1. **Network Construction:** Tessellate the target area into census-based subpopulations connected by a mobility graph weighted by daily commuter and airline flows. 2. **Compartmental Definition:** Implement an age-structured SLIR (Susceptible-Latent-Infectious-Removed) model within each subpopulation. 3. **Stochastic Simulation:** Initialize agents and simulate time steps using a Monte Carlo approach: update agent locations based on mobility fluxes, then update infection status based on local prevalence and contact matrices. 4. **Ensemble Generation:** Execute independent simulation runs to generate a distribution of possible outcomes. 5. **Aggregation:** Sum agent counts to the state level and compute quantiles across the simulations.” NEU_ISI-AdaptiveEnsemble “A meta-ensemble that filters and aggregates long-term scenario projections by rejecting trajectories that diverge significantly from recent ground-truth data. 1. **Ingestion:** Load the set of long-term projection trajectories from the Scenario Modeling Hub. 2. **Validation:** Calculate the likelihood of the recent observed data given each scenario trajectory. 3. **Filtering/Pruning:** Apply an adaptive rejection scheme to discard scenario runs with likelihood scores below a defined threshold (or keep the top matches). 4. **Aggregation:** Construct the final forecast distribution by aggregating the remaining, empirically validated scenario trajectories.” UMass-gbqr “A gradient boosting quantile regression model that utilizes extensive feature engineering, including lagged values and time-index covariates, to minimize the pinball loss function. 1. **Variance Stabilization:** Apply a fourth-root transformation () to the raw hospitalization counts. 2. **Feature Engineering:** Construct a feature matrix including autoregressive lags (), demographic constants, cyclical calendar features (sin/cos of day-of-year), and the integer forecast horizon. 3. **Ensemble Training:** Initialize a Gradient Boosting regressor (e.g., LightGBM) with the ‘quantile‘ objective. Train separate estimators for each required quantile by minimizing the Pinball Loss. 4. **Prediction:** Generate transformed quantile forecasts
65
Prospective multi-pathogen disease forecasting using autonomous LLM-guided tree search
for the horizon. 5. **Inverse Transformation:** Raise predicted values to the fourth power () to return to the natural count scale.” CADPH-CovidCAT_Ensemble “A multi-method ensemble that aggregates forecasts from linear autoregressive processes, exponential smoothing, and machine learning regressors trained on multiple surveillance data streams. 1. **Data Stream Processing:** Ingest and clean historical time series from three distinct sources: Clinical Laboratory data, Public Health Laboratory data, and Hospitalization census. 2. **Time-Series Modeling:** For each data stream, fit an Autoregressive Integrated Moving Average (ARIMA) process to capture stochastic trends and an Exponential Smoothing (ETS) state-space model to capture error, trend, and seasonality components. 3. **Feature Construction:** Create a feature matrix for predictive modeling containing current values, lagged observations (e.g., ), and interaction terms between the different surveillance streams. 4. **Predictive Modeling:** Train a Linear Regression model to map features to future hospitalizations using linear coefficients. Simultaneously, train an ensemble of decision trees (Random Forest) to capture non-linear dependencies and feature interactions. 5. **Ensemble Aggregation:** Align the forecast trajectories from all component models (ARIMA, ETS, Linear, Tree Ensemble) on a temporal grid. 6. **Consensus Generation:** Compute the central tendency (median) of the predictions at each time step to produce the final point forecast and derive quantiles from the ensemble distribution.” UMass-ar6_pooled “A pooled autoregressive model of order 6 where coefficients are shared globally across all locations, while variance is estimated locally. 1. **Transformation:** Apply a variance-stabilizing fourth-root transformation to all data. 2. **Global Estimation:** Stack data from all locations and estimate a single set of AR(6) coefficients minimizing the global squared error. 3. **Local Variance Estimation:** Calculate residuals for each specific location using the global coefficients and compute location-specific error variance (). 4. **Recursive Forecasting:** Iterate the AR equation forward, adding noise sampled from for each step. 5. **Inversion:** Apply the inverse transformation to the generated trajectories.” CFA-EpiAutoGP “A non-parametric statistical model using Gaussian Process regression to infer the underlying structure of the epidemic curve without rigid mechanistic assumptions. 1. **Kernel Specification:** Define a covariance kernel function (e.g., Matérn or Spectral Mixture) that dictates the smoothness and periodicity of the modeled function. 2. **Hyperparameter Optimization:** Maximize the marginal likelihood of the Gaussian Process given the historical hospitalization data to solve for optimal kernel hyperparameters (length-scale, amplitude). 3. **Posterior Conditioning:** Condition the GP prior on the training data to derive the posterior distribution of functions. 4. **Function Sampling:** Draw function samples from the posterior predictive distribution for the forecast horizon. 5. **Transformation:** If data was log-transformed, exponentiate the samples to return to the count scale and compute quantiles.”
66
Prospective multi-pathogen disease forecasting using autonomous LLM-guided tree search
Metaculus-cp “An aggregation model that compiles probabilistic forecasts submitted by a crowd of human forecasters and computes a consensus distribution. 1. **Elicitation:** Distribute forecasting questions (targets and horizons) to the forecaster community. 2. **Collection:** Aggregate individual user submissions, which are full probability distributions. 3. **Consensus Algorithm:** Apply a weighting algorithm (e.g., RecencyWeighted Median or accuracy-weighted "Metaculus Prediction") to combine individual distributions into a single consensus PDF. 4. **Extraction:** Derive the specific reporting quantiles directly from the consensus Probability Density Function.” UGA_flucast-INFLAenza “A spatial time-series model using Integrated Nested Laplace Approximation (INLA) to perform Bayesian inference on a latent field with spatial (BYM) and temporal (Random Walk) components. 1. **Graph Definition:** Construct the adjacency matrix representing the neighbor relationships between all US states. 2. **Model Specification:** Define the additive linear predictor, including a Random Walk component for temporal trends, a Besag-York-Mollié (BYM) component for spatial smoothing, and harmonic terms for seasonality. 3. **Bayesian Inference:** Use the ‘R-INLA‘ package to deterministically approximate the posterior marginals of the latent field and hyperparameters. 4. **Sampling:** Draw joint samples from the approximated posterior distribution. 5. **Forecast Construction:** Project the sampled latent field components forward, apply the inverse link function (exponentiation), and compute quantiles.” CMU-TimeSeries “A CDF-space-averaged ensemble combining a population-normalized autoregressive model, a covariate-assisted quantile regression model, and a regime-switching forecaster based on epidemic direction. 1. **Component 1 (Point AR):** Normalize hospitalization data by population. Fit an autoregressive model on a rolling 28-day window. Generate distributions by resampling from studentized residuals of the fit. 2. **Component 2 (Quantile Regression):** Construct a feature set using lagged hospitalization counts and confirmed case counts. Train separate linear quantile regression models for each target quantile ($ \tau \in {0.025,..., 0.975} $) on a 28-day count-scale window. Apply post-hoc sorting to ensure quantile monotonicity. 3. **Component 3 (Directional Stratification):** Calculate trend metrics to classify the epidemic trajectory into discrete states (e.g., "up", "steady", "down"). Select a pre-trained quantile regression model specific to the identified directional class to generate conditional forecasts. 4. **Ensemble Integration:** Convert the output quantiles of all three components into Cumulative Distribution Functions (CDFs). 5. **Probabilistic Averaging:** Compute the arithmetic mean of the probabilities across the CDFs (linear pooling). Invert the averaged CDF to extract the final submission quantiles.” JHU_CSSE-CSSE_Ensemble “A hierarchical ensemble blending state-level ARIMA models with regional and national-level Long Short-Term Memory (LSTM) networks. 1. **Hierarchical Aggregation:** Aggregate state-level data to form Regional (HHS regions) and National time series. 2. **Local Modeling:** Fit independent ARIMA models to each state’s series to generate local base forecasts. 3. **Regional/National Modeling:** Train LSTM networks on the regional and national datasets to capture broader spatial dynamics. 4. **Disaggregation:** Map the regional and national forecasts down to the state level (e.g., via
67
Prospective multi-pathogen disease forecasting using autonomous LLM-guided tree search
population weighting or historical ratios). 5. **Ensemble Weighting:** Compute the final forecast for each state by taking a weighted combination of its specific ARIMA result and the disaggregated LSTM predictions.” CMU-climate_baseline “A climatological baseline model that projects future trajectories based on the historical mean and variance of the target variable for the specific time of year. 1. **Temporal Alignment:** Aggregate all available historical hospitalization data and align observations by a seasonal index (e.g., Day of Year or Epidemiological Week). 2. **Statistical Profiling:** For each seasonal index in the forecast horizon, calculate the central tendency (mean/median) and dispersion (standard deviation/interquartile range) of the historical observations. 3. **Projection:** Populate the forecast trajectory using the calculated historical central tendencies for the corresponding future dates. 4. **Uncertainty Construction:** Generate prediction intervals by applying the historical dispersion parameters to the point forecast, assuming a distribution shape consistent with the historical variance.” CEPH-Rtrend_covid “A semi-mechanistic renewal model that estimates the posterior distribution of the effective reproduction number (Rt) via MCMC and projects it forward based on recent trends. 1. **Signal Processing:** Apply a low-pass filter to smooth the hospitalization time series and interpolate to a daily resolution to match the generation interval time-step. 2. **Bayesian Inference:** Define a likelihood function dependent on the convolution of past incidence, the generation interval, and . Use the Metropolis-Hastings MCMC algorithm to sample the posterior distribution of given the data and an informed prior. 3. **Trend Estimation:** Extract estimates from the most recent 21-day window. Fit a trend line (e.g., linear regression on log space) to this window to characterize the current trajectory. 4. **Forward Simulation:** Project into the forecast horizon by continuing the identified trend. 5. **Incidence Calculation:** Input the projected trajectory into the renewal equation to iteratively simulate future daily hospitalizations. 6. **Quantile Extraction:** Aggregate the simulated trajectories to compute the required probability quantiles.”
Source: Prior Season (Double Adapted / Deep Research) DeepResearch_CounterfactualSimulation “Implement an unconditional uncertainty quantification module using Monte Carlo simulation. **Scenario Definition**: Define a set of key uncertain future drivers (e.g., ’New_Variant_Emerge’, ’Holiday_Mobility_Change’) and assign a probability distribution to each. **Simulation Loop**: Execute N (e.g., 10,000) simulation runs. In each run, sample a complete future trajectory for all defined drivers, creating one plausible future scenario. Feed this specific vector of future covariates into the pre-trained Causal Spatio-Temporal Graph Transformer (Pillar 1) to generate a single median forecast trajectory conditioned on that scenario. **Aggregation**: Collect the N median forecast trajectories. The final submitted forecast quantiles (0.025, 0.5, 0.975) must be the empirical quantiles calculated across this aggregated distribution of N trajectories for each forecast date.”
68
Prospective multi-pathogen disease forecasting using autonomous LLM-guided tree search
CMU_TimeSeries-UMass_gbqr “We have up until now done experiments with two major types of codes, that are described in detail below. PLEASE CREATE AN ALGORITHM THAT USES THE BEST PARTS OF BOTH STRATEGIES TO CREATE A HYBRID STRATEGY THAT IS TRULY WONDERFUL AND SCORES HIGHER THAN EITHER OF THE INDIVIDUAL STRATEGIES. These two code solutions tackle the problem of forecasting COVID-19 hospitalizations using distinct methodological approaches: Code 1 employs a classical statistical time-series method (Quantile Regression with an Autoregressive structure), while Code 2 leverages a modern machine learning ensemble technique (Gradient Boosting Quantile Regression with LightGBM). Here’s a comparison highlighting their main differing principles: — ### CODE 1: Quantile Regression with AR Lags (Statsmodels ‘QuantReg‘) **Main Principles:** 1. **Classical Statistical Time Series (AR Model):** * **Autoregressive (AR) Structure:** The core idea is that future values are a linear function of past values of the series itself. The model explicitly creates lagged features (‘lag_1‘, ‘lag_2‘, etc.) directly from the target variable’s history. * **Iterative (Recursive) Forecasting:** For multi-step-ahead predictions (horizons 0, 1, 2, 3), the model uses its *own median prediction* for the previous time step as an input for the next. For example, to predict week 2, it uses the prediction for week 1 as a feature. This can lead to compounding errors over longer horizons. * **Linear Model (Quantile Regression):** ‘statsmodels.regression.quantile_regression.QuantReg‘ is used. This fits a linear relationship between the lagged features and the target, but estimates different coefficients for each quantile. This is a direct approach to quantile regression, estimating the conditional quantile function. 2. **Explicit Data Normalization and Smoothing:** * **Population Normalization:** The ’Total COVID-19 Admissions’ are explicitly divided by the ‘population‘ to create a ‘per-capita‘ rate (‘admissions_norm‘). This is done to make the series comparable across different-sized states. The model is trained on these normalized rates, and predictions are denormalized at the end. * **Rolling Mean Smoothing:** A ‘smoothing_window_weeks‘ (e.g., 1 week) rolling mean is applied to the normalized admissions. This reduces noise in the time series, making the AR model’s coefficients more stable and the underlying trend more apparent, which can be beneficial for linear models that are sensitive to noise. 3. **Minimal Feature Set:** * The features for the ‘QuantReg‘ model are primarily the **lagged values of the smoothed, normalized admissions** and a constant (intercept) term. It does not explicitly incorporate other features like calendar dates (year, month, week, day of year), or interaction terms. The model relies almost entirely on the historical patterns of the admissions data itself. 4. **Training Data Windowing:** * The model explicitly filters the training data to a recent ‘train_window_days‘ (e.g., 28 days = 4 weeks). This assumes that only recent historical data is relevant for predicting future trends, and it helps manage computational complexity by limiting the size of the dataset for ‘QuantReg‘. 5. **Robustness and Post-processing:** * Includes ‘try-except‘ blocks for model fitting, handling cases where specific quantile models might fail. * Applies critical post-processing steps: * **Non-negativity:** Ensures predicted counts are not less than zero. * **Monotonicity (Quantile Sorting):** Uses ‘np.maximum.accumulate‘ to ensure that predicted quantiles are strictly non-decreasing (e.g., ‘quantile_0.01 <= quantile_0.025 <= ...‘). This is essential for valid probabilistic forecasts. * **Rounding:** Rounds predictions to the nearest integer, as admissions are counts. — ### CODE 2: Gradient Boosting Quantile Regression (LightGBM) **Main Principles:**
69
Prospective multi-pathogen disease forecasting using autonomous LLM-guided tree search
1. **Machine Learning Ensemble (Gradient Boosting Trees):** * **LightGBM (LGBMRegressor):** This is an implementation of Gradient Boosting Decision Trees. It’s a powerful non-linear, non-parametric ensemble method. Instead of fitting a single linear model, it builds a sequence of weak prediction models (decision trees), with each new tree correcting the errors of the previous ones. * **Direct Quantile Optimization:** LightGBM supports a ’quantile’ objective function, meaning each model directly optimizes the specified quantile loss, similar to ‘QuantReg‘ but through an ensemble tree structure. * **Direct Multi-Horizon Forecasting:** Unlike the iterative AR approach, this code frames the problem as predicting a future value given its ‘reference_date‘, ‘target_end_date‘, ‘location‘, and a specific ‘horizon‘. The ‘horizon‘ itself is included as an explicit feature, allowing the model to learn different patterns or biases for 1-week ahead vs. 4-week ahead forecasts independently, without compounding errors from previous predictions. 2. **Rich Feature Engineering:** * **Datetime Features:** Extensively extracts features from ‘target_end_date‘ and ‘reference_date‘, including year, month, week of year, day of year, and cyclical ‘sin/cos‘ transformations for week of year. This allows the model to capture complex seasonal patterns and trends that are not directly captured by simple lags. * **Population as a Covariate:** ‘log_population‘ is used as a direct feature in the model, along with an ‘interaction term‘ (‘log_population_horizon_interaction‘). This allows the model to learn how admission counts scale with population and how this scaling might change with the forecast horizon. * **Lagged Features of Raw Data:** Lags are also included (‘lag_admissions_Xw‘), but these are typically applied to the **raw, unsmoothed** admission counts. LightGBM’s tree-based nature is more robust to noise and can learn complex relationships directly from the raw data, often reducing the need for explicit smoothing. 3. **Data Structuring for ML Forecasting:** * **Extended Training Data:** A key difference is how the training data is prepared. For each historical observation, multiple training samples are generated by associating it with different ‘reference_date‘ and ‘horizon‘ pairs that *could have produced* that observation. This converts the time-series problem into a tabular supervised learning problem where each row is a (location, reference_date, horizon) tuple with its corresponding ‘Total COVID-19 Admissions‘. This enables the direct multi-horizon forecasting approach. 4. **Categorical Feature Handling:** * Explicitly identifies and converts features like ’location’ and ’horizon’ to ’category’ dtype. LightGBM is highly optimized to handle categorical features directly, often outperforming one-hot encoding for tree-based models. 5. **Robustness and Post-processing:** * Similar to Code 1, it implements crucial postprocessing: * **Non-negativity:** Ensures predictions are not less than zero. * **Rounding:** Rounds predictions to the nearest integer. * **Monotonicity (Quantile Sorting):** Uses ‘np.sort(..., axis=1)‘ to enforce that predicted quantiles are non-decreasing for each forecast. * Includes a ‘_return_empty_predictions‘ helper for robust error handling, returning zero-filled DataFrames if data issues prevent training or prediction. — ### Summary of Main Principles Differences: | Feature/Principle | CODE 1: QuantReg AR Model | CODE 2: LightGBM Quantile Regression | | :———————- | :——————————————————– | :—————————— ————————– | | **Core Algorithm** | **Statistical (Quantile Regression)**; linear, parametric. | **Machine Learning (Gradient Boosting Trees)**; non-linear, non-parametric ensemble. | | **Forecasting Method** | **Iterative/Recursive AR**: Predictions fed back as lags. | **Direct Multi-Horizon**: ‘horizon‘ is a feature; each forecast independent. | | **Target Preprocessing**| **Normalized** by population, then **smoothed**. | Uses **raw admissions**; population as a covariate. | | **Lagged Features** | Lags of ‘smoothed, normalized‘ target variable. | Lags of ‘raw‘ target variable (and many other features). | | **Feature Engineering** | Minimalist: primarily ‘lags‘ + intercept. | Rich: ‘lags‘, ‘calendar features‘, ‘population‘, ‘interaction terms‘. | | **Data Framing** | Traditional time series. | Transforms time series into tabular supervised learning via **extended training data**.
70
Prospective multi-pathogen disease forecasting using autonomous LLM-guided tree search
| | **Flexibility** | Lower (linear assumptions). | Higher (captures complex, non-linear relationships and interactions). | | **Interpretability** | Higher (linear coefficients). | Lower (ensemble "black-box"). | | **Computation** | Potentially slower due to ‘QuantReg‘ optimization and iterative forecasting. | Generally faster and more scalable for large datasets due to optimized GBDT implementation. | In essence, Code 1 is a more traditional, parsimonious time-series model well-suited for understanding the direct linear impact of past values on future ones. Code 2 is a more modern, flexible machine learning approach that can capture complex, non-linear relationships by integrating a wide array of features and learning patterns across different forecast horizons directly. The choice between them often depends on the specific problem characteristics, required performance, and interpretability needs. For competitive forecasting, models like LightGBM (Code 2) often achieve superior performance due to their flexibility.” DeepResearch_RegimeSwitchingDetection “Implement an automated regime-switching detection and response system. **Detection Module**: Apply a Bayesian Change Point (BCP) detection algorithm continuously to a sliding window of key indicator time series (e.g., hospitalization growth rate, wastewater momentum). **Trigger Mechanism**: If the BCP algorithm identifies a change point with a posterior probability exceeding a predefined threshold (e.g., 0.95), an automated ’regime switch’ alert must be triggered. **Response Protocol**: This alert must programmatically initiate a cascade of four actions: 1. force immediate retraining of core predictive models; 2. apply a temporal decay weight to the loss function during retraining to down-weight pre-switch data; 3. signal the meta-ensemble to adjust weights; 4. apply a multiplicative factor to widen final prediction intervals.” CEPH_Rtrend_covid-CMU_climate_baseline “We have up until now done experiments with two major types of codes, that are described in detail below. PLEASE CREATE AN ALGORITHM THAT USES THE BEST PARTS OF BOTH STRATEGIES TO CREATE A HYBRID STRATEGY THAT IS TRULY WONDERFUL AND SCORES HIGHER THAN EITHER OF THE INDIVIDUAL STRATEGIES. These two code solutions represent fundamentally different approaches to time series forecasting, particularly in the context of epidemiological modeling. Here’s a comparison explaining the main principles that differ between them: — ### CODE 1: Mechanistic/Epidemiological Renewal Equation Model with MCMC **Main Principles:** 1. **Mechanistic Modeling (Renewal Equation):** * **Principle:** This code explicitly models the underlying epidemiological process of disease transmission using the renewal equation. It assumes that the number of new hospitalizations at time ‘t‘ (‘I_t‘) is a function of the effective reproduction number (‘R_t‘) and the effective infectious pressure from past hospitalizations (‘sum(I_{t-s} * g_s)‘), where ‘g_s‘ is the generation interval distribution. * **Implementation:** It calculates a ‘generation_interval_pmf‘ (probability mass function), which describes how infections are distributed in time from a primary case to secondary cases. It then uses this PMF to relate past cases to current infections. 2. **Bayesian Inference for ‘R_t‘ (MCMC):** * **Principle:** The core of the model is to estimate the *effective reproduction number (R_t)*, which represents the average number of secondary cases generated by one infected individual. It does this using a Bayesian approach, specifically Metropolis-Hastings Markov Chain Monte Carlo (MCMC). MCMC allows for the sampling from the posterior distribution of ‘R_t‘, given the observed hospitalization data
71
Prospective multi-pathogen disease forecasting using autonomous LLM-guided tree search
and prior beliefs about ‘R_t‘. * **Implementation:** The ‘run_mcmc_rt_estimation‘ function defines a likelihood (Poisson distribution for observed cases given expected cases) and a prior (Gamma distribution for ‘R_t‘). The MCMC algorithm then iteratively proposes new ‘R_t‘ values, accepts or rejects them based on the Metropolis-Hastings criterion, thereby generating samples from ‘R_t‘’s posterior distribution. This provides a measure of uncertainty in ‘R_t‘. 3. **Stochastic Simulation for Forecasting:** * **Principle:** Once a distribution of ‘R_t‘ values is obtained, future hospitalizations are forecasted by running multiple Monte Carlo simulations. Each simulation uses a sampled ‘R_t‘ value (which itself can evolve stochastically over the forecast horizon) and the renewal equation to project daily hospitalizations forward in time. Randomness (e.g., Poisson noise for new cases, Gaussian noise for ‘R_t‘ drift) is incorporated at each step to capture inherent variability. * **Implementation:** The ‘forecast_hospitalizations‘ function takes ‘rt_samples‘ and the ‘gi_pmf‘ to simulate many trajectories. It includes parameters for ‘rt_forecast_noise_sd‘ and ‘drift_alpha‘ (Rt tending towards 1.0, representing an endemic state) to model ‘R_t‘’s future evolution. 4. **Signal Processing for Data Smoothing:** * **Principle:** Real-world epidemiological data can be noisy. To extract the underlying trend and make ‘R_t‘ estimation more robust, the code applies a lowpass filter. * **Implementation:** The ‘lowpass_filter_data‘ function uses a Butterworth filter (‘scipy.signal.butter‘, ‘filtfilt‘) to smooth the weekly hospitalization time series, removing high-frequency noise while preserving the general shape of the epidemic curve. 5. **Daily Resolution Interpolation:** * **Principle:** The renewal equation typically operates on a daily time step, as infections and their spread are continuous processes. Weekly data needs to be disaggregated to fit this model. * **Implementation:** The ‘interpolate_to_daily‘ function converts weekly sums into average daily values for each week, providing a daily time series for the renewal equation and ‘R_t‘ estimation. — ### CODE 2: Climatological/Historical Averaging Model **Main Principles:** 1. **Empirical/Climatological Forecasting (Historical Analogs):** * **Principle:** This code operates on the assumption that future patterns for a given week of the year will resemble past patterns for that same week of the year. It does not attempt to model disease transmission dynamics explicitly. Instead, it directly uses historical observations from similar time periods to form its predictions. * **Implementation:** The model identifies the "week of year" for each target forecast date and then looks back at all historical data (from ‘min_year‘ onwards) for that specific week (and a surrounding ‘window_size‘ of weeks). 2. **Ensemble of Quantiles (Geo-Specific and Geo-Aggregated):** * **Principle:** To improve robustness and incorporate different levels of historical information, the model generates two sets of quantiles and then averages them. * **Geo-Specific:** Uses historical data *only from the target location* for the relevant weeks. This captures local patterns. * **GeoAggregated:** Uses historical data *from all locations* for the relevant weeks, normalized by population (‘case_rate‘). This helps provide a more robust estimate, especially for locations with sparse data, by leveraging broader trends. * **Implementation:** It calculates ‘geo_specific_quantiles‘ and ‘geo_aggregated_quantiles‘ independently and then combines them (e.g., simple average). This implicitly handles uncertainty by considering the range of past observations. 3. **Direct Quantile Estimation:** * **Principle:** Instead of running stochastic simulations to generate a distribution of outcomes, this model directly calculates empirical quantiles from the aggregated historical data points. * **Implementation:** It collects all relevant historical hospitalization counts (or rates) within the defined time window and then uses ‘numpy.percentile‘ to find the values corresponding to the required quantile levels. 4. **Population Normalization:** * **Principle:** When combining data across different locations (for geo-aggregated quantiles), it’s important to account for varying population
72
Prospective multi-pathogen disease forecasting using autonomous LLM-guided tree search
sizes to ensure fair comparison and scaling. * **Implementation:** It calculates a ‘case_rate‘ (hospitalizations per 100,000 population) in the training data for geo-aggregated analysis. When converting aggregated rate quantiles back to counts for the target location, it scales them by the target location’s population. 5. **Smoothing Factor:** * **Principle:** For highly volatile or sparse data (e.g., low incidence periods), direct historical quantiles might be noisy or result in overly high predictions if a few historical outliers exist. A smoothing factor allows for blending the historical estimate with a default value (usually 0) to pull predictions towards zero, making the model more conservative in low-incidence scenarios. * **Implementation:** ‘final_quantiles = final_quantiles * (1 smoothing_factor) + (0 * smoothing_factor)‘. — ### Summary of Key Differences: | Feature | CODE 1 (Mechanistic/Renewal Equation) | CODE 2 (Climatological/Historical Averaging) | | :———————– | :—————————————————————————| :—————————————————————————— | | **Core Model Type** | **Mechanistic/Epidemiological**: Models disease transmission dynamics (Rt, GI). | **Empirical/Climatological**: Relies on historical patterns and seasonality. | | **Forecasting Mechanism**| Estimates ‘R_t‘ and simulates future cases using the renewal equation. | Directly computes quantiles from historical observations in similar time periods. | | **Temporal Focus** | Primarily uses **recent history** to infer current ‘R_t‘ and projects forward. | Uses **all relevant historical data** from the *same week of the year* (seasonality). | | **Uncertainty Gen.** | Achieved via **Monte Carlo simulations** (sampling ‘R_t‘, Poisson noise). | Derived from the **empirical distribution** of historical observations. | | **Data Preprocessing** | Extensive: Lowpass filtering, daily interpolation. | Simpler: Date parsing, week/year extraction, population normalization. | | **Adaptability** | Can adapt to *recent changes* in transmission if reflected in recent data. | Less responsive to *novel dynamics* not seen historically; robust to sparse recent data. | | **Interpretability** | Provides insights into epidemiological parameters like ‘R_t‘. | Simpler, "what happened last year (or similar years) at this time?". | | **Computational Cost** | Higher due to MCMC and Monte Carlo simulations. | Lower, mainly data aggregation and percentile calculations. | | **Key Assumption** | Renewal equation holds, ‘R_t‘ dynamics are predictable. | Future resembles past seasonality, historical variability represents future uncertainty. | In essence, Code 1 attempts to understand and predict *how* the disease spreads, making it a more "scientific" or "mechanistic" model. Code 2, on the other hand, is a more "statistical" or "data-driven" model that simply observes and projects what has happened before, assuming historical patterns will repeat.” CMU_climate_baseline-UMass_ar6_pooled “We have up until now done experiments with two major types of codes, that are described in detail below. We also have created a hybrid model that performs well on a short period of time. PLEASE MODIFY THE ALGORITHM THAT USES THE BEST PARTS OF BOTH STRATEGIES TO CREATE A HYBRID STRATEGY THAT IS TRULY WONDERFUL AND SCORES HIGHER THAN EITHER OF THE INDIVIDUAL STRATEGIES ON THE ENTIRE SEASON. The two code solutions offer fundamentally different approaches to forecasting COVID-19 hospital admissions, despite addressing the same problem of generating probabilistic quantile forecasts for multiple US states. Here’s a comparison highlighting their main differing principles: — ### Code 1: Climatological Model **Main Principles:** 1. **Modeling Paradigm: Non-Parametric, Historical Quantile Lookup** * This model is a *climatological* or *analog-based* approach. It does not fit a traditional statistical model
73
Prospective multi-pathogen disease forecasting using autonomous LLM-guided tree search
(like regression or ARIMA). Instead, it predicts future values by looking at actual historical observations from similar time periods (i.e., the same week of the year). * Quantiles are derived directly by calculating percentiles of these historical observations. This makes it non-parametric, meaning it doesn’t assume a specific underlying distribution for the data. 2. **Temporal Dependencies & Seasonality: Explicit Seasonal Focus** * The core assumption is strong *seasonality*. Forecasts for a given ‘target_week_num‘ are based on data from that same week (and a surrounding ‘window_size‘ of weeks) across *all past years* included in the ‘min_year‘ filter. * It explicitly handles year wrap-around for week numbers (e.g., week 53 to week 1). * It does *not* explicitly model short-term autocorrelation (how last week’s value directly influences this week’s value, independent of season). 3. **Cross-Location Information Sharing (Borrowing Strength): Ensemble Approach** * It uses an *ensemble* of two types of historical data: * **Geo-specific:** Uses data *only from the target location* for the relevant weeks. This captures local patterns. * **Geoaggregated:** Uses data from *all locations* (normalized by population as ’case_rate’) for the relevant weeks. This "borrows strength" across states, assuming general disease trends are somewhat shared, especially for locations with sparse data. * The final prediction is a simple average of these two sets of quantiles, balancing local and general patterns. 4. **Data Transformation: Normalization (Rate-based)** * It normalizes ‘Total COVID-19 Admissions‘ to ‘case_rate‘ (per 100,000 population) *only for the geo-aggregated component*. This is a normalization step to make data comparable across different population sizes, not a statistical transformation to improve model linearity or normality. * The geo-specific quantiles are calculated directly on raw admission counts. 5. **Robustness & Fallbacks:** * Uses ‘min_samples_geo_specific‘ and ‘min_samples_geo_aggregated‘ to ensure a minimum number of historical observations before calculating quantiles. If not met, it sets quantiles to NaN, which are then handled by the ensemble logic (favoring available quantiles, or falling back to 0). * A ‘smoothing_factor‘ is applied to pull predictions towards 0, which can be useful for sparse or declining case counts, acting as a form of regularization. * Enforces non-negativity and monotonicity of quantiles. — ### Code 2: Autoregressive (AR) Model **Main Principles:** 1. **Modeling Paradigm: Parametric, Time-Series Forecasting** * This model uses a *parametric statistical time-series* approach, specifically an Autoregressive (AR) model of ‘AR_ORDER‘. It explicitly models the current value as a linear function of its own past values (lags). * Quantiles are derived by assuming a Gaussian (Normal) distribution for the transformed data (or its errors) around the predicted mean, scaled by a standard deviation. 2. **Temporal Dependencies & Seasonality: Autocorrelation Focus** * The core mechanism for capturing temporal patterns is ‘AR_ORDER‘ lags. This directly models *autocorrelation* – how past values influence current values. It’s excellent for capturing short-to-medium-term persistence, trends, and oscillations inherent in the time series itself. * It does *not* explicitly model seasonality (e.g., there’s no specific 52-week lag or seasonal differencing). Any seasonal patterns would have to be implicitly captured by the short AR lags, or might be missed. * It includes a robust multi-step forecasting loop, where predicted values are fed back into the lag buffer to forecast further into the future. It also explicitly handles "retrospective" forecasts (target date is within the training set). 3. **Cross-Location Information Sharing (Borrowing Strength): Shared Parameters** * It employs a *shared parameter* approach: * **Shared AR Coefficients:** A single set of AR coefficients and an intercept are estimated using Ordinary Least Squares (OLS) from the *entire pooled training dataset across all locations*. This is a strong form of borrowing strength, assuming a common underlying dynamic for admissions (after transformation) across states. * **Location-Specific Variances:** While the AR dynamics are shared, the
74
Prospective multi-pathogen disease forecasting using autonomous LLM-guided tree search
variance of the model’s errors (residuals) is estimated *separately for each location*. This allows the model to capture differing levels of uncertainty or volatility unique to each state. 4. **Data Transformation: Fourth-Root Transform** * A significant difference is the application of a **fourth-root transformation** ‘(x + epsilon_transform)^(1/4)‘ to the raw admission counts. This is done to: * Stabilize the variance of the data (make it less dependent on the mean). * Make the data distribution more symmetric and closer to Gaussian, which is an assumption underlying OLS and subsequent normal quantile generation. * Handle zero values by adding a small ‘epsilon‘. * An inverse transformation ‘(y_transformed)^4 epsilon_transform‘ is applied before outputting the final integer counts. 5. **Robustness & Fallbacks:** * Uses ‘epsilon_transform‘ and ‘epsilon_std‘ for numerical stability and to prevent zero-width prediction intervals. * Employs global mean/variance as fallback parameters if OLS fitting fails (e.g., due to insufficient data for a specific location or rank deficiency in the design matrix), or for locations entirely new to the training data. * Pads the lag buffer with the shared intercept if a location has insufficient historical data to form the required lags. * Enforces non-negativity and monotonicity of quantiles after inverse transformation. — ### Summary of Key Differences: | Feature | Code 1: Climatological Model | Code 2: AR Model | | :———————— | :——— ————————————————— | :——————————————————————| | **Modeling Paradigm** | Non-parametric, historical lookup/averaging | Parametric, statistical time-series (ARIMA subclass) | | **Primary Driver** | Explicit seasonality (week of year) | Autocorrelation (lags), capturing short/medium-term trends/persistence | | **Data Transformation** | Normalization to "rate" for aggregation; no statistical transform of target variable | Fourth-root transformation of target variable to stabilize variance and normalize errors | | **Info Sharing (States)** | Ensemble of geo-specific and geo-aggregated patterns | Shared AR coefficients across states; location-specific error variances | | **Quantile Method** | Direct ‘np.percentile‘ from historical data | Assumed Gaussian distribution; ‘norm.ppf‘ on transformed data | | **Multi-step Forecasting**| Indirect (looks up historical data for target week) | Explicit recursive prediction (feeding forecasts back into lags) | | **Assumption on Data** | Patterns recur yearly for a given week | Linear relationship with past values in transformed space; Gaussian errors | | **Handling New Locations**| Primarily falls back to geo-aggregated (population-normalized) data, or 0 if truly sparse | Uses shared AR coefficients and global/average variance as strong fallback | In essence, Code 1 is a "memory-based" model that looks back to what happened *at this time of year* in the past. Code 2 is a "momentum-based" model that looks at the *recent trend and past values* to predict the next steps, assuming a consistent underlying dynamic across states. The choice between them (or combining them) often depends on which characteristics of the target variable (COVID-19 admissions) are most dominant and predictable.”
75
Prospective multi-pathogen disease forecasting using autonomous LLM-guided tree search
D.3. Method Descriptions for RSV (RSVHub)
76
Prospective multi-pathogen disease forecasting using autonomous LLM-guided tree search
References [1] Eser Aygün, Anastasiya Belyaeva, Gheorghe Comanici, Marc Coram, Hao Cui, Jake Garrison, Renee Johnston Anton Kast, Cory Y. McLean, Peter Norgaard, Zahra Shamsi, David Smalling, James Thompson, Subhashini Venugopalan, Brian P. Williams, Chujun He, Sarah Martinson, Martyna Plomecka, Lai Wei, Yuchen Zhou, Qian-Ze Zhu, Matthew Abraham, Erica Brand, Anna Bulanova, Jeffrey A. Cardille, Chris Co, Scott Ellsworth, Grace Joseph, Malcolm Kane, Ryan Krueger, Johan Kartiwa, Dan Liebling, Jan-Matthis Lueckmann, Paul Raccuglia, Xuefei, Wang, Katherine Chou, James Manyika, Yossi Matias, John C. Platt, Lizzie Dorfman, Shibl Mourad, and Michael P. Brenner. An AI system to help scientists write expert-level empirical software, September 2025. arXiv:2509.06503 [cs] version: 1. [2] Cécile Viboud, Kaiyuan Sun, Robert Gaffey, Marco Ajelli, Laura Fumanelli, Stefano Merler, Qian Zhang, Gerardo Chowell, Lone Simonsen, Alessandro Vespignani, and RAPIDD Ebola Forecasting Challenge group. The RAPIDD ebola forecasting challenge: Synthesis and lessons learnt. Epidemics, 22:13–21, March 2018. [3] Sara Y. Del Valle, Benjamin H. McMahon, Jason Asher, Richard Hatchett, Joceline C. Lega, Heidi E. Brown, Mark E. Leany, Yannis Pantazis, David J. Roberts, Sean Moore, A Townsend Peterson, Luis E. Escobar, Huijie Qiao, Nicholas W. Hengartner, and Harshini Mukundan. Summary results of the 2014-2015 DARPA Chikungunya challenge. BMC Infectious Diseases, 18(1):245, May 2018. [4] Michael A. Johansson, Karyn M. Apfeldorf, Scott Dobson, Jason Devita, Anna L. Buczak, Benjamin Baugher, Linda J. Moniz, Thomas Bagley, Steven M. Babin, Erhan Guven, Teresa K. Yamana, Jeffrey Shaman, Terry Moschou, Nick Lothian, Aaron Lane, Grant Osborne, Gao Jiang, Logan C. Brooks, David C. Farrow, Sangwon Hyun, Ryan J. Tibshirani, Roni Rosenfeld, Justin Lessler, Nicholas G. Reich, Derek A. T. Cummings, Stephen A. Lauer, Sean M. Moore, Hannah E. Clapham, Rachel Lowe, Trevor C. Bailey, Markel García-Díez, Marilia Sá Carvalho, Xavier Rodó, Tridip Sardar, Richard Paul, Evan L. Ray, Krzysztof Sakrejda, Alexandria C. Brown, Xi Meng, Osonde Osoba, Raffaele Vardavas, David Manheim, Melinda Moore, Dhananjai M. Rao, Travis C. Porco, Sarah Ackley, Fengchen Liu, Lee Worden, Matteo Convertino, Yang Liu, Abraham Reddy, Eloy Ortiz, Jorge Rivero, Humberto Brito, Alicia Juarrero, Leah R. Johnson, Robert B. Gramacy, Jeremy M. Cohen, Erin A. Mordecai, Courtney C. Murdock, Jason R. Rohr, Sadie J. Ryan, Anna M. Stewart-Ibarra, Daniel P. Weikel, Antarpreet Jutla, Rakibul Khan, Marissa Poultney, Rita R. Colwell, Brenda Rivera-García, Christopher M. Barker, Jesse E. Bell, Matthew Biggerstaff, David Swerdlow, Luis Mier-y Teran-Romero, Brett M. Forshey, Juli Trtanj, Jason Asher, Matt Clay, Harold S. Margolis, Andrew M. Hebbeler, Dylan George, and Jean-Paul Chretien. An open challenge to advance probabilistic forecasting for dengue epidemics. Proceedings of the National Academy of Sciences of the United States of America, 116(48):24268–24274, November 2019. [5] Chelsea S. Lutz, Mimi P. Huynh, Monica Schroeder, Sophia Anyatonwu, F. Scott Dahlgren, Gregory Danyluk, Danielle Fernandez, Sharon K. Greene, Nodar Kipshidze, Leann Liu, Osaro Mgbere, Lisa A. McHugh, Jennifer F. Myers, Alan Siniscalchi, Amy D. Sullivan, Nicole West, Michael A. Johansson, and Matthew Biggerstaff. Applying infectious disease forecasting to public health: a path forward using influenza forecasting examples. BMC Public Health, 19(1):1659, December 2019. [6] Nicholas G. Reich, Justin Lessler, Sebastian Funk, Cecile Viboud, Alessandro Vespignani, Ryan J. Tibshirani, Katriona Shea, Melanie Schienle, Michael C. Runge, Roni Rosenfeld, Evan L. Ray, Rene Niehus, Helen C. Johnson, Michael A. Johansson, Harry Hochheiser, Lauren Gardner, Johannes
77
Prospective multi-pathogen disease forecasting using autonomous LLM-guided tree search
Bracher, Rebecca K. Borchering, and Matthew Biggerstaff. Collaborative Hubs: Making the Most of Predictive Epidemic Modeling. American Journal of Public Health, 112(6):839–842, June 2022. [7] Jonathon Mellor, Maria L Tang, Owen Jones, Thomas Ward, Steven Riley, and Sarah R Deeny. Forecasting COVID-19, influenza, and RSV hospitalizations over winter 2023–4 in England. International Journal of Epidemiology, 54(3):dyaf066, June 2025. [8] Consortium of Infectious Disease Modeling Hubs, Melissa Kerr, Rebecca Borchering, Alvaro Castro Rivadeneira, Lucie Contamin, Sebastian Funk, Harry Hochheiser, Emily Howerton, Anna Krystalli, Li Shandross, and Nicholas G. Reich. Coordinating collaborative infectious disease modeling projects with the hubverse, April 2026. ISSN: 3067-2007 Pages: 2025.10.03.25337284. [9] Craig J. McGowan, Matthew Biggerstaff, Michael Johansson, Karyn M. Apfeldorf, Michal Ben-Nun, Logan Brooks, Matteo Convertino, Madhav Erraguntla, David C. Farrow, John Freeze, Saurav Ghosh, Sangwon Hyun, Sasikiran Kandula, Joceline Lega, Yang Liu, Nicholas Michaud, Haruka Morita, Jarad Niemi, Naren Ramakrishnan, Evan L. Ray, Nicholas G. Reich, Pete Riley, Jeffrey Shaman, Ryan Tibshirani, Alessandro Vespignani, Qian Zhang, and Carrie Reed. Collaborative efforts to forecast seasonal influenza in the United States, 2015–2016. Scientific Reports, 9(1):683, January 2019. Bandiera_abtest: a Cc_license_type: cc_by Cg_type: Nature Research Journals Number: 1 Primary_atype: Research Subject_term: Computational models;Epidemiology;Influenza virus Subject_term_id: computational-models;epidemiology;influenza-virus. [10] Nicholas G. Reich, Logan C. Brooks, Spencer J. Fox, Sasikiran Kandula, Craig J. McGowan, Evan Moore, Dave Osthus, Evan L. Ray, Abhinav Tushar, Teresa K. Yamana, Matthew Biggerstaff, Michael A. Johansson, Roni Rosenfeld, and Jeffrey Shaman. A collaborative multiyear, multimodel assessment of seasonal influenza forecasting in the United States. Proceedings of the National Academy of Sciences, 116(8):3146–3154, February 2019. [11] Estee Y. Cramer, Evan L. Ray, Velma K. Lopez, Johannes Bracher, Andrea Brennen, Alvaro J. Castro Rivadeneira, Aaron Gerding, Tilmann Gneiting, Katie H. House, Yuxin Huang, Dasuni Jayawardena, Abdul H. Kanji, Ayush Khandelwal, Khoa Le, Anja Mühlemann, Jarad Niemi, Apurv Shah, Ariane Stark, Yijin Wang, Nutcha Wattanachit, Martha W. Zorn, Youyang Gu, Sansiddh Jain, Nayana Bannur, Ayush Deva, Mihir Kulkarni, Srujana Merugu, Alpan Raval, Siddhant Shingi, Avtansh Tiwari, Jerome White, Neil F. Abernethy, Spencer Woody, Maytal Dahan, Spencer Fox, Kelly Gaither, Michael Lachmann, Lauren Ancel Meyers, James G. Scott, Mauricio Tec, Ajitesh Srivastava, Glover E. George, Jeffrey C. Cegan, Ian D. Dettwiller, William P. England, Matthew W. Farthing, Robert H. Hunter, Brandon Lafferty, Igor Linkov, Michael L. Mayo, Matthew D. Parno, Michael A. Rowland, Benjamin D. Trump, Yanli Zhang-James, Samuel Chen, Stephen V. Faraone, Jonathan Hess, Christopher P. Morley, Asif Salekin, Dongliang Wang, Sabrina M. Corsetti, Thomas M. Baer, Marisa C. Eisenberg, Karl Falb, Yitao Huang, Emily T. Martin, Ella McCauley, Robert L. Myers, Tom Schwarz, Daniel Sheldon, Graham Casey Gibson, Rose Yu, Liyao Gao, Yian Ma, Dongxia Wu, Xifeng Yan, Xiaoyong Jin, Yu-Xiang Wang, YangQuan Chen, Lihong Guo, Yanting Zhao, Quanquan Gu, Jinghui Chen, Lingxiao Wang, Pan Xu, Weitong Zhang, Difan Zou, Hannah Biegel, Joceline Lega, Steve McConnell, V. P. Nagraj, Stephanie L. Guertin, Christopher Hulme-Lowe, Stephen D. Turner, Yunfeng Shi, Xuegang Ban, Robert Walraven, Qi-Jun Hong, Stanley Kong, Axel van de Walle, James A. Turtle, Michal Ben-Nun, Steven Riley, Pete Riley, Ugur Koyluoglu, David DesRoches, Pedro Forli, Bruce Hamory, Christina Kyriakides, Helen Leis, John Milliken, Michael Moloney, James Morgan, Ninad Nirgudkar, Gokce Ozcan, Noah Piwonka, Matt Ravi, Chris Schrader, Elizabeth Shakhnovich, Daniel Siegel,
78
Prospective multi-pathogen disease forecasting using autonomous LLM-guided tree search
Ryan Spatz, Chris Stiefeling, Barrie Wilkinson, Alexander Wong, Sean Cavany, Guido España, Sean Moore, Rachel Oidtman, Alex Perkins, David Kraus, Andrea Kraus, Zhifeng Gao, Jiang Bian, Wei Cao, Juan Lavista Ferres, Chaozhuo Li, Tie-Yan Liu, Xing Xie, Shun Zhang, Shun Zheng, Alessandro Vespignani, Matteo Chinazzi, Jessica T. Davis, Kunpeng Mu, Ana Pastore y Piontti, Xinyue Xiong, Andrew Zheng, Jackie Baek, Vivek Farias, Andreea Georgescu, Retsef Levi, Deeksha Sinha, Joshua Wilde, Georgia Perakis, Mohammed Amine Bennouna, David Nze-Ndong, Divya Singhvi, Ioannis Spantidakis, Leann Thayaparan, Asterios Tsiourvas, Arnab Sarker, Ali Jadbabaie, Devavrat Shah, Nicolas Della Penna, Leo A. Celi, Saketh Sundar, Russ Wolfinger, Dave Osthus, Lauren Castro, Geoffrey Fairchild, Isaac Michaud, Dean Karlen, Matt Kinsey, Luke C. Mullany, Kaitlin Rainwater-Lovett, Lauren Shin, Katharine Tallaksen, Shelby Wilson, Elizabeth C. Lee, Juan Dent, Kyra H. Grantz, Alison L. Hill, Joshua Kaminsky, Kathryn Kaminsky, Lindsay T. Keegan, Stephen A. Lauer, Joseph C. Lemaitre, Justin Lessler, Hannah R. Meredith, Javier Perez-Saez, Sam Shah, Claire P. Smith, Shaun A. Truelove, Josh Wills, Maximilian Marshall, Lauren Gardner, Kristen Nixon, John C. Burant, Lily Wang, Lei Gao, Zhiling Gu, Myungjin Kim, Xinyi Li, Guannan Wang, Yueying Wang, Shan Yu, Robert C. Reiner, Ryan Barber, Emmanuela Gakidou, Simon I. Hay, Steve Lim, Chris Murray, David Pigott, Heidi L. Gurung, Prasith Baccam, Steven A. Stage, Bradley T. Suchoski, B. Aditya Prakash, Bijaya Adhikari, Jiaming Cui, Alexander Rodríguez, Anika Tabassum, Jiajia Xie, Pinar Keskinocak, John Asplund, Arden Baxter, Buse Eylul Oruc, Nicoleta Serban, Sercan O. Arik, Mike Dusenberry, Arkady Epshteyn, Elli Kanal, Long T. Le, Chun-Liang Li, Tomas Pfister, Dario Sava, Rajarishi Sinha, Thomas Tsai, Nate Yoder, Jinsung Yoon, Leyou Zhang, Sam Abbott, Nikos I. Bosse, Sebastian Funk, Joel Hellewell, Sophie R. Meakin, Katharine Sherratt, Mingyuan Zhou, Rahi Kalantari, Teresa K. Yamana, Sen Pei, Jeffrey Shaman, Michael L. Li, Dimitris Bertsimas, Omar Skali Lami, Saksham Soni, Hamza Tazi Bouardi, Turgay Ayer, Madeline Adee, Jagpreet Chhatwal, Ozden O. Dalgic, Mary A. Ladd, Benjamin P. Linas, Peter Mueller, Jade Xiao, Yuanjia Wang, Qinxia Wang, Shanghong Xie, Donglin Zeng, Alden Green, Jacob Bien, Logan Brooks, Addison J. Hu, Maria Jahja, Daniel McDonald, Balasubramanian Narasimhan, Collin Politsch, Samyak Rajanala, Aaron Rumack, Noah Simon, Ryan J. Tibshirani, Rob Tibshirani, Valerie Ventura, Larry Wasserman, Eamon B. O’Dea, John M. Drake, Robert Pagano, Quoc T. Tran, Lam Si Tung Ho, Huong Huynh, Jo W. Walker, Rachel B. Slayton, Michael A. Johansson, Matthew Biggerstaff, and Nicholas G. Reich. Evaluation of individual and ensemble probabilistic forecasts of COVID-19 mortality in the United States. Proceedings of the National Academy of Sciences, 119(15):e2113561119, April 2022. [12] Juliette Paireau, Alessio Andronico, Nathanaël Hozé, Maylis Layan, Pascal Crépey, Alix Roumagnac, Marc Lavielle, Pierre-Yves Boëlle, and Simon Cauchemez. An ensemble model based on early predictors to forecast COVID-19 health care demand in France. Proceedings of the National Academy of Sciences, 119(18):e2103302119, May 2022. [13] Katharine Sherratt, Hugo Gruson, Rok Grah, Helen Johnson, Rene Niehus, Bastian Prasse, Frank Sandmann, Jannik Deuschel, Daniel Wolffram, Sam Abbott, Alexander Ullrich, Graham Gibson, Evan L Ray, Nicholas G Reich, Daniel Sheldon, Yijin Wang, Nutcha Wattanachit, Lijing Wang, Jan Trnka, Guillaume Obozinski, Tao Sun, Dorina Thanou, Loic Pottier, Ekaterina Krymova, Jan H Meinke, Maria Vittoria Barbarossa, Neele Leithauser, Jan Mohring, Johanna Schneider, Jaroslaw Wlazlo, Jan Fuhrmann, Berit Lange, Isti Rodiah, Prasith Baccam, Heidi Gurung, Steven Stage, Bradley Suchoski, Jozef Budzinski, Robert Walraven, Inmaculada Villanueva, Vit Tucek, Martin Smid, Milan Zajicek, Cesar Perez Alvarez, Borja Reina, Nikos I Bosse, Sophie R Meakin, Lauren Castro, Geoffrey Fairchild, Isaac Michaud, Dave Osthus, Pierfrancesco Alaimo Di Loro, Antonello Maruotti, Veronika Eclerova, Andrea Kraus, David Kraus, Lenka Pribylova, Bertsimas Dimitris, Michael Lingzhi Li, Soni Saksham, Jonas Dehning, Sebastian Mohr, Viola Priesemann, Grzegorz
79
Prospective multi-pathogen disease forecasting using autonomous LLM-guided tree search
Redlarski, Benjamin Bejar, Giovanni Ardenghi, Nicola Parolini, Giovanni Ziarelli, Wolfgang Bock, Stefan Heyder, Thomas Hotz, David E Singh, Miguel Guzman-Merino, Jose L Aznarte, David Morina, Sergio Alonso, Enric Alvarez, Daniel Lopez, Clara Prats, Jan Pablo Burgard, Arne Rodloff, Tom Zimmermann, Alexander Kuhlmann, Janez Zibert, Fulvia Pennoni, Fabio Divino, Marti Catala, Gianfranco Lovison, Paolo Giudici, Barbara Tarantino, Francesco Bartolucci, Giovanna Jona Lasinio, Marco Mingione, Alessio Farcomeni, Ajitesh Srivastava, Pablo MonteroManso, Aniruddha Adiga, Benjamin Hurt, Bryan Lewis, Madhav Marathe, Przemyslaw Porebski, Srinivasan Venkatramanan, Rafal P Bartczuk, Filip Dreger, Anna Gambin, Krzysztof Gogolewski, Magdalena Gruziel-Slomka, Bartosz Krupa, Antoni Moszyński, Karol Niedzielewski, Jedrzej Nowosielski, Maciej Radwan, Franciszek Rakowski, Marcin Semeniuk, Ewa Szczurek, Jakub Zielinski, Jan Kisielewski, Barbara Pabjan, Kirsten Holger, Yuri Kheifetz, Markus Scholz, Biecek Przemyslaw, Marcin Bodych, Maciej Filinski, Radoslaw Idzikowski, Tyll Krueger, Tomasz Ozanski, Johannes Bracher, and Sebastian Funk. Predictive performance of multi-model ensemble forecasts of COVID-19 across European nations. eLife, 12:e81916, April 2023. [14] Nicholas G. Reich, Craig J. McGowan, Teresa K. Yamana, Abhinav Tushar, Evan L. Ray, Dave Osthus, Sasikiran Kandula, Logan C. Brooks, Willow Crawford-Crudell, Graham Casey Gibson, Evan Moore, Rebecca Silva, Matthew Biggerstaff, Michael A. Johansson, Roni Rosenfeld, and Jeffrey Shaman. Accuracy of real-time multi-model ensemble forecasts for seasonal influenza in the U.S. PLOS Computational Biology, 15(11):e1007486, November 2019. [15] Velma K. Lopez, Estee Y. Cramer, Robert Pagano, John M. Drake, Eamon B. O’Dea, Madeline Adee, Turgay Ayer, Jagpreet Chhatwal, Ozden O. Dalgic, Mary A. Ladd, Benjamin P. Linas, Peter P. Mueller, Jade Xiao, Johannes Bracher, Alvaro J. Castro Rivadeneira, Aaron Gerding, Tilmann Gneiting, Yuxin Huang, Dasuni Jayawardena, Abdul H. Kanji, Khoa Le, Anja Mühlemann, Jarad Niemi, Evan L. Ray, Ariane Stark, Yijin Wang, Nutcha Wattanachit, Martha W. Zorn, Sen Pei, Jeffrey Shaman, Teresa K. Yamana, Samuel R. Tarasewicz, Daniel J. Wilson, Sid Baccam, Heidi Gurung, Steve Stage, Brad Suchoski, Lei Gao, Zhiling Gu, Myungjin Kim, Xinyi Li, Guannan Wang, Lily Wang, Yueying Wang, Shan Yu, Lauren Gardner, Sonia Jindal, Maximilian Marshall, Kristen Nixon, Juan Dent, Alison L. Hill, Joshua Kaminsky, Elizabeth C. Lee, Joseph C. Lemaitre, Justin Lessler, Claire P. Smith, Shaun Truelove, Matt Kinsey, Luke C. Mullany, Kaitlin Rainwater-Lovett, Lauren Shin, Katharine Tallaksen, Shelby Wilson, Dean Karlen, Lauren Castro, Geoffrey Fairchild, Isaac Michaud, Dave Osthus, Jiang Bian, Wei Cao, Zhifeng Gao, Juan Lavista Ferres, Chaozhuo Li, Tie-Yan Liu, Xing Xie, Shun Zhang, Shun Zheng, Matteo Chinazzi, Jessica T. Davis, Kunpeng Mu, Ana Pastore y Piontti, Alessandro Vespignani, Xinyue Xiong, Robert Walraven, Jinghui Chen, Quanquan Gu, Lingxiao Wang, Pan Xu, Weitong Zhang, Difan Zou, Graham Casey Gibson, Daniel Sheldon, Ajitesh Srivastava, Aniruddha Adiga, Benjamin Hurt, Gursharn Kaur, Bryan Lewis, Madhav Marathe, Akhil Sai Peddireddy, Przemyslaw Porebski, Srinivasan Venkatramanan, Lijing Wang, Pragati V. Prasad, Jo W. Walker, Alexander E. Webber, Rachel B. Slayton, Matthew Biggerstaff, Nicholas G. Reich, and Michael A. Johansson. Challenges of COVID-19 Case Forecasting in the US, 2020–2021. PLOS Computational Biology, 20(5):e1011200, May 2024. [16] Sarabeth M. Mathis, Alexander E. Webber, Tomás M. León, Erin L. Murray, Monica Sun, Lauren A. White, Logan C. Brooks, Alden Green, Addison J. Hu, Roni Rosenfeld, Dmitry Shemetov, Ryan J. Tibshirani, Daniel J. McDonald, Sasikiran Kandula, Sen Pei, Rami Yaari, Teresa K. Yamana, Jeffrey Shaman, Pulak Agarwal, Srikar Balusu, Gautham Gururajan, Harshavardhan Kamarthi, B. Aditya Prakash, Rishi Raman, Zhiyuan Zhao, Alexander Rodríguez, Akilan Meiyappan, Shalina Omar, Prasith Baccam, Heidi L. Gurung, Brad T. Suchoski, Steve A. Stage, Marco Ajelli, Allisandra G. Kummer, Maria Litvinova, Paulo C. Ventura, Spencer Wadsworth, Jarad Niemi, Erica Carcelen, Alison L. Hill, Sara L. Loo, Clifton D. McKee, Koji Sato, Claire Smith, Shaun 80
Prospective multi-pathogen disease forecasting using autonomous LLM-guided tree search
Truelove, Sung-mok Jung, Joseph C. Lemaitre, Justin Lessler, Thomas McAndrew, Wenxuan Ye, Nikos Bosse, William S. Hlavacek, Yen Ting Lin, Abhishek Mallela, Graham C. Gibson, Ye Chen, Shelby M. Lamm, Jaechoul Lee, Richard G. Posner, Amanda C. Perofsky, Cécile Viboud, Leonardo Clemente, Fred Lu, Austin G. Meyer, Mauricio Santillana, Matteo Chinazzi, Jessica T. Davis, Kunpeng Mu, Ana Pastore y Piontti, Alessandro Vespignani, Xinyue Xiong, Michal Ben-Nun, Pete Riley, James Turtle, Chis Hulme-Lowe, Shakeel Jessa, V. P. Nagraj, Stephen D. Turner, Desiree Williams, Avranil Basu, John M. Drake, Spencer J. Fox, Ehsan Suez, Monica G. Cojocaru, Edward W. Thommes, Estee Y. Cramer, Aaron Gerding, Ariane Stark, Evan L. Ray, Nicholas G. Reich, Li Shandross, Nutcha Wattanachit, Yijin Wang, Martha W. Zorn, Majd Al Aawar, Ajitesh Srivastava, Lauren A. Meyers, Aniruddha Adiga, Benjamin Hurt, Gursharn Kaur, Bryan L. Lewis, Madhav Marathe, Srinivasan Venkatramanan, Patrick Butler, Andrew Farabow, Naren Ramakrishnan, Nikhil Muralidhar, Carrie Reed, Matthew Biggerstaff, and Rebecca K. Borchering. Evaluation of FluSight influenza forecasting in the 2021–22 and 2022–23 seasons with a new target laboratory-confirmed influenza hospitalizations. Nature Communications, 15(1):6289, July 2024. [17] David Silver, Julian Schrittwieser, Karen Simonyan, Ioannis Antonoglou, Aja Huang, Arthur Guez, Thomas Hubert, Lucas Baker, Matthew Lai, Adrian Bolton, et al. Mastering the game of Go without human knowledge. Nature, 550(7676):354–359, 2017. [18] Johannes Bracher, Evan L. Ray, Tilmann Gneiting, and Nicholas G. Reich. Evaluating epidemic forecasts in an interval format. PLOS Computational Biology, 17(2):e1008618, 2021. [19] Nikos I. Bosse, Sam Abbott, Anne Cori, Edwin van Leeuwen, Johannes Bracher, and Sebastian Funk. Scoring epidemiological forecasts on transformed scales. PLOS Computational Biology, 19(8):e1011393, August 2023. [20] Dave Osthus, James Gattiker, Reid Priedhorsky, and Sara Y. Del Valle. Dynamic Bayesian Influenza Forecasting in the United States with Hierarchical Discrepancy (with Discussion). Bayesian Analysis, 14(1):261–312, March 2019. [21] Dave Osthus. Fast and accurate influenza forecasting in the United States with Inferno. PLOS Computational Biology, 18(1):e1008651, January 2022. [22] Tilmann Gneiting and Adrian E. Raftery. Strictly proper scoring rules, prediction, and estimation. Journal of the American Statistical Association, 102(477):359–378, 2007. [23] Xin He, Kai Zhao, and Xiaowen Chu. Automl: A survey of the state-of-the-art. Knowledge-Based Systems, 212:106622, 2021. [24] Alyssa M. Bilinski, Joshua A. Salomon, and Laura A. Hatfield. Adaptive metrics for an evolving pandemic: A dynamic approach to area-level COVID-19 risk designations. Proceedings of the National Academy of Sciences, 120(32):e2302528120, August 2023. [25] Aaron Gerding, Nicholas G Reich, Benjamin Rogers, and Evan L Ray. Evaluating infectious disease forecasts with allocation scoring rules. Journal of the Royal Statistical Society Series A: Statistics in Society, 188(4):1299–1325, October 2025. [26] Cathal Mills, Nicholas J. Irons, Joseph L.-H. Tsui, Sarah Sparrow, Luiz M. Carvalho, Adam J. Kucharski, Oliver Ratmann, Ben Lambert, Christl A. Donnelly, and Moritz U. G. Kraemer. From metric to action: The decision value of infectious disease forecasts, March 2026. ISSN: 3067-2007 Pages: 2025.07.20.25331802.
81
Prospective multi-pathogen disease forecasting using autonomous LLM-guided tree search
[27] Gerda Claeskens, Jan R. Magnus, Andrey L. Vasnev, and Wendun Wang. The forecast combination puzzle: A simple theoretical explanation. International Journal of Forecasting, 32(3):754–762, July 2016. [28] Evan L. Ray, Logan C. Brooks, Jacob Bien, Matthew Biggerstaff, Nikos I. Bosse, Johannes Bracher, Estee Y. Cramer, Sebastian Funk, Aaron Gerding, Michael A. Johansson, Aaron Rumack, Yijin Wang, Martha Zorn, Ryan J. Tibshirani, and Nicholas G. Reich. Comparing trained and untrained probabilistic ensemble forecasts of COVID-19 cases and deaths in the United States. International Journal of Forecasting, 39(3):1366–1383, July 2023. [29] Centers for Disease Control and Prevention. Weekly hospital respiratory data (HRD) metrics by jurisdiction, National Healthcare Safety Network (NHSN), 2025. Accessed 2026. [30] Simon Pollett, Michael A Johansson, Nicholas Giangreco, Olivia Collignon, Julia J Morgan, Arthur M Hersh, Robert B McQueen, and Cécile Viboud. Recommended reporting items for epidemic forecasting and prediction research: The EPIFORGE 2020 guidelines. PLOS Medicine, 18(10):e1003793, 2021.
82