Seasonality of agricultural exposure as an important predictor of seasonal yellow fever spillover in Brazil

Yellow fever (YF) is maintained in South America primarily via a sylvatic cycle involving NHPs and forest mosquitoes (Haemagogus and Sabethes), with humans as incidental hosts. Urban transmission can occur via Aedes aegypti, leading to explosive epidemics. In Brazil, since 1942, cases have been attributed to the sylvatic cycle, but since 1998 the risk area has expanded, with large outbreaks in 2016–2018 affecting densely populated southeastern states. While YF seasonality is recognized and often attributed to climatic drivers that affect vector populations and transmission suitability, a critical counterpart is seasonality of exposure. Approximately 45% of human YF cases in Brazil occur in individuals engaged in agriculture or extractivism—activities with strong seasonal patterns. The study investigates whether seasonal agricultural activities (as proxies for exposure to sylvatic habitats) contribute to seasonal spillover. The research question is whether incorporating agricultural seasonality alongside climate/vegetation and demographic covariates improves prediction of monthly YF occurrence in humans and NHPs across Brazil, and which specific agricultural activities are most associated with increased YF reporting.

Previous work has characterized YF seasonality largely through climate-related mechanisms influencing vectors and transmission suitability, enabling seasonal forecasting. Broader literature links agriculture and land-use to infectious disease via landscape change effects on exposure, vector ecology, and reservoir host behavior, but there is limited research on how the seasonality of agriculture itself modulates spillover risk for vector-borne zoonoses. Historical studies document seasonal associations between rainfall, vector abundance (e.g., Haemagogus), and YF incidence. Research on spillover determinants emphasizes human–animal contact frequency and ecological interfaces, yet agricultural seasonality as a driver of contact and exposure has been underexplored. This study aims to fill that gap by explicitly testing agricultural calendar effects on YF reporting seasonality in Brazil.

Study design: A machine-learning classification framework (Random Forests) was used to predict monthly occurrences of YF reports at the municipality level in Brazil for humans, NHPs, and municipalities with both, aggregated over 2003–2018.

Outcomes: For each municipality and each calendar month (July–June surveillance year), outcomes were coded as binary presence/absence of a YF report, separately for: (i) human report, (ii) NHP report, and (iii) both human and NHP report.

Data sources and processing:

• Human YF reports: 2,423 initial cases (2003–2018) from the Brazilian Ministry of Health; after exclusions (missing date n=10, not geolocatable n=18), 2,395 cases remained, aggregated into 694 monthly occurrences across 434 unique municipalities.
• NHP epizootic reports: 3,209 confirmed events; after excluding 10 prior to 2003, 3,199 remained, aggregated into 771 monthly occurrences across 409 municipalities.
• Host demographics (H): Number of NHP species per municipality from NatureServe distributions; human population and proportion working in agriculture from IBGE; logarithm of human population included.
• Agricultural seasonality (A): Monthly binary indicators (0/1) for planting and harvesting states (first administrative level) from CONAB/Mapa agricultural calendars. Fifteen crops were catalogued; eight with substantial production were selected: peanuts, rice, common bean, castor beans, corn, soya, sorghum, and wheat. This yielded 16 binary covariates (planting and harvesting for each of 8 crops), mapped to municipalities by state.
• Agricultural output (O): Number of farms producing each of the eight crop types at municipality level from the 2017 Agricultural, Forestry and Aquaculture Census (IBGE).
• Climate/vegetation seasonality (C): Monthly multi-year averages (2003–2016) of rainfall, day and night temperature, diurnal–nocturnal temperature range, and Enhanced Vegetation Index (EVI), with 1- and 2-month lags for each. Variables were spatially aggregated to municipalities using population-weighted means based on LandScan 2015.

Covariate groupings and models: Covariates were standardized (zero mean, unit variance). Four covariate groups were defined: O (agricultural output), H (host demographics), A (agricultural seasonality), and C (climate/vegetation seasonality). All 15 non-empty combinations of these groups were fit as separate Random Forest models for each classification target, resulting in 15 model formulations.

Random Forest implementation: Models were implemented in R 3.5.1 using the ranger package. RFs modeled non-linearities and interactions. For interpretability, permutation variable importance and partial dependence plots were computed.

Validation and performance metrics: A spatial-block bootstrapping cross-validation was used to assess out-of-sample performance and mitigate spatial autocorrelation. A 5°×5° latitude/longitude grid assigned municipalities to blocks; with replacement, blocks were sampled to form training sets comprising ~60–70% of municipalities, with the remaining 30–40% used for validation. This procedure was repeated 200 times. Performance was evaluated using area under the ROC curve (AUC) for each classification outcome and Brier score for overall probabilistic accuracy. Within-sample fit was also compared via absolute total deviations between observed and predicted monthly totals and by Pearson correlations of monthly predictions versus data.

Collinearity considerations: RFs are robust to correlated predictors through random subspace selection. Empirically, seasonally varying agricultural covariates showed only moderate correlations with climate/vegetation variables, suggesting they capture partly independent processes (supported by a covariate correlation heatmap).

