Improvement of quantitative microbiological risk assessment (QMRA) methodology through integration with gaenetic data

Abstract

Quantitative microbiological risk assessment (QMRA) methodology aims to estimate and describe the transmission of pathogenic microorganisms from animals and food to humans. In microbiological literature, the availability of whole genome sequencing (WGS) data is rapidly increasing, and incorporating this data into QMRA has the potential to enhance the reliability of risk estimates. This study provides insight into which are the key pathogen properties for incorporating WGS data to enhance risk estimation, through examination of example risk assessments for important foodborne pathogens: Listeria monocytogenes (Lm), Salmonella, Campylobacter and Shiga toxin-producing Escherichia coli. By investigating the relationship between phenotypic pathogen properties and genetic traits, a better understanding was gained regarding their impact on risk assessment. Virulence of Lm was identified as a promising property for associating different symptoms observed in humans with specific genotypes. Data from a genome-wide association study were used to correlate lineages, serotypes, sequence types, clonal complexes and the presence or absence of virulence genes of each strain with patient's symptoms. We also investigated the effect of incorporating WGS data into a QMRA model including relevant genomic traits of Lm, focusing on the dose–response phase of the risk assessment model, as described with the case/exposure ratio. The results highlighted that WGS studies which include phenotypic information must be encouraged, so as to enhance the accuracy of QMRA models. This study also underscores the importance of executing more risk assessments that consider the ongoing advancements in OMICS technologies, thus allowing for a closer investigation of different bacterial subtypes relevant to human health.

Summary

The quantitative microbiological risk assessment (QMRA) methodology aims to estimate and describe the transmission of pathogenic microorganisms from animals and food to humans. QMRA can enhance the understanding of the impact of interventions in reducing the public health burden, as well as identify the most significant transmission routes for specific pathogens. However, there are important limitations associated with QMRA. Estimated numbers of human cases often greatly exceed what is found by epidemiological studies suggest. Further, subtypes of pathogens are not distinguished, even though certain subtypes may dominate the food chain but are rarely found in humans, or the other way around. In microbiological literature, the availability of whole genome sequencing (WGS) data is rapidly increasing, and incorporating this data into QMRA has the potential to enhance the reliability of risk estimates. This study aimed to improve the QMRA methodology by using WGS data for different subtypes of a food-borne pathogen. The initial phase of the study focused on investigating the variability of pathogen properties among strains and their impact on risk. Four farm-to-fork QMRAs were studied, for Listeria monocytogenes (Lm), Salmonella, Campylobacter and Shiga toxin-producing Escherichia coli. The goal was to identify the most critical pathogen properties (e.g. growth, virulence) and related parameters that could be used to improve QMRA through the inclusion of genotype information. Therefore, a comprehensive list of properties and parameters used throughout the food chain was made, and the sensitivity of each parameter in every QMRA was evaluated to assess the magnitude of their impact on model outcomes. Furthermore, the magnitude of variability of each parameter was examined through a review of relevant literature studies, while aiming also to establish connections between the properties/parameters of interest and the genome variability observed in different strains. Studies investigating the relationships between subtypes, genes and the variation in pathogen properties by examining the phenotypes of different subtypes were considered crucial for improving QMRA. Therefore, the set of identified pathogen properties considered valuable for QMRA improvement was compared with the available literature data, so as to gain a comprehensive understanding of the WGS data gaps from the perspective of QMRA. This knowledge is essential for identifying the key aspects that warrant attention in future research. The virulence of Lm, which impacts the dose–response phase of the QMRA model, was identified as a promising candidate for associating different symptoms observed in humans with specific Lm genotypes. This was investigated further using a genome-wide association study to correlate lineages, serotypes, sequence types, clonal complexes and the presence or absence of virulence genes of each strain with disease outcomes.

Finally, the study aimed at investigating the size and nature of the effect of incorporating WGS data into the QMRA using a QMRA model that integrated the genomic traits into the virulence property of Lm. The results highlighted that additional studies focusing on these aspects are needed to enhance the knowledge and refinement of the QMRA methodology.

This study highlighted the importance of developing improved risk assessment methodology in response to the continuous advancements in OMICS technologies, which enables a more detailed examination of various bacterial subtypes relevant to human health.

Introduction

