Abstract
The emergence of H5N1 clade 2.3.4.4b in dairy cattle raised concerns over the safety of fluid milk. We developed two stochastic quantitative risk assessment models to represent the United States raw and pasteurized fluid milk supply chains and employed these models in baseline, sensitivity, and scenario analysis. Median (5th, 95th percentiles) probabilities of infection per 240-mL pasteurized and farmstore-purchased raw milk serving were 5.68E-15 (1.766E-16, 2.98E-13) and 1.13E-03 (5.16E-06, 3.82E-02), respectively. This metric was highly sensitive to the viral titers in infected cows’ milk. Pasteurization is highly effective at reducing this risk. Bulk tank milk PCR testing is more effective at reducing the probability of infection per raw milk serving than improving the identification and diversion of infected cows’ milk at harvest. These findings emphasize the importance of pasteurization and dairy cow disease surveillance in reducing the risk of H5N1 infection associated with fluid milk consumption.
Summary In this epidemiologic study, the public health risk of human H5N1 infection from drinking raw and pasteurized milk is investigated by use of quantitative risk assessment models.
Biosketch Katherine J. Koebel earned her Bachelor’s in Animal Science from Michigan State University and her Doctor of Veterinary Medicine from Cornell University. Her research draws from skills and experience in quantitative analysis, epidemiologic modeling, and dairy science. She is currently a PhD student in Epidemiology at Cornell University, researching One Health in the context of dairy cattle population medicine.
Introduction
Avian influenza (AI) threatens both avian and mammalian species. AI is an illness caused by multiple strains of the influenza A virus (IAV).(1,2) Of particular interest are the viruses of the H5N1 clade 2.3.4.4b. This clade emerged in 2020 and was subsequently detected in North America in late 2021.(3,4) Epidemiologic data highlight the predominance of this clade in reports of influenza A(H5) samples from both humans and animals.(5,6).
In March 2024, a spillover event resulted in the infection of dairy cattle in Texas with H5N1 clade 2.3.4.4b.(7,8) The initial infections have since spread to more than 800 herds nationwide as of 19 December 2024.(9) Infected cows present predominantly with clinical mastitis, usually in mid- to late-lactation, with additional signs such as anorexia, rumen hypomotility, and nasal discharge. Milk from these animals contains extremely high concentrations of virus (4.0-8.8 log10TCID50/mL (note: hereafter “logTCID50”)) and is infectious in a murine oral inoculation model.(7,10,11)
As of 19 December 2024, 37 human H5N1 cases with cattle exposure have been confirmed.(12,13) These cases were attributed to workplace exposure and no mortalities nor documented cases of foodborne transmission have been reported. Human cases are generally mild, including conjunctivitis and upper respiratory symptoms, and no human-to-human transmission has been proven.(14) Exposure in 2 cases is unknown (13). Additionally, on 24 November 2024, a voluntary recall was issued after a lot of retail raw milk from a California dairy herd tested positive for H5N1.(15) No resultant illnesses have been reported, but the California Department of Public Health acknowledges the potential risk of human infection from contact with or ingestion of the product.
Many questions remain about the risk of human infection with H5N1 from fluid dairy milk. Recent studies demonstrate high efficacy of heat inactivation of the virus in milk (10,16); for example, high-temperature short-time (HTST) pasteurization is posited to yield a minimum 12-log reduction in viral loads.(17) Despite the importance of pasteurization in reducing the risk of foodborne illness, consumption of raw milk remains popular in the U.S.; 4.4% of U.S. adults consume raw milk at least once per year.(18–20)
Timely research into the human infection risk associated with fluid milk consumption is essential to inform actions aimed at reducing said risk. The objectives of this study were to (i) estimate the public health risk of human H5N1 infection from consumption of raw or pasteurized fluid milk via quantitative risk assessment (QRA), (ii) assess intervention strategies, and (iii) identify key knowledge gaps.
Methods
Due to differences in the pasteurized and raw milk supply chains, two separate QRA models were developed in @RISK v.8.6.1 (Lumivero, Denver, CO), an add-on to Microsoft Excel 365 v.2406 (Microsoft Corp., Redmond, WA). Additionally, two raw milk purchasing pathways (farmstore and retail) were differentiated. Both models consisted of three main, and several additional, components (Figure 1); parameters are defined in Table 1.
We assume uniform spatial dispersion of the virus in fluid milk; ignorance of cumulative risk from consumption of multiple servings over time; a fixed number of herds shipping to a given processing plant; and unspecified fluid milk fat percentage (i.e. skim, 2%, 1%, or whole). The implications thereof are elaborated upon in Discussion.
Pasteurized Milk Model
This model represents a given fluid milk processing plant receiving, processing, packaging, and distributing milk from a fixed number of herds. Fluid milk is followed from the plant through retail until its consumption by a consumer. The plant receives milk from 100 herds (n), where the number of affected herds (k) shipping to that plant is calculated with Eqs. S1-3.
The herds are of average size HZ, with a healthy cow producing YH gal/d.(21,22) Infected cows comprise pherd proportion of the herd and give YI relative yield with milk titers VM.(7,23) A proportion D of these animals is identified and their milk is discarded at the time of harvest (i.e. does not enter the bulk tank). Diversion of this milk means identification of the animal as unfit to contribute to the milk supply, including but not limited to recognition of clinical illness or a laboratory-confirmed H5N1 infection.
Milk is pooled at the plant and held at 4°C for 24 hours, during which the viral titers experience log decay VD of 0.1 logTCID50/d.(16) Temperature-dependent log decay is calculated with Eq. S4. This method was applied to all storage/transportation stages in both models.
The milk is HTST pasteurized, which is estimated to yield >12-log reduction of H5N1.(17) Pasteurization log reduction (L) is therefore set at a baseline of 12 with additional evaluation during scenario analysis. Milk is packaged into gallon (3,785 mL) containers and the partitioning of virions modeled according to Nauta.(24) Gallons proceed through the cold chain: storage at the plant (UP-Z, TP-Z), transport to retail (UTR-Z, TTR-Z), storage at retail (USR-Z, TSR-Z), transport to the consumer residence (UTC-Z, TTC-Z), and storage at the residence (USC-Z, and TSC-Z).(25) Food waste of gallons (PDisc-Gal-Z) and servings (PDisc-Serv-Z) is accounted for.(26) The model culminates in the consumption of an 8-ounce (240 mL) serving. An exponential dose-response curve with parameter r is employed. Parameter r was calibrated using a mustelid model of ingestion of H5N1-contaminated meat (Eq. S5).(27) To predict the number of H5N1 cases nationally from drinking pasteurized fluid dairy milk, per annum U.S. consumption statistics were obtained from an FDA QRA.(28)
Validation of the pasteurized milk model (its exposure assessment) was performed using nationwide dairy food surveillance data which reported the logEID50/mL (50% egg infectious doses; considered synonymous with logTCID50 for our purposes) in pasteurized fluid milk at retail calculated via quantitative real-time RT-PCR (qrRT-PCR).(29) We reproduced this data by setting K=26 (see Eqs. 2-3), running the model with L=0, and recording the concentration of virus (logTCID50/mL) in fluid milk at retail. The output (logTCID50/mL) can be interpreted as the equivalent amount of viral genetic material remaining in the product post-pasteurization, approximating the approach in Spackman et al (compared to their results, we combined all fluid milk types as our model does not differentiate between product fat percentages). We compared our predicted logTCID50/mL with Spackman et al’s logEID50/mL graphically.
Raw Milk Model
This model represents a herd of size HW selling raw milk through either farmstore or retail pathways. No pooling of milk from multiple herds nor pasteurization occurs. Parameters pnat’l, pherd, YH, YI, VM, VD, D, and r are shared with the pasteurized milk model. In the farmstore pathway, gallon packages are sold directly to consumers after storage at the farmstore (UF, TF). In the retail pathway, on-farm storage lag (UP-W, TP-w), retail transport (UTR-W, TTR-W), and retail storage (USR, TTR) are modeled.(30) In both models packages are transported to (UTC-W, TTC-W, TH,) and stored at (USC-W, and TSC-W) the consumer residence. The milk temperature during transport to the residence (TTC-W) is calculated with Eq. S6. Accounting for raw milk perishability, the probability of spoilage and discard of the product is adapted from Crotta et al (Eq. S7).(31)
Analysis
The primary risk output reported is the “probability of infection per serving” (hereafter “p(infection)”) for both models. Sensitivity analysis was conducted with partial rank correlation coefficients (⍴) achieving statistical significance after Bonferroni correction at ɑ=0.05/(number of stochastic parameters in model). Parameters with |⍴| ≤ 0.1 were considered negligible (and of low practical significance) and were omitted from further discussion.(32) Coefficient strength ranges were adapted from this paper as well, where 0.1 < |⍴| < 0.40 is weak, 0.40 ≤ |⍴| < 0.70 is moderate, and 0.70 ≤ |⍴| is strong.
Scenario analyses investigated the effects of four parameters on p(infection) (Table 1): in the pasteurized milk model, (i) pasteurization (L); in both models, (ii) dose-response parameter r; and in the raw milk model, (iii) improved infected animal diversion (D) and (iv) bulk tank PCR testing. In the pasteurized milk model, L was tested at levels 6-, 8-, 10-, 12-(baseline), and 14-log reduction. Similarly, in both models, the dose-response parameter r was tested at levels 1E-6, 1E-8 (baseline), and 1E-10. Lastly, full factorial analysis was performed to examine two interventions aimed at reducing risk from consumption of raw milk. Parameters sensitivity (Se), specificity (Sp), and limit of detection (LoD) are employed regarding PCR testing;(33,34) all positive tanks (true or false) are discarded. All negatives (true or false) proceed to distribution.
Simulations are comprised of 50,000 iterations generated with Latin hypercube sampling and Mersenne twister pseudorandom number generation. Convergence testing with 5% tolerance and 95% confidence confirmed the sufficiency of the iteration number. Figures were produced in R v4.3.2 (R Foundation for Statistical Computing, Vienna, Austria).
Results
QRA models were developed for the raw and pasteurized milk supply chains (Figure 1). Projections for the concentration of viral material in retail gallon packages conform well to the results of previously conducted surveillance testing, (29) supporting the validity of the QRA models’ exposure assessment (Figure 2). Baseline p(infection) per serving indicated similarity between the two raw milk purchasing pathways (Figure 3 and Table S1). Paired with the increased access to farmstore-sold over retail raw milk in terms of state legality,(35) further description of results is limited to the farmstore pathway, except in sensitivity analysis.
Baseline scenario
Risk metrics in the pasteurized milk model were very low (Figure 3). 5th and 95th percentiles for p(infection) were 1.77E-16 and 2.98E-13, respectively, with a median value of 5.68E-15 (Table S1). The quantity of virus per contaminated serving was similarly low, ranging from 0 to 0.30 logTCID50 with a mean value of 0.019 logTCID50. A value of 0 logTCID50 corresponds to a single TCID50 in a serving and 0.30 logTCID50 corresponds to 2 TCID50 in a serving. Prevalence of serving contamination took percentiles of 2.43E-08 (5th) and 4.03E-05 (95th) with a median of 7.67E-07. Note that this refers to the presence of live infectious virus and not viral material alone. Out of 50,000 iterations, n=290 (0.58%) projected milk-borne human H5N1 infections with 1 (n=281), 2 (n=8), or 3 (n=1) infections projected in these iterations. An iteration represents a full calendar year of average U.S. fluid milk consumption under the simulated conditions.
In the raw milk model, risk metrics were considerably higher (12-log difference in medians) than in the pasteurized milk model (Figure 3 and Table S1). Median p(infection) was 1.13E-03, with percentiles of 5.16E-06 (5th) to 3.82E-02 (95th). Percentiles for the quantity of virus per contaminated serving ranged from 5.22 (5th) to 8.49 (95th) with a median value of 6.88 logTCID50. Median prevalence of serving contamination was 3.13E-02 (5th percentile: 3.21E-03; 95th percentile: 1.07E-01).
Sensitivity analysis
Tornado plots are given in Figures 4a-c. Parameter notation is available in Table 1. In both models, VM was very strongly correlated with the p(infection) (pasteurized: ⍴=0.94; raw, farmstore: ⍴=0.88, raw, retail: ⍴=0.83). This metric was also sensitive to pnat’l in both models (pasteurized: ⍴=0.80; raw, farmstore: ⍴=0.72; raw, retail: ⍴=0.63).
In the pasteurized milk model, parameter pherd demonstrated moderate positive correlation (⍴=0.52). YI was weakly correlated (⍴=0.15). HZ (⍴=-0.10) yielded negligible correlation. In the raw milk model, weak positive correlation was observed in pherd (farmstore: ⍴=0.29; retail: ⍴=0.24). For the farmstore purchase pathway, moderate negative correlations were observed in UF (⍴=-0.62), TF (⍴=-0.54), USC (⍴=-0.48), and TSC (⍴=-0.43). Similarly, in the retail purchase pathway, USR (⍴=-0.54) and TSR (⍴=-0.42) display moderate negative correlation. YI (farmtore: ⍴=0.09; retail: ⍴=0.08) was negligible in both purchasing pathways. Compared against their counterparts in the pasteurized milk model, the strengths and directions of statistically significant parameters were comparable with the exceptions of HW (farmstore: ⍴=-0.12; retail: ⍴=-0.10), which was considered weak in the farmstore pathway and negligible in retail, and parameters TH (farmstore: ⍴=-0.01) and UTC (farmstore: ⍴=-0.01), which failed to achieve significance in the retail pathway.
Scenario analysis
The median p(infection) was equal to 3.71E-10 and 3.49E-16 in the 6- and 14-log pasteurization scenarios, respectively (Figure 5 and Table S2). Increases in the pasteurization efficacy parameter L produced a progressively lower median p(infection) at each level tested.
A 100-fold increase in the dose-response parameter r (i.e. from a baseline of 1E-8 to 1E-6), meaning a decrease in the quantity of virus necessary to result in infection, proportionately increases (by 1-2 logs) p(infection), and vice versa (Figure 6 and Table S3). For example, for pasteurized milk, median p(infection) increases from 5.68E-15 to 5.68E-13 when r increases by 2-log. This effect was consistent in both raw milk purchasing pathways (Figure S1).
In the raw milk supply chain, bulk tank PCR testing was more effective than improved animal diversion in reducing p(infection); for all three levels of diversion, the addition of bulk tank PCR testing reduced median p(infection) by 1-2 log (Figures 7 and S2). In the farmstore purchase pathway (Table S4), without bulk tank PCR, increasing infected cow diversion over a 25% baseline to 50% reduced the median p(infection) from 1.13E-03 to 7.90E-04; with the addition of bulk tank PCR, this value at 25% diversion was 1.81E-05, then 1.26E-05 at 50% diversion. Results for the retail purchase pathway were similar (Table S5).
Discussion
Our model predicts the risk of human H5N1 infection from consumption of pasteurized fluid milk is extremely low (Figure 3 and Table S1), supporting claims about the safety of the domestic pasteurized fluid milk supply chain. According to research surrounding H5N1 and HTST pasteurization, 12-log reduction is highly efficacious at reducing the public health risk.(17) Even if a small quantity of infectious virions survive pasteurization, partitioning of these units during packaging further reduces infection risk. In iterations producing a non-zero amount of live virus in a serving, the dose-response parameter r is a significant driver of infection risk, as even 1-2 TCID50 can cause infection under the exponential dose-response model, albeit with a very low probability. With the importance of r demonstrated in scenario analysis (Figure 6 and Table S3), further research into H5N1 ingestion dose-response is needed. As seen in sensitivity analysis (Figure 4), the predicted low infection risk is also due in small part to the dilution effect of pooling from multiple herds and the reduced milk yield of infected cows that are not diverted from the supply chain. Also in sensitivity analysis, the decreased H5N1 infection risk in raw milk associated with increased storage temperatures/times is a function of viral decay; excessive or inappropriate storage of dairy products is not advisable due to the risk of illness from spoilage pathogens.
The human public health burden from pasteurized milk is predicted to be low, in that <1% of iterations predicted a maximum of 3 cases attributable to its consumption. Note that this output is calculated with per annum consumption data and thus represents the projected annual number of infections extrapolating from epidemiologic data collected between March and November of this EID situation; the implications behind this must be considered. Expansion of the virus into California has caused rapid increases in herd infections; as H5N1 continues to spread, by a means still not well understood, it is unknown where in the epidemic curve the U.S. national herd is currently. The model addresses this with stochasticity in the monthly national herd-level prevalence parameter pnat’l, but projections of risk may be over- or underestimated.
While a >12-log reduction is posited (17), demonstration of pasteurization efficacy is fundamentally limited by the viral concentration of the initial sample.(16) In order to definitively demonstrate efficacy, assays with supraphysiological concentrations of the virus must be conducted. Pasteurization remains the most effective method for reducing risk amidst this H5N1 outbreak, and per scenario analysis (Table S2), should its efficacy be overestimated or a pasteurization failure go undetected, the risk to public health would be substantial.
Unknowns surrounding the clinical manifestations of human H5N1 infection by ingestion must be considered. The two cases of Michigan dairy workers cite potential exposure from splashing milk and direct contact with oronasal secretions during animal care procedures.(36) It is unknown what symptoms an infection from ingestion would produce in a human. In laboratory animal studies, clinical signs after H5N1 ingestion ranged from none to weight loss and lethargy, up to and including mortality.(27,10,11,37) It is thought that viral contact with the oropharynx during deglutition is the route of entry sufficient to establish infection in these models. Perhaps presentation with conjunctivitis in humans is a function of the route of exposure, including the splash of contaminated material into the mucous membranes of the face. Other routes of exposure may produce a different array of symptoms. Therefore, the absence of a documented case attributable to milk consumption may be a result of asymptomaticity or non-reporting. Exposure in two H5N1 cases remains unattributed.(13)
Wastewater surveillance has proven applications in the context of this outbreak. Detection of H5N1 in Texas wastewater coincided with the emergence of the disease syndrome in dairy cattle and predates the official announcement of the causative agent.(38) While genetic analysis of wastewater samples indicates primarily bovine contribution (likely milk effluent from farms or processing plants), viral shedding in human sewage cannot be ruled out.(39) As such, wastewater surveillance in areas without cattle could be used to identify human cases even if said cases are asymptomatic.
Ingestion of raw milk carries inherent risk of foodborne illness, primarily bacterial.(18–20) Our model demonstrates that, without pasteurization, the p(infection) is dramatically higher. Given the popularity of raw milk amongst U.S. consumers, in light of this EID there is an urgent need for new technological, educational, and policy solutions to protect public health. Here we have demonstrated the efficacy of two interventions. PCR testing of raw milk herd bulk tanks is highly efficacious in reducing p(infection), more so than improving the ability to identify and divert milk from infected cows. However, despite its comparatively weaker effect, diverting infected cows is still important, as this relates to sensitivity in both models to viral titers in infected cow milk (VM). The U.S.D.A. is currently implementing a bulk tank milk testing surveillance plan at the state and locoregional levels.(40)
The assumptions and limitations of any QRA must be acknowledged for appropriate interpretation of its results. We assume uniform spatial dispersion of virus in milk; if this is incorrect and clustering is present, the resultant distribution for the quantity of virus per serving, and therefore the predicted p(infection), is incorrect as a smaller proportion of packages will contain more virions paired with a lower prevalence of serving contamination. This may also be true for different milkfat levels, if virions distribute differently in fat globules versus the liquid fraction. As Spackman et al report slight differences in the concentration of viral material amongst different milkfat levels,(29) continued surveillance is required to determine if this is a function of processing or sample size. Our model does not differentiate between milkfat levels; if sufficient data become available to allow for inclusion of milkfat as a parameter, projections of risk may shift in either direction. Next, a fixed number of herds contribute to the modeled plant. Due to the dilution effect demonstrated previously, deviation from a fixed n in real life may increase or decrease risk. Lastly, note that we report risk p(infection) as the probability of infection per single, 8-ounce (240 mL) serving. We assume that p(infection) of a given serving is independent of the “next” consumed serving. In reality, there may be a cumulative risk from consumption of multiple subinfectious doses over time, therefore underestimating true risk.
Additional information is needed to make a QRA accounting for all dairy products. As research characterizing H5N1 in other dairy foods becomes available, the framework of our model may be adapted. As of 6 December 2024, a new Federal Order has been mandated requiring processor-level raw milk H5N1 surveillance. Implementation of this initiative is supported by the projections of our models, and results from this heightened surveillance will prove useful in future model developments.(41) For example, academics and industry stakeholders have expressed interest in studying H5N1 in raw milk cheese. The expansion of our models to fill this knowledge gap is being pursued by the authors.
In conclusion, based on the results of this QRA and the most up-to-date scientific literature, the risk of human H5N1 infection from consumption of pasteurized fluid milk is extremely low. Consumption of raw milk has historically been discouraged due to the risk of foodborne illness; in the context of this EID, the emergence of a novel potential milkborne pathogen emphasizes both the risk of this practice and the need for solutions to effectively protect public health. This risk may be reduced with implementation of risk mitigation interventions, including thermal treatment and bulk tank PCR testing. As new information becomes available, the models will be refined and new applications considered.
Data Availability
All data produced in this work are available upon reasonable request to the corresponding author.
Acknowledgements
Research reported in this publication was supported by the Office of the Director, National Institutes of Health of the National Institutions of Health (NIH) under Award Number T32ODO011000. The content is solely the responsibility of the authors and does not necessarily represent the official views of the National Center for Research Resources or the NIH. Additionally, this research was partially supported by grants to RI from the National Institute of Food and Agriculture, USDA, Hatch under Accession Number 7000433, as well as Multistate Research Funds Accession Number 1016738.