• Seasonality patterns: YF reports are highly seasonal. Nationally, 79.8% of all reports (human, NHP, or both) occurred January–March; the minimum monthly total was 1 report in October and the maximum 255 in January. NHP reporting occurred year-round at higher background levels than human reports (minimum 12 in August; maximum 181 in January).
• Model performance: Including agricultural seasonality markedly improved prediction of seasonal patterns, particularly for human reports. Out-of-sample AUC for human reports ranged from 0.80 (0.73–0.87) in the simplest model to 0.93 (0.90–0.96) in the full model with all covariate groups (OHAC). For NHP reports, AUCs ranged from 0.78 (0.75–0.82) to 0.92 (0.90–0.94) in models including both agriculture and climate (HAC/OHAC). For municipalities with both human and NHP reports, AUCs ranged from 0.73 (0.69–0.77) to 0.84 (0.81–0.87) in OHAC.
• Monthly fit metrics: The best model including both seasonality types (OHAC) had the lowest absolute total monthly deviations (Human 278.7; NHP 260.0; Both 71.2) and the highest Pearson monthly correlations (Human 0.99; NHP 0.95; Both 0.97). Models with only agricultural seasonality (OHA) also performed well (Human 0.99; NHP 0.95; Both 0.94), while climate-only seasonality (OHC) underperformed in reproducing the seasonal amplitude (Human 0.80; NHP 0.83; Both 0.80).
• Seasonal prediction behavior: Climate-only models overpredicted during low-season months (June–December) and underpredicted peaks (January–April). Models including agricultural seasonality better captured the timing and magnitude of peaks, though all models tended to underpredict the absolute epidemic peak.
• Spatial patterns: Predicted spatial distributions from the best model (OHAC) matched observed hotspots: southeastern Atlantic states, western and Amazonian areas for human reports; NHP reports were more widespread, including states without human cases (e.g., Bahia, Tocantins). Some overprediction occurred in parts of the Amazon where human cases were not reported in the study period.
• Variable importance: In OHAC, climate variables (especially night temperature and its 1–2 month lags, rainfall and its lags, temperature range, and EVI) ranked highest. Among agricultural output variables, numbers of bean, corn, and soya farms were influential. Among agricultural seasonality variables, rice harvesting and peanut planting had the greatest importance.
• Validation: Out-of-sample AUCs from spatial-block bootstrapping were similar to within-sample AUCs for human and NHP reports, with slightly larger discrepancies for the “both” classification, suggesting some overfitting for that category.

Findings demonstrate that seasonality of exposure, proxied by agricultural planting and harvesting cycles, is a key determinant of seasonal YF spillover risk, complementing climate-driven transmission suitability. Incorporating agricultural seasonality significantly enhanced prediction of human YF reports and improved reproduction of observed seasonal amplitude, indicating that human interactions with sylvatic habitats during agricultural activities drive seasonal spillover. Identification of specific crop-related activities (e.g., rice harvesting, peanut planting) associated with higher YF reporting provides actionable insights for targeted public health interventions. While climate remains essential for enabling sylvatic transmission, agricultural seasonality introduces an “anthropogenic seasonality” that modulates spillover timing and magnitude. Public health implications include refining surveillance timing and prioritizing vaccination and monitoring in high-exposure periods and locales. However, differences in surveillance sensitivity—especially for NHP epizootics—may influence observed patterns and model performance, and models predicting the “both” category may not fully capture distinct human versus NHP exposure dynamics.

The study shows that agricultural seasonality substantially improves prediction of seasonal YF spillover in Brazil beyond climate alone, underscoring the role of exposure dynamics linked to human agricultural activities. By pinpointing crop types and activities associated with elevated risk, the approach enables more precise targeting of vaccination and surveillance to periods and regions of highest spillover risk—crucial under resource constraints and amid changing YF epidemiology and land use. Future research should incorporate finer-resolution, quantitative measures of agricultural activity across more crop types, improved geo-localization of reports, and mechanistic studies to elucidate how specific agricultural practices elevate exposure to the sylvatic cycle.

• Surveillance heterogeneity: NHP reporting varies across regions, with some areas underreporting epizootics; intensified NHP surveillance often follows human case detection, introducing bias. Human reports are more consistent but still subject to under-ascertainment.
• Outcome definition: Presence/absence per month aggregates across years and does not reflect true transmission intensity; asymptomatic and mild cases are underdetected.
• Agricultural data granularity: Agricultural seasonality indicators are binary and at state-level, not quantitative or municipality-specific; limited to eight crop types.
• Spatial resolution and mapping: Some misalignment between state-level agricultural calendars and municipality-level outcomes; potential overprediction in certain Amazonian municipalities.
• Model concerns: Slight overfitting indicated for the “both human and NHP” category; RFs provide limited mechanistic insight despite strong predictive performance.
• Confounding and correlation: Although correlations between agriculture and climate covariates are moderate, residual confounding between climate-driven suitability and agricultural timing cannot be fully excluded.