The quantitative microbiological risk assessment (QMRA) is recognised as a robust mathematical modelling approach for estimating human exposure risks to pathogenic microorganisms through food, environment or direct contact with infected animals. QMRA provides valuable insights into the effectiveness of interventions aimed at reducing public health burden and helps to identify the most important transmission routes for specific pathogens (Evers and Bouwknegt, 2016). In the context of food-borne transmission, farm-to-fork models estimate the number of illness cases resulting from pathogen ingestion following the consumption of contaminated food. These models provide a comprehensive description of the entire food chain, starting from the farm phase and ending with food preparation by consumer (Chardon and Evers, 2017). A critical component of QMRA models is the estimation of the number of human cases for the specific disease under consideration. Therefore, there is growing interest in improving the precision of these estimates since QMRA values often exceed epidemiological estimates. In QMRA studies, there is an increasing awareness to differentiate pathogen subtypes based on phenotypic characteristics that are related to genotypic traits. This distinction allows for the identification of strains that may be dominant in the food chain but rarely found in humans, or strains that prevail in humans compared to their presence in food sources. The availability of whole genome sequencing (WGS) data in the microbiological literature is expanding rapidly, and its incorporation into QMRA enables the development of more targeted risk assessments and management strategies (Nielsen et al., 2017; Fritsch et al., 2018, 2019). Therefore, the aim of this work was to improve the QMRA methodology by using WGS data to discriminate between different subtypes of a food-borne pathogen. This subtype distinction, based on the analysis of genomes, aims to provide a deeper understanding of the potential of integrating genotypic information to improve the reliability of risk estimates.

Data and methodologies Study setting

The study of four farm-to-fork QMRAs on four important food-borne pathogens (L. monocytogenes, Salmonella, Campylobacter and Shiga toxin-producing E. coli – STEC) was performed to determine which pathogen properties and related parameters are the most important candidates for QMRA improvement through the integration of genomic data. A list of properties and parameters used throughout the food chain was made, evaluating the magnitude of the effect of changing each parameter on the model outcome (sensitivity) in every QMRA studied. Moreover, the relation between the presence of particular genomic traits and phenotypic differences (variability) and the magnitude of this variability was evaluated studying literature. For each property with both a high sensitivity and variability, the WGS data gaps were assessed from the viewpoint of QMRA, to understand which aspects are important to focus on in future research.

An in-dept study was performed on the virulence property of Lm in the dose–response (or case/exposure) phase of the QMRA model using a RIVM (Rijksinstituut voor Volksgezondheid en Milieu, the Netherlands) Lm database. The database included Lm strains collected in humans and in food in the Netherlands. For each strain were reported the lineage, serotype, sequence type, clonal complex (CC), presence or absence of virulence genes (as reported by the ‘Virulence factor database’, http://www.mgc.ac.cn/VFs) and the allelic number for each virulence gene, obtained analysing the complete genome (Coipan et al., 2023). Concerning human strains, the patient's characteristics were included in the database (such as age, gender, decease, pregnancy status, underlined disorders, systematic use of antacids or immunosuppressive). Moreover, when available, the symptoms following the Lm infection were reported, as sepsis, meningitis, gastrointestinal disease, lung inflammation, encephalitis, endocarditis. Concerning food strains, for each isolate was reported the food source, which included vegetable, pork, fish, beef and veal, chicken, turkey, mutton/lamb and game meat.

Methodologies Random forest

The random forest method was applied to human isolates from the Lm database. The goal was to predict specific symptoms during Lm infection based on isolate characteristics, including CC, patient status (e.g. deceased, pregnant), underlying disorders, immunosuppressive or antacid use, gender, age and presence or absence of certain virulence genes. The random forest method is a machine learning algorithm able to combine the outputs of multiple decision trees, considering the isolate characteristics, to predict the occurrence of specific symptoms resulting from the ingestion of Lm-contaminated food. Additionally, a heatmap was generated and correlation coefficients were calculated, considering virulence genes and symptoms as variables, to better understand the most significant characteristics to focus on. The analysis was performed using Rstudio software (version 4.2.1, http://www.rstudio.org).

Swift quantitative microbiological risk assessment

To assess the impact of incorporating WGS data into QMRA, a QMRA model was expanded to include relevant genomic features obtained from the analysis of the Lm database. Further, various scenarios were examined using the sQMRA (swift QMRA) tool, to study the implications of dividing exposures based on the prevalence of each CC in chicken products (Chardon and Evers, 2017). The sQMRA tool is a risk assessment model that provides a comprehensive description of the transmission and propagation of a pathogen in or on a food portion throughout part of the food chain, and its impact on consumers. Uncertainty of prevalence was included in our calculations using @Risk (version 8.3.2, https://www.palisade.com/risk).

Results

A comprehensive study of farm-to-fork QMRAs was performed for four important food-borne pathogens: Lm, Salmonella, Campylobacter and STEC. The risk assessment reports for each pathogen and each phase of the assessment were thoroughly examined to identify the key pathogen properties and related parameters that could be targeted for QMRA improvement. A list of properties and parameters used throughout the food chain was made, and the sensitivity of each parameter in the QMRAs was evaluated to determine the extent of the impact of parameter changes on the model outcomes. Moreover, the variability of each parameter among different strains, reflecting differences in genotype, was assessed based on available literature studies. The pathogen properties with significant variability and sensitivity were identified, highlighting the need to prioritise them in future QMRA improvements.

Concerning Campylobacter, a farm-to-fork QMRA model focusing on Campylobacter in broiler meat was studied (RIVM, 2005). Within the processing phase of the model, the invasiveness and colonisation potential of Campylobacter was found to have an important impact on the risk outcomes, and this can also be influenced by genotype variations among different strains (Chaloner et al., 2014). Moreover, during the storage phase at home, the survival or inactivation of Campylobacter due to refrigerated storage conditions was identified as an important property affecting the overall risk. Several studies reported variation in viability profiles of different strains at various temperatures, potentially linked to differences in their genotypes, but further studies are needed to clarify this aspect (Murphy et al., 2006). In the ingestion phase, the ability of Campylobacter to cause illness was found an important property; studies indicated that variations in the expression of virulence genes can lead to differences in adherence, invasion, intracellular survival capacities and toxin production (Janssen et al., 2008). Other noteworthy properties included environmental persistence, biofilm formation capacities and Campylobacter inactivation through heating.

The study of Salmonella was performed using a farm-to-consumption model on Salmonella in slaughter and breeder pigs (Hill et al., 2010). The assessment of Salmonella transmission during the process determined some promising properties and parameters involving differences in virulence, pathogenesis and intestinal persistence based on the genomic traits of the strains (Cui et al., 2021). Additionally, literature reported a significant diversity in stress resistance among various Salmonella strains due to the expression of different stress-induced genes (Wang et al., 2020). However, further studies are needed to understand whether the model risk outcome is sensitive to these properties. Finally, the ability to survive cleaning processes due to biofilm formation capacity can be influenced by differences in genotype characteristics that lead to better persist in the environment, but it was shown to have a minimal impact on the model risk outcome (Dantas et al., 2020). This study did not produce important pathogen properties for our aim and based on our criteria.

The study of STEC involved a farm-to-fork QMRA focused on steak tartare (Nauta, 2001). The STEC intestinal colonisation potential can affect the excretion rate, resulting in a possible contamination of the carcass in the slaughterhouse phase of the model, but it can also affect the dose–response phase because it can be important for the STEC virulence. This property had a high sensitivity in relation to the model outcome, and some studies reported that the adherence and gut colonisation capacities may be associated with differences in the expression of virulence and adherence genes (Pielaat et al., 2015; Barth et al., 2020). The STEC ability to cluster and their susceptibility to inactivation during storage were identified as important properties due to their high sensitivity in relation to the model outcome, but further studies are needed to determine if these properties are associated with specific genetic characteristics.

The study of Lm involved the analysis of two risk assessments which focused on estimating the risks of listeriosis associated with the consumption of ready-to-eat (RTE) foods (FDA and FSIS, 2003; Interagency Workgroup, 2013). This provided valuable insights into the properties and parameters of Lm important for enhancing the accuracy of QMRA models and risk estimation. In particular, the ability of Lm to grow at storage temperature showed high variability attributed to genotypic traits, even if its sensitivity was not quantified in terms of its impact on the outcome variation (Cordero et al., 2016). The thermal inactivation during food preparation showed differences based on the lineage of diverse Lm strains, but also in this case, its sensitivity remains unknown (Liu et al., 2021). The ability of Lm to persist in the environment and form biofilms was important during food preparation, with high sensitivity related to the risk outcome, and numerous studies suggested that persistent Lm strains have distinct genotypes characterised by the presence or absence of specific genes (Ferreira et al., 2014; Mazza et al., 2015). In the dose–response phase, variation of the ability of Lm to cause infection was identified as an important property, although the impact of varying this property on the risk outcome remains unknown. Recent studies highlighted differences in virulence among Lm strains based on the presence or absence of a truncated inlA gene, suggesting that the severity of illness should be considered when improving QMRA models using WGS data (Pouillot et al., 2015; Chen et al., 2020). Moreover, several studies reported differences in virulence between strains based on variations in lineages, serotypes and CCs (Nielsen et al., 2017; Fritsch et al., 2019; Maury et al., 2019; Quereda et al., 2021; Cardenas-Alvarez et al., 2022; Muchaamba et al., 2022).

In conclusion, considering the four pathogens studied, further research is required to establish possible connections between specific phenotypic and genotypic features and to gain a better understanding of how variations in these properties impact the overall risk outcome.

Given that the Lm virulence property affecting the risk characterisation phase of the model was found as a promising property, it was further investigated using the RIVM Lm database to find a relationship between the symptoms experienced by patients and specific genotype features. To identify the most promising symptom to be predicted, the data were visualised using a heat map, and correlation coefficients were calculated. The random forest method was then performed to predict the symptoms based on the presence or absence of specific virulence genes, but the results were not favourable. The error rate was determined to be 65%, indicating a high level of misclassification. The model accuracy (which equals 100 minus the error rate), representing the proportion of correctly predicted outcomes, was thus only 35%, which suggested that the predictive power of the gene presence/absence for symptoms in the database was limited. Further, random forest analysis was performed to predict symptoms based on isolate characteristics such as age, gender, decease, pregnancy status, underlined disorders, systematic use of antacids or immunosuppressive and an analysis using only CCs as predictors was also conducted. However, in both cases, the error rate remained high, and the model accuracy was low. These findings indicated that the presence or absence of specific virulence genes, isolate characteristics and CCs have no predictive power for symptoms in this database. Further studies and possible evaluation of additional factors are necessary to improve the prediction of symptoms to be used possibly in QMRA models.

The study used prevalence data from the Lm database, including both human and food isolates, to investigate the impact of splitting up exposures over CCs, using a QMRA model. To estimate the total exposure for each CC, the fraction of positive food isolates per CC and food in the RIVM database was used. A stochastic approach was used to calculate the uncertainty of these fractions. The number of human cases per CC per year in the Netherlands was determined using the Lm incidence in humans in 2021 (Benincà et al., 2022). By calculating the case/exposure ratio, it was possible to identify the CCs that exhibited higher virulence compared to others, despite the high uncertainties calculated. It proved that our results were consistent with existing literature reporting certain CCs as primarily found in food environments and considered hypo-virulent (Quereda et al., 2021). On the other hand, some CCs were predominantly associated with human cases and characterised as hyper-virulent due to their enhanced ability to colonise the intestine and invade the intestinal mucosa (Quereda et al., 2021).

Calculations with a baseline sQMRA model provided a reference point for comparing the model with and without genetic information, useful for a better understanding of the specific contributions of genomic characteristics to the risk estimates. Investigation of different scenarios for some CC's of interest were performed, to provide insights into the specific risk outcomes associated with different genomic characteristics of Lm strains.

Conclusion

This study highlighted that by incorporating WGS data into the QMRA model, a more comprehensive and accurate assessment and understanding of the risks associated with pathogen contamination throughout the food chain may be achieved. This integration will allow for a better understanding of the potential differences in risk estimates due to the specific genomic traits of pathogen strains. Such a model can assist in developing more targeted risk management strategies and interventions to mitigate the risks associated with food-borne pathogen contamination. In this study the most important pathogen properties to focus on were determined, and first attempts were made to combine WGS with QMRA.

Recommendations

The findings of this study emphasised the importance of including WGS data in QMRA. However, further research and in-dept studies that combine genomic data with phenotype analysis are indispensable to make progress in this. Such studies will allow for incorporation of pathogen subtypes based on genomic traits into QMRA models, so that more accurate risk assessments and risk management strategies can be achieved. Otherwise this will not be possible.