microorganisms Article Prevalence and Genomic Diversity of Salmonella enterica Recovered from River Water in a Major Agricultural Region in Northwestern Mexico Irvin González-López 1 , José Andrés Medrano-Félix 2 , Nohelia Castro-del Campo 1 , Osvaldo López-Cuevas 1 , Jean Pierre González-Gómez 1 , José Benigno Valdez-Torres 1 , José Roberto Aguirre-Sánchez 1 , Jaime Martínez-Urtaza 3 , Bruno Gómez-Gil 4 , Bertram G. Lee 5 , Beatriz Quiñones 5 and Cristóbal Chaidez 1, * 1 Centro de Investigación en Alimentación y Desarrollo A.C. (CIAD), Coordinación Regional Culiacán, Laboratorio Nacional para la Investigación en Inocuidad Alimentaria, Culiacán 80110, Sinaloa, Mexico; irvin.gonzalez@estudiantes.ciad.mx (I.G.-L.); ncastro@ciad.mx (N.C.-d.C.); osvaldo.lopez@ciad.mx (O.L.-C.); jgonzalez.219@estudiantes.ciad.mx (J.P.G.-G.); jvaldez@ciad.mx (J.B.V.-T.); jose.aguirre.dc18@estudiantes.ciad.mx (J.R.A.-S.) 2 Investigadoras e Investigadores por México, Centro de Investigación en Alimentación y Desarrollo A.C. (CIAD), Coordinación Regional Culiacán, Laboratorio Nacional para la Investigación en Inocuidad Alimentaria, Culiacán 80110, Sinaloa, Mexico; jose.medrano@ciad.mx 3 Department of Genetics and Microbiology, Universitat Autờnoma de Barcelona, 08193 Bellaterra, Spain; jaime.martinez.urtaza@uab.cat 4 Centro de Investigación en Alimentación y Desarrollo A.C. (CIAD), Coordinación Regional Mazatlán, Citation: González-López, I.; Medrano-Félix, J.A.; Castro-del Acuicultura y Manejo Ambiental, Mazatlán 82100, Sinaloa, Mexico; bruno@ciad.mx 5 U.S. Department of Agriculture, Agricultural Research Service, Western Regional Research Center, Campo, N.; López-Cuevas, O.; Produce Safety and Microbiology Research Unit, Albany, CA 94710, USA; bertram.lee@usda.gov (B.G.L.); González-Gómez, J.P.; Valdez-Torres, beatriz.quinones@usda.gov (B.Q.) J.B.; Aguirre-Sánchez, J.R.; Martínez- * Correspondence: chaqui@ciad.mx; Tel.: +52-(667)-480-6950 Urtaza, J.; Gómez-Gil, B.; Lee, B.G.; et al. Prevalence and Genomic Abstract: Salmonella enterica is a leading cause of human gastrointestinal disease worldwide. Given Diversity of Salmonella enterica that Salmonella is persistent in aquatic environments, this study examined the prevalence, levels and Recovered from River Water in a Major Agricultural Region in genotypic diversity of Salmonella isolates recovered from major rivers in an important agricultural Northwestern Mexico. region in northwestern Mexico. During a 13-month period, a total of 143 river water samples Microorganisms 2022, 10, 1214. were collected and subjected to size-exclusion ultrafiltration, followed by enrichment, and selective https://doi.org/10.3390/ media for Salmonella isolation and quantitation. The recovered Salmonella isolates were examined microorganisms10061214 by next-generation sequencing for genome characterization. Salmonella prevalence in river water was lower in the winter months (0.65 MPN/100 mL) and significantly higher in the summer months Academic Editors: Javier Garaizar and Lorena Laorden (13.98 MPN/100 mL), and a Poisson regression model indicated a negative effect of pH and salinity and a positive effect of river water temperature (p = 0.00) on Salmonella levels. Molecular subtyping Received: 24 May 2022 revealed Oranienburg, Anatum and Saintpaul were the most predominant Salmonella serovars. Single Accepted: 12 June 2022 nucleotide polymorphism (SNP)-based phylogeny revealed that the detected 27 distinct serovars Published: 14 June 2022 from river water clustered in two major clades. Multiple nonsynonymous SNPs were detected in Publisher’s Note: MDPI stays neutral stiA, sivH, and ratA, genes required for Salmonella fitness and survival, and these findings identified with regard to jurisdictional claims in relevant markers to potentially develop improved methods for characterizing this pathogen. published maps and institutional affil- iations. Keywords: environmental microbiology; food safety; foodborne pathogen; genomics; river water; Salmonella; serovars; single nucleotide polymorphisms; ultrafiltration method Copyright: © 2022 by the authors. Licensee MDPI, Basel, Switzerland. 1. Introduction This article is an open access article distributed under the terms and Salmonella enterica is recognized as a significant causative agent of human gastrointesti- conditions of the Creative Commons nal infections worldwide from foodborne and waterborne sources, and non-typhoid human Attribution (CC BY) license (https:// salmonellosis is responsible for millions of human illnesses in a typical year [1–3]. Among creativecommons.org/licenses/by/ the bacterial foodborne pathogens, S. enterica is one of most commonly reported in Mex- 4.0/). ico, according to the National Epidemiological Surveillance System (Sistema Nacional de Microorganisms 2022, 10, 1214. https://doi.org/10.3390/microorganisms10061214 https://www.mdpi.com/journal/microorganisms Microorganisms 2022, 10, 1214 2 of 17 Vigilancia Epidemiológica de México (SINAVE)) [4]. In particular, human salmonellosis is associated with the consumption of a diverse variety of food commodities, such as various fruits, vine vegetables, leafy greens, spices, tree nuts, poultry, beef, sprouts, flour, and ready- to-eat foods [4]. The health benefits associated with eating a well-balanced diet, including fresh fruits and vegetables, have contributed to a significant increase in the consumption of fresh produce in Mexico and many other countries during the past decades [5–7]. How- ever, this higher consumption of fresh and minimally processed ready-to-eat produce has been associated with an increased incidence in foodborne related outbreaks and sporadic illnesses [8]. The identification of relevant environmental reservoirs as potential sources of Salmonella transmission to fresh produce is thus needed as a preventive measure to ensure a safer food supply. Although Salmonella commonly resides in the gastrointestinal tract of animals, this pathogen is highly persistent in aquatic environments due to its ability to prevail for long periods in soils and sediments and to successfully adapt to stressful conditions, such as fluc- tuating temperature and pH, desiccation, and nutrient deprivation [9–11]. In the subtropical environments found in the Culiacan Valley, surface waters such as rivers and streams have been reported as relevant reservoirs for Salmonella due to optimal growth temperatures and the availability of nutrients (glucose, N-acetyl-D-glucosamine and D-glucosamine acid), enabling the prevalence of Salmonella in these aquatic environments [12,13]. Additionally, several studies indicated that organic matter transported into rivers by stormwater runoff promotes an increase in Salmonella populations [14–16]. Approximately 2500 S. enterica serovars have been currently identified worldwide based on cell surface H- and O-antigens [17]. In Mexico, the highest incidence in human infections has been associated with S. Typhimurium, S. Enteritidis, S. Derby, S. Agona and S. Anatum, which are also the most predominant serovars recovered from diverse food sources [4,18]. However, these serovars, causing the majority of human illnesses, rarely coincide with those with the highest prevalence in aquatic environments, including S. Oranienburg, S. Saintpaul, S. Minnesota, S. Give, and S. Infantis [19–21] recovered from river water belonging to the Culiacan Valley in northeastern Mexico. The low environmental prevalence of clinically relevant serovars provides a challenge for accurately determining the origin of the source of contamination and the identification of conducive conditions promoting the prevalence of these serovars in aquatic environments [9,18,22]. To aid in identifying sources of contamination, the use of high-resolution sequencing methods, such as whole genome sequencing, in conjunction with bioinformatic approaches can assist with the identification of genetic markers distinguishing environmental and clinical isolates, allowing a better understanding of the genotypic profiles of Salmonella strains recovered from aquatic environments in Mexico [23]. To implement strategies that allow an improved detection sensitivity for Salmonella in agricultural environments, the present study employed an ultrafiltration system [24,25] as an effective recovery method for Salmonella from samples collected from the major rivers in the Culiacan Valley, an important region in Mexico for the production of fresh produce commodities for domestic and international consumption [5–7]. The effects of several water physicochemical parameters on Salmonella recovery levels were also further investigated, and genome sequencing and bioinformatics tools were subsequently employed to deter- mine the phylogenetic relationships of the S. enterica serovars recovered from irrigation rivers in this important agricultural region in northwestern Mexico. 2. Materials and Methods 2.1. Collection of River Water at Locations in the Culiacan Valley, Northwestern Mexico The Culiacan Valley in the state of Sinaloa is a subtropical region located in northwest- ern Mexico, and the main activity of this region is agriculture, which is supplied with water mainly by three rivers. Originating at the Sierra Madre Occidental, the Humaya River and Tamazula River are supplied with water by the Adolfo López Mateos Dam and the Sanalona Dam, respectively, and converge in the city of Culiacan to form the Culiacan River, Microorganisms 2022, 10, x FOR PEER REVIEW 3 of 17 Microorganisms 2022, 10, 1214 water mainly by three rivers. Originating at the Sierra Madre Occidental, the Humaya 3 of 17 River and Tamazula River are supplied with water by the Adolfo López Mateos Dam and the Sanalona Dam, respectively, and converge in the city of Culiacan to form the Culiacan River, which which thenfrom then flows flowseast from to east westtotowest empty to empty into theintoGulf theofGulf of California California (Figure(Figure 1). These 1). Theserivers three threecomprise rivers comprise the Culiacanthe Culiacan River basin, River which basin, which atiscoordinates is located 105◦ 410 and located at coordinates 108 ◦ 40 W 105°41′ and and 24◦ 24W 108°4′ 0 and and 26 ◦ 240 N 24°24′ and 26°24′ [26]. Eleven N [26]. riverEleven sites inriver thesites Humayain theRiver, Humaya River, Tamazula Tamazula River, River, and and Culiacan Culiacan River River the throughout throughout Culiacanthe Culiacan Valley Valley were were sampled duringsampled dur- a 13-month period, betweenperiod, ing a 13-month June 2018 and November between June 2018 and 2019.November The sampled 2019.sites Thewere named sampled sites sites wereA, B, and C for the Humaya River; sites D, E, and F for the Tamazula named sites A, B, and C for the Humaya River; sites D, E, and F for the Tamazula River; River; and sites G, H, I, andJ, and sitesKG, forH,the I, J,Culiacan and K for River the (Table Culiacan 1). River The criteria (Tableemployed for the 1). The criteria selectionfor employed of the sampling selection of points was based the sampling on the points wasfollowing: based on the (i) following: sites in agricultural (i) sites in areas where areas agricultural rural populations were present;were where rural populations (ii) sites where present; (ii)livestock sites where andlivestock agriculture and were practiced agriculture not prac- were more than ticed1not kmmore awaythanfrom the away 1 km samplingfromcollection the sampling site; and (iii) sites collection where site; there and (iii) was sites somethere where type of pluvial discharge into the river. was some type of pluvial discharge into the river. Figure 1. Schematic diagram diagram of the sampling sampling sites sites in in the the Culiacan Culiacan Valley Valley in northwestern northwestern Mexico. Mexico. Sampling sites A, B, and C were located along the Humaya River, and sampling sites D, E, and F were located along the Tamazula River. River. Site Site G G is located at the point of convergence for both rivers and along with sites H, I, J and K, they belong and along with sites H, I, J and K, they belong to to the the Culiacan Culiacan River. River. TableA1.total Latitude of 10and longitude L of coordinates river water of each sampling was collected in sterilesite for the recovery polypropylene of water from containers per the Humaya River, Tamazula River, and Culiacan River in northwestern Mexico. sampling site each month. Each sample was collected at a depth of 20 cm and at a distance of 50 cm from the edge of the river [27], and all samples were Sampling subjected toSite a collection Samplingschedule, Site River which Sampling Site was performed Location during the last Name days of each month. After collection, samples were Latitude Longitude immediately transported to the Microbiology Laboratories of CIAD’s National Laboratory A Adolfo López Mateos Dam 25°02′46.9″ −107°23′50.5″ for Food Safety [Laboratorio Nacional para la Investigación en Inocuidad Alimentaria at Humaya B Agua Caliente 24°55′44.2″ the Centro de Investigación en Alimentación y Desarrollo (CIAD), Culiacan, −107°23′14.9″ Mexico)] in refrigerated coolersCat a temperature La Guásima of 4–8 25°52′10.2″ ◦ C. The microbiological analysis−107°24′37.1″ of each sample D was performed routinely within 24Sanalona Dam h after collection. 25°41′54.8″ −108°38′41.3″ Tamazula E Imala 24°51′11.7″ −107°13′17.2″ F Las Peñitas 24°51′41.8″ −107°15′21.1″ Microorganisms 2022, 10, 1214 4 of 17 Table 1. Latitude and longitude coordinates of each sampling site for the recovery of water from the Humaya River, Tamazula River, and Culiacan River in northwestern Mexico. Sampling Site Sampling Site River Sampling Site Location Name Latitude Longitude Adolfo López Mateos A 25◦ 020 46.900 −107◦ 230 50.500 Dam Humaya B Agua Caliente 24◦ 550 44.200 −107◦ 230 14.900 C La Guásima 25◦ 520 10.200 −107◦ 240 37.100 D Sanalona Dam 25◦ 410 54.800 −108◦ 380 41.300 Tamazula E Imala 24◦ 510 11.700 −107◦ 130 17.200 F Las Peñitas 24◦ 510 41.800 −107◦ 150 21.100 G Puente Negro 24◦ 480 24.300 −107◦ 240 34.400 H San Pedro 24◦ 470 06.200 −107◦ 330 31.700 Culiacan I Cofradía de San Pedro 24◦ 460 47.300 −107◦ 350 39.800 J Iraguato 24◦ 370 39.200 −107◦ 390 39.200 K El Castillo 24◦ 320 39.700 −107◦ 420 22.800 2.2. Evaluation of Environmental and Water Parameters Several physicochemical parameters were measured at each sampling site in the Culiacan Valley. In particular, water temperature (◦ C), pH, electrical conductivity, total dissolved solids, and salinity were recorded for each sampling site during each collection of river water with a TRACER PockeTester™ monitoring device (Model1766; LaMotte, Chestertown, MD, USA). Additionally, rainfall (mm/m2 ) and relative humidity percentage data from the closest weather station for each sampling point were also recorded using the information provided by the website https://www.tutiempo.net/ (accessed on 1 May 2018), which reports meteorological data from different stations. For the analysis of environmental data, the averages of the measurements during each month of sampling were collected. 2.3. Ultrafiltration Process for River Water Samples To increase the recovery of Salmonella from river water, an ultrafiltration system, based on a size-exclusion method, was employed to efficiently concentrate water volumes from liters to milliliters, resulting in a final suspension of the targeted pathogen in the retentive volume [27]. Briefly, pre-sterilized Masterflex L/S® high-performance precision pump tubing (Cole-Parmer, Vernon Hills, IL, USA) and a Masterflex L/S® Variable-Speed pump system (Model 7528-10; Cole-Parmer) were used for all ultrafiltration experiments with settings to generate a flow rate of 2900 mL/min and to maintain pressure at 5 to 10 psi. Single-use polysulfone dialysis filters (Model Optiflux® F180NR; Fresenius Medical Care, Waltham, MA, USA) were used for each sample after blocking non-specific binding by recirculating 1 L of sterile 0.1% (w/v) sodium hexametaphosphate (NaPP) solution (Sigma- Aldrich, St Louis, MO, USA) for 10 min. After this incubation time, the permeate port was opened to remove the blocking solution and initiate the sample filtration process. Filtration was carried out until a final volume of 400 mL was obtained. The concentrated sample was removed from the ultrafiltration system, and 500 mL of a backflushing solution (0.01% (w/v) NaPP, 0.01% (w/v) Tween 80 and 0.001% (w/v) Antifoam Y-30 emulsion) (Sigma-Aldrich) was added to elute microorganisms until 100 mL of solution was obtained. Finally, the concentrated sample (400 mL) and the elution solution (100 mL) were mixed to obtain 500 mL of final retentate that was used for the recovery of Salmonella isolates. 2.4. Recovery and Characterization of S. enterica Isolates by Culturing, Biochemical and Molecular Methods A standardized analytical protocol was employed for the identification and confirma- tion of non-typhoidal Salmonella isolates using selective media, followed by biochemical and molecular confirmation [28]. For the recovery of presumptive Salmonella isolates, various amounts of each collected surface water sample (undiluted, 1 mL, 10 mL, 20 mL) Microorganisms 2022, 10, 1214 5 of 17 were subjected to an enrichment step in tryptic soy broth (TSB; Becton Dickinson Bioxon, Mexico City, Mexico) by incubating at 37 ◦ C for 24 h. After incubation, the TSB cultures were used for inoculating modified semi-solid Rappaport Vassiliadis (MSRV) agar plates with a total of six discrete 30 µL drops from each TSB enrichment, followed by incubation at 42 ◦ C for 18 h. After incubation on MSRV agar, presumptive Salmonella with a white halo around the colony were selected with a sterile platinum loop and re-streaked on xylose lysine desoxycholate (XLD) agar for further incubation at 37 ◦ C for 24 h. Isolates with a characteristic Salmonella colony morphology (consisting of black colonies with a whitish halo) were selected for subsequent confirmatory biochemical and molecular tests. Biochemical tests on the presumptive Salmonella colonies were subsequently performed using triple sugar iron agar, lysine iron agar, and urea broth and further incubated under aerobic conditions at 37 ◦ C for 24 h, following standardized procedures [28]. The results from the multiple-tube biochemical assays were then employed to obtain an estimation of the bacterial concentration based on the most probable number (MPN) technique and expressed as MPN per 100 mL of sample, according to the protocol for quantitative sample analyses in the EPA Method 1200 [28]. Presumptive Salmonella colonies were subjected to a confirmatory molecular test consisting of polymerase chain reaction (PCR) amplification of the invA gene, which en- codes an invasion protein in Salmonella spp. [29], by using the primers Sal1598 Forward (50 -AACGTGTTTCCGTGCGTAAT-30 ) and Sal1859 Reverse (50 -TCCATCAAATTAGCGGAG GC-30 ) [30]. All oligonucleotides were purchased from Sigma-Aldrich. For preparation of crude lysates for each isolate, TSB cultures were inoculated with Salmonella presumptive colonies grown on XLD agar and incubated for 24 h at 37 ◦ C under aerobic conditions. One milliliter of the overnight TSB culture was subjected to centrifugation at 11,300× g for 5 min, and the cell pellet was washed twice and resuspended in 200 µL of molecular biology-grade water (UltraPure, ThermoFisher, Waltham, MA, USA) for further incubation at 98 ◦ C for 10 min in a ThermoMixer (Eppendorf Latin America, Mexico City, Mexico), followed by 10 min incubation on ice. The cell lysates were subjected to subsequent centrifugation at 6700× g for 5 min, and the supernatants (100–150 µL) were stored at −20 ◦ C until further use as templates for PCR. PCR amplification reactions consisted of 2×GoTaq® Green Master Mix (Promega Corporation, Madison, WI, USA), Sal1598 and Sal1859 primers (1 µM each), and 2 µL of bacterial lysate in a final volume of 25 µL for incubation in a Mastercycler® thermal cycler (Eppendorf Latin America) with the following parameters: 1 cycle at 95 ◦ C for 5 min, followed by 26 cycles of 1 min at 95 ◦ C, 1 min at 63 ◦ C, 1 min at 72 ◦ C, and a final extension cycle at 72 ◦ C for 5 min [30]. The PCR amplicons were separated on a 1% agarose gel and visualized by staining with 1 µL of 10 mg/mL Gel Red® Nucleic Acid Stain (Sigma-Aldrich). Positive samples were selected based on the amplification of a 244 bp fragment of the invA gene [30]. All PCR-positive Salmonella isolates were selected and subsequently preserved at −80 ◦ C in a 50/50 mixture of glycerol and TSB. 2.5. Genome Sequencing and Assembly of the Recovered S. enterica Isolates For whole genome sequencing, the recovered Salmonella isolates were grown for 24 h in TSB at 37 ◦ C under aerobic conditions, and the genomic DNA was extracted using a DNeasy Blood and Tissue kit (QIAGEN, Mexico City, Mexico) following the manufacturer’s instruc- tion, that was initially assessed and quantified with a NanoDrop 2000c Spectrophotometer (Thermo Fisher Scientific, Waltham, MA, USA). For the sequencing reactions, genomic DNA from each isolate was then quantified with a Qubit™ 2.0 Fluorometer (Invitrogen, Carlsbad, CA, USA), adjusted to 0.2 ng/µL, and 1 ng aliquot was used for preparing ge- nomic DNA libraries with the Nextera XT DNA Library Preparation Kit (Illumina Inc., San Diego, CA, USA). Subsequently, the prepared genomic libraries from the Salmonella isolates were sequenced using a MiSeq™ Reagent Kit v2 (300-cycle format) to obtain 2 × 150 bp paired-end reads with a MiSeq™ System (Illumina, Inc.) at the Earlham Institute (Norwich Research Park, Norwich, UK) or a MiniSeq™ System (Illumina, Inc.) at the Molecular Diagnostics Laboratory (Laboratorio de Diagnóstico Molecular, Centro de Investigación Microorganisms 2022, 10, 1214 6 of 17 en Alimentación y Desarrollo, Mazatlán, Mexico) (see Supplementary Table S1). Prior to genome assembly, fastp 0.20.1 [31] was used to assess and filter the quality of the Salmonella genomic sequences and ensure the removal of the trimmed adapter sequences from the 30 -end of the reads. The genome assembly pipelines A5-miseq [32] and SPAdes v3.13.0 [33] were used with default parameters to assemble the genomes for isolates sequenced by the MiSeq™ and MiniSeq™ Systems, respectively. 2.6. In Silico Serotyping and Phylogenetic Relationships of the Recovered S. enterica Isolates The Salmonella in silico typing resource (SISTR) v1.0.2 [34] was employed to perform serovar predictions of the assembled genome sequences for each Salmonella isolate. The default parameters were specified for using a combination of serogroup-specific probes and core genome multilocus sequence typing (cgMLST) analysis for serovar prediction. Analy- sis of the diversity in the detected serovars was assessed by calculating Simpson’s diversity index, where high diversity in the population is indicated by diversity index values close to 1 [35,36]. Furthermore, the relationships of the isolates were investigated through the identification of single nucleotide polymorphisms (SNPs) in the core genome by employing Parsnp v1.2, available in the Harvest software suite [37], where Salmonella enterica subsp. enterica serovar Typhimurium strain LT2 (GenBank Accession No. AE006468.2) was used as the reference strain. A maximum-likelihood phylogenetic tree was constructed using the Randomized Axelerated Maximum Likelihood (RAxML) v8 program [38] with the General- Time-Reversible (GTR) model plus gamma distribution as the nucleotide substitution rate model [39,40] using a bootstrap number of 100. The constructed phylogenetic tree was then visualized and annotated using the Interactive Tree Of Life (iTOL) v6 web-based tool [41]. For conducting the subsequent analysis based on single nucleotide polymorphisms (SNP) phylogeny and multilocus sequence typing (MLST), the assembled genomes were submit- ted to the pipeline CSIPhylogeny version 1.4 [42] and the pipeline multilocus sequence typing (MLST) version 2.0.4 to obtain sequence types (ST) [43] (both available on the web- site for the Center for Genomic Epidemiology http://www.genomicepidemiology.org/ (accessed on 1 March 2022) by using default parameters and S. Typhimurium strain LT2, listed above, as the reference strain. To further characterize the variable and high-quality single nucleotide polymorphisms (SNPs) identified in the recovered S. enterica sequenced genomes, an alignment file consisting of 131,259 SNPs was generated and subsequently visualized with GrapeTree [44] software to generate a minimum spanning tree. The nodes on the minimum spanning tree were numbered by using GrapeTree software and the list of the variable SNPs contributing to the branches are listed in the Supplementary Table S2. Branch lengths of less than 400 SNPs in the minimum spanning tree were collapsed to create a single node per serovar identified in the recovered S. enterica isolates. 2.7. Statistical Analysis All statistical analyses were performed with Minitab® 18 statistical software (Minitab, LLC, State College, PA, USA) to determine the Pearson linear correlation coefficients [14] and Poisson regression analyses [45] for describing the relationships between water and environmental parameters and the prevalence of recovered Salmonella levels during each sampling month. For all measures of association, p values <0.05 were considered statistically significant. 3. Results 3.1. Salmonella Recovery from River Sampling Locations in the Culiacan Valley, Northwestern Mexico The present study employed a size-exclusion ultrafiltration method [24] followed by an enrichment step with selective media [28] to enable improved detection of Salmonella when present at low concentrations in water samples. Presumptive Salmonella spp. isolates were recovered from a total of 54.5% (78/143) of the river water samples collected from 11 river sampling sites in the Culiacan Valley from June 2018 to November 2019 (Figure 1 3.1. Salmonella Recovery from River Sampling Locations in the Culiacan Valley, Northwestern Mexico The present study employed a size-exclusion ultrafiltration method [24] followed by an enrichment step with selective media [28] to enable improved detection of Salmonella Microorganisms 2022, 10, 1214 when present at low concentrations in water samples. Presumptive Salmonella spp.7 of iso- 17 lates were recovered from a total of 54.5% (78/143) of the river water samples collected from 11 river sampling sites in the Culiacan Valley from June 2018 to November 2019 (Figure and Table 1 and Table 1). The 1). The highest highest percentage percentage of positive of positive water water samples wassamples obtained was obtained from sites B from sites B and D with 69.2% (9/13), followed by sites A, C, E, G, and I with and D with 69.2% (9/13), followed by sites A, C, E, G, and I with 61.5% (8/13) of positive 61.5% (8/13) of positiveSites samples. samples. J andSites J and K had theKlowest had theprevalence lowest prevalence with 30.8%with(4/13) 30.8%and(4/13) and (3/13) 23.1% 23.1% (3/13) positivepositive samples,samples, respectively. respectively. BasedBased on theonnumber the number of combined of combined positive positive samples samples per per river, river, the highest the highest number number of positive of positive samples samples was was identified identified for sites for sites along along the Hu- the Humaya maya River64.1% River with with (25/39), 64.1% (25/39), followedfollowed by the by the Tamazula Tamazula River River with with(24/39), 61.5% 61.5% (24/39), and thenandby then by the Culiacan River with 44.6% (29/65) the Culiacan River with 44.6% (29/65) of the samples. of the samples. As shown in Figure 2, quantitation of Salmonella in the collected river water samples was determined by employing the most probable number technique [28]. A spatial exami- exam- ination nation of of the the recovered recovered Salmonella Salmonella levels revealed higher than average recovery levels for for sites B and H on the Humaya sites Humaya and Culiacan Culiacan rivers, rivers, followed followed byby sites sites A, A, D, D, and and FF on on the Humaya and Humaya and Tamazula rivers (Figure 2a). The lowest lowest detected detected levels Salmonella were levels of Salmonella observed for Site K on the Culiacan Culiacan River River prior prior toto discharging into the Pacific Ocean (Sea of Cortez). of Cortez). Furthermore, Furthermore, temporal temporal analysis analysis of of the the recovered Salmonella showed recovered Salmonella showed thatthat the detected levels detected levels were highly highly variable variable but significantly higher during the summer months when compared to other months of sampling (p == 0.000) 0.000) (Figure (Figure 2b). 2b). Throughout Throughout the the entire entire sampling period, sampling period, the the detectable detectable levels levels of Salmonella in of Salmonella in the river water samples had mean values ranging values ranging from from <0.64 <0.64 MPN/100 MPN/100mL mLto to13.98 13.98MPN/100 MPN/100mL mL(Figure (Figure2b). 2b). (a) (b) Figure Figure 2.2. Estimation EstimationofofSalmonella Salmonellalevels recovered levels from recovered river from water river samples. water A standardized samples. ana- A standardized lytical protocol employing culturing, biochemical, and molecular confirmation [28] was used to de- analytical protocol employing culturing, biochemical, and molecular confirmation [28] was used to termine the Salmonella levels based on the most probable number (MPN) technique per 100 mL of determine the Salmonella levels based on the most probable number (MPN) technique per 100 mL sample for the various sites (a) and dates (b). The red dashed line represents the overall mean of of sample for Salmonella thefor levels various sites (a)sites all sampling and and datesdates. (b). The red dashed line represents the overall mean of Salmonella levels for all sampling sites and dates. In particular, the highest levels of Salmonella were detected during the summer In particular, the highest levels of Salmonella were detected during the summer months of June, July and August of 2018 with average value of 13.8 MPN/100 mL, 7.88 months of June, July and August of 2018 with average value of 13.8 MPN/100 mL, MPN/100 mL and 6.75 MPN/100 mL, respectively (Figure 2b). The lowest levels were ob- 7.88 MPN/100 mL and 6.75 MPN/100 mL, respectively (Figure 2b). The lowest levels served from December 2018 to March 2019 and September 2019 to November 2019, corre- were observed from December 2018 to March 2019 and September 2019 to November 2019, sponding to thetofall corresponding theand fallwinter seasons and winter in thein seasons Culiacan ValleyValley the Culiacan with levels ranging with levels from ranging 0.65 from 0.65 MPN/100 mL to 0.93 MPN/100 mL (Figure 2b). During the summer of 2019,Sal- MPN/100 mL to 0.93 MPN/100 mL (Figure 2b). During the summer of 2019, the the monella levels were 2.41 and 2.33 MPN/100 mL in July and August, respectively. Salmonella levels were 2.41 and 2.33 MPN/100 mL in July and August, respectively. These These recovery recovery levels levels were were three-fold three-fold lower lower when when compared compared to to those those detected detected in in the the summer summer months of the previous year but at least two-fold higher than the levels detected during the winter and fall of 2019. 3.2. Effect of Environmental and Physicochemical Water Parameters on Salmonella Detection Levels The effects of several environmental parameters on the recovery levels of Salmonella from the river water were further evaluated. During the sampling period in the present study, the seasonal values for rainfall ranged from 3.0 to 194.6 mm/m2 , and the recorded relative humidity fluctuated from 51.09 to 75.5% (Supplementary Table S3). The lowest values for these environmental parameters were recorded in the fall and winter months Microorganisms 2022, 10, 1214 8 of 17 (October to April), whereas the highest values were registered in summer (June to Septem- ber). It is noteworthy that rainfall levels exceeded 214 mm/m2 in September 2018 due to Tropical Depression 19-E, which caused flooding throughout the Culiacan Valley [46]. Further analysis of physicochemical water parameters revealed that the river water had temperatures ranging from 21.8 ◦ C to 36.4 ◦ C (Supplementary Table S3). Additional ex- amination of other physicochemical water parameters for sites A to I revealed similar values throughout the year for electrical conductivity (<613.46 µS/cm), total dissolved solids (<387 ppm), and salinity (<0.29 PSU). As expected, the farthest site, K, located on the Culiacan River and in close proximity to the Pacific Ocean, had the highest mean values for electrical conductivity (6477 µS/cm), total dissolved solids (6066 ppm), and salinity (9 PSU) (Supplementary Table S3). When compared to all other sampling sites, site J on the Culiacan River showed the greatest fluctuations in the recorded levels for electrical conductivity (4130–6400 µS/cm), total dissolved solids (2510–4950 ppm), and salinity (1.9–2.76 PSU) during the months of January to July. According to Pearson’s analysis, the physicochemical parameters, such as total dissolved solids, electrical conductivity, and salinity showed a high positive correlation, indicating a statistically significant (p < 0.05) linear relationship with each other (Table 2). Moreover, rainfall and relative humidity also had a significant positive correlation (r = 0.762, p = 0.000). Further analysis of the correlations for the various environmental parameters with the detected Salmonella levels revealed that only the river water temperature showed a significant positive correlation (r = 0.204, p = 0.017) (Table 2): therefore, at higher or lower temperatures, Salmonella levels increased or decreased, respec- tively. By contrast, the river water pH was found to be negatively correlated with the levels of Salmonella in the Culiacan Valley (r = −0.231, p = 0.007). To further quantify any significant correlations between the physicochemical param- eters and the enumerated Salmonella spp. from river water, a Poisson regression model was performed (Table 3). The results indicated that only the environmental data recorded during June, July, and August were significantly associated with an increase in Salmonella levels (p = 0.000). In particular, the month of June in the first year of sampling showed the highest levels of Salmonella, and when compared to the other months of sampling, the levels of Salmonella in June were observed to have a significant log increase of e2.845 = 17.2 (95% confidence interval (CI): 12.94–22.87) (Table 3). In the month of July for both sampling years, a log increase in Salmonella concentrations was calculated as e1.371 = 3.94 (95% CI: 2.8–5.53), followed by a gradual decrease in the month of August with an average level of e1.303 = 3.68 (95% CI: 2.61–5.21). Furthermore, the results indicated that for each unit increase of pH and salinity, the average of concentration of Salmonella spp. decreased by e−0.397 = 0.67 (95% CI: 0.487–0.923) and e−0.115 = 0.89 (95% CI: 0.835–0.951), respectively (Table 3). Based on these quantitative analyses, only two physicochemical water parame- ters, pH and salinity, were found to have a significant quantitative negative effect on the detected levels of Salmonella recovered from river water in the Culiacan Valley. Microorganisms 2022, 10, 1214 9 of 17 Table 2. Pearson’s correlation coefficients (r) and p-values (p) determined for the environmental parameters and Salmonella detection levels for all sampling sites throughout the sampling period. Environmental Salmonella River Water Total Dissolved Electrical pH Salinity Rainfall Parameter Levels a Temperature Solids Conductivity r b p r p r p r p r p r p r p pH −0.231 0.007 River water temperature 0.204 0.017 0.173 0.044 Total dissolved solids −0.143 0.097 0.219 0.010 0.100 0.246 Salinity −0.133 0.123 −0.236 0.006 0.106 0.221 0.840 0.000 Electrical conductivity −0.143 0.097 0.215 0.012 0.095 0.269 0.985 0.000 0.861 0.000 Rainfall 0.056 0.518 −0.012 0.890 0.767 0.000 −0.030 0.731 −0.027 0.753 −0.035 0.685 Relative humidity −0.049 0.572 0.069 0.426 0.561 0.000 0.045 0.605 0.048 0.582 0.036 0.681 0.762 0.000 a Estimation of Salmonella levels based on the most probable number (MPN) technique per 100 mL of sample [28]. b Pearson’s correlation coefficients displayed in bold are statistically significant (p-value < 0.05). Table 3. Parameters examined with Poisson regression analysis for Salmonella detection levels recovered from river water in northwestern Mexico. Regression Environmental Regression Standard Error Variance Inflation Coefficient 95% Z-Value p-Value Parameter Coefficient (ex ) Coefficient Factor Confidence Interval (CI) Constant 2.954 1.264 (0.48, 5.43) 2.34 0.019 Salmonella levels in June a 2.845 0.147 (2.56, 3.13) 19.40 0.000 1.80 Salmonella levels in July 1.371 0.175 (1.03, 1.71) 7.86 0.000 1.53 Salmonella levels in August 1.303 0.177 (0.96, 1.65) 7.35 0.000 1.52 pH −0.397 0.163 (−0.72, −0.08) −2.43 0.015 1.11 Salinity −0.115 0.035 (−0.18, −0.05) −3.29 0.001 1.07 a Estimation of Salmonella levels based on the most probable number (MPN) technique per 100 mL of sample [28]. Microorganisms 2022, 10, 1214 10 of 17 Microorganisms 2022, 10, x FOR PEER REVIEW 10 of 17 3.3. Serovar Identification and Phylogenetic Relationships among the Recovered Salmonella Isolates 3.3. fromSerovar Identification River Water and Phylogenetic in the Culiacan Valley Relationships among the Recovered Salmonella Isolates from River Water in the Culiacan Valley To perform a genotypic classification of the Salmonella isolates recovered from the river To sampling perform sites in athe genotypic Culiacanclassification Valley, a totalofofthe Salmonellaisolates 73 Salmonella isolateswererecovered fromtothe subjected whole river sampling sites in the Culiacan Valley, a total of 73 Salmonella genome sequencing (Supplementary Tables S1 and S4). The sequencing results revealed isolates were subjected to whole genome sequencing (Supplementary Tables S1 and S4). The sequencing results that 100% (73/73) of the recovered isolates were classified as S. enterica subspecies enterica. revealed that 100% (73/73) of the recovered isolates were classified as S. enterica subspecies Subsequent in silico identification of serovar-specific genes in S. enterica subspecies enterica enterica. Subsequent in silico identification of serovar-specific genes in S. enterica subspe- revealed 27 distinct serovars, and the isolates belonging to serovar Oranienburg were cies enterica revealed 27 distinct serovars, and the isolates belonging to serovar Oranien- the most prevalent, representing 23% (17/73) of the recovered isolates. Moreover, lower burg were the most prevalent, representing 23% (17/73) of the recovered isolates. Moreo- prevalences were detected for S. Anatum (8% (6/73)), S. Saintpaul (7% (5/73)), S. Pomona ver, lower prevalences were detected for S. Anatum (8% (6/73)), S. Saintpaul (7% (5/73)), (5% (4/73)), and S. Sandiego (5% (4/73)). S. Oranienburg had the highest distribution S. Pomona (5% (4/73)), and S. Sandiego (5% (4/73)). S. Oranienburg had the highest distri- and was detected in 90.9% (10/11) of the sampling sites, followed by S. Anatum, and bution and was detected in 90.9% (10/11) of the sampling sites, followed by S. Anatum, S. Sandiego, which were detected in 4 of the 11 sampling sites (Table S3). The highest and S. Sandiego, which were detected in 4 of the 11 sampling sites (Table S3). The highest diversity in serovar distribution was detected in samples collected from the Tamazula and diversity in serovar distribution was detected in samples collected from the Tamazula and Culiacan Rivers, with a total of 16 serovars identified. The lowest diversity was detected Culiacan Rivers, with a total of 16 serovars identified. The lowest diversity was detected in the S. enterica isolates recovered from the Humaya River, with 11 distinct serovars. The in the S. enterica isolates recovered from the Humaya River, with 11 distinct serovars. The diversity in diversity in the the recovered recovered serovars serovars was was examined examined by by calculating calculating Simpson’s Simpson’s diversity diversity in- index, and the dex, and results from the results thisthis from quantification quantification revealed revealed high highdiversity diversityindex indexvalues values of of 0.92 0.92 for the Humaya River, 0.95 for the Tamazula River, and 0.90 for the for the Humaya River, 0.95 for the Tamazula River, and 0.90 for the Culiacan River. Fur-Culiacan River. Further quantification ther quantificationof Simpson’s of Simpson’sdiversity index diversity for the index serovars for the serovarsrecovered recovered from fromall all three rivers three rivers revealed a diversity index of 0.91 during the summer months of June, July and Au- as revealed a diversity index of 0.91 during the summer months of June, July and August well as well gust as a highas diversity index of a high diversity 0.95offor index thefor 0.95 remaining non-summer the remaining non-summer months. months. To further analyze the genetic relatedness of the isolates, the S. To further analyze the genetic relatedness of the isolates, the S. enterica genomes enterica genomes werewere subjected to MLST as well as a phylogenetic analysis based on the identification of eacheach subjected to MLST as well as a phylogenetic analysis based on the identification of isolate’s core isolate’s coregenome genomeofofSNPs SNPsasasa amethod method forfor achieving achieving high high resolution resolution genotyping genotyping andand for classifying for classifyingisolates isolatesbased basedonondifferent differentpolymorphisms polymorphisms (Figure (Figure 3).3). Figure Figure 3. 3. Phylogenetic Phylogenetictree treeofofS.S.enterica serovars enterica isolated serovars from isolated rivers from in the rivers Culiacan in the Valley, Culiacan north- Valley, north- western Mexico. western Mexico.AAmaximum-likelihood maximum-likelihoodphylogenetic phylogenetic tree was tree constructed was from constructed corecore from genome SNPs genome SNPs using using RAxML RAxMLv8 v8[38] [38]under underthe theGeneral-Time-Reversible General-Time-Reversible plus plusgamma gamma distribution model distribution withwith model a a boot strap value of 100. The constructed phylogenetic tree was visualized and annotated using the boot strap value of 100. The constructed phylogenetic tree was visualized and annotated using the Interactive Tree of Life (iTOL) v6 web-based tool [41], and two clades are indicated by the red and Interactive blue Tree of Life (iTOL) v6 web-based tool [41], and two clades are indicated by the red and branches. blue branches. Microorganisms Microorganisms 2022, 10, x1214 FOR PEER REVIEW 11 11 of of 17 17 The The typing typing analysis analysis by by MLST MLST demonstrated demonstrated that that the the ST ST ofof each each isolate isolate corresponded corresponded to those to those belonging belonging to to the the same same serovar serovar except except for for S. Minnesota, Minnesota, wherewhere one isolate isolate had had anan ST of 285 ST of 285 andand the other isolate other isolate had had an ST of 548 (Table S2). Moreover, the SNP analysis analysis revealed that revealed that the the recovered enterica isolates recovered S. enterica isolates from from river river water water could could bebe classified classified in in two two major clades major clades(Figure (Figure3). 3).Most Mostof of thethe recovered recovered S. enterica S. enterica isolates isolates clustered clustered in clade in clade 1, 1, sig- significantly nificantly represented represented by S. S. Oranienburg, byOranienburg, the predominant the predominant serovar serovar identified identified in theinpre- the present sent study.study. CladeClade 2 consisted 2 consisted of only of only 38% 38%of allofisolates, all isolates, mostlymostly recovered recovered fromfromthe Hu-the Humaya River, where Anatum and Saintpaul were the most maya River, where Anatum and Saintpaul were the most common serovars (Figure 3). common serovars (Figure 3). Sequence analysis of the core genome revealed S. enterica isolates Sequence analysis of the core genome revealed S. enterica isolates belonging to the same belonging to the same serotype and serotype andrecovered recoveredfrom fromdifferent differentrivers rivers were were closely closely related related andandonlyonly differed differed by 400by 400 SNPs. SNPs. By contrast, By contrast, isolates isolates belonging belonging to to different different serovarswithin serovars withinthethesame samecladeclade were were found to found to vary vary byby approximately approximately 15,000 15,000 toto 20,000 20,000 SNPs. SNPs. OneOne notable notable exception exception concerns concerns the S. Bovismorbificans the S. Bovismorbificansisolates,isolates,belonging belongingtoto clade clade 1, which 1, which werewere foundfound to differ to differ by ap-by approximately 30,000 SNPs from those isolates belonging proximately 30,000 SNPs from those isolates belonging to both clades. to both clades. As shown As shown in in Figure Figure4,4,a asubsequent subsequent validation validation of of all all SNPsSNPsdetected in the detected in core genome the core ge- for all sequenced S. enterica isolates was conducted with a novel nome for all sequenced S. enterica isolates was conducted with a novel SNP analysis SNP analysis method employing method the CSIPhylogeny employing pipeline (Kaas the CSIPhylogeny 2014). pipeline In agreement (Kaas 2014). In with the corewith agreement genome SNP the core analysis SNP genome (Figure 3), results analysis (Figurefrom3),the SNP-based results from thephylogeny SNP-based revealed phylogenythat revealed both synonymous that both and nonsynonymous SNPs in the sequenced S. enterica synonymous and nonsynonymous SNPs in the sequenced S. enterica isolates isolates clustered in two separate clustered in two separate clades (Figure 4). Serovars belonging to clade 1 were more closelyof clades (Figure 4). Serovars belonging to clade 1 were more closely related than those related in clade 2 (Figure 4a). All three sampled rivers were found to be well represented in both than of those in clade 2 (Figure 4a). All three sampled rivers were found to be well repre- clades (Figure 4b). sented in both clades (Figure 4b). (a) (b) Figure Figure 4. 4. Minimum spanning tree Minimum spanning treeof ofSNPs SNPsidentified identifiedininthe theS.S. enterica enterica isolates isolates from from river river water water sam- samples ples in the Culiacan Valley, northwestern Mexico. An alignment of highly variable SNPs identified in the Culiacan Valley, northwestern Mexico. An alignment of highly variable SNPs identified in the in the genomes of the recovered S. enterica isolates was generated using CSIPhylogeny version 1.4 genomes of the recovered S. enterica isolates was generated using CSIPhylogeny version 1.4 [42]. The [42]. The alignment was visualized using GrapeTree software with the MSTree V2 algorithim [44] alignment to generatewas visualized a minimum using GrapeTree spanning tree. SNPssoftware with the clusters were MSTree colored V2 algorithim according [44] to generate to the detected serovar a minimum spanning tree. SNPs clusters were colored according to the detected serovar (a) or the sampled river (b) in the Culiacan Valley. Nodes less than 400 SNPs apart were collapsed (a) or the sampled to form a river single(b) in the SNP Culiacan cluster, Valley. and the nodeNodes size is less than 400 indicative of SNPs apart were the number collapsed to form a of isolates. single SNP cluster, and the node size is indicative of the number of isolates. From the analysis of 600 SNPs that distinguished clade 1 from clade 2, a total of 62 From the analysis nonsynonymous SNPs,ofaffecting 600 SNPsprotein that distinguished clade function, were 1 from clade identified 2, a total of 62 (Supplementary non- Table synonymous SNPs, affecting protein function, were identified (Supplementary S2). Further analysis of the detected nonsynonymous SNPs revealed that the vast major- Table S2). Further analysis of the detected nonsynonymous SNPs revealed that the vast majority, ity, corresponding to 48% (30/62) of the polymorphisms, reside in genes conferring viru- cor- responding to 48% (30/62) of the polymorphisms, reside in genes conferring virulence traits Microorganisms 2022, 10, 1214 12 of 17 associated with adhesion, motility, and colonization of the host, as well as survival from nutritional and oxidative stress (Supplementary Table S3). Other SNPs were found in genes responsible for metabolic traits, motility, and membrane function (Supplementary Table S3). Interestingly, five of the detected SNPs affected stiA, a starvation-inducible loci associated with survival during prolonged nutritional starvation due to deficiencies of phosphate, carbon, and nitrogen [47]. Other multiple SNPs were identified in the virulence genes ratA and sivH, both located on the Salmonella pathogenicity island CS54 and required for intestinal colonization in the host [48]. 4. Discussion The present study examined the prevalence and diversity of S. enterica serovars recov- ered from major rivers in the Culiacan Valley, a major fresh produce production region in northwestern Mexico. Subsequent analysis evaluated physicochemical water parameters and their influence on the prevalence of Salmonella in rivers in the Culiacan Valley. This fundamental information will aid the identification of potential contamination sources for the development of focused and efficient intervention strategies for mitigating outbreaks caused by this foodborne pathogen. In Mexico, the National Epidemiological Surveillance System (Sistema Nacional de Vigilancia Epidemiológica, SINAVE) demonstrated a consid- erable increase in cases of salmonellosis from 1984 to 2017 (from 31,943 to 104,471), and the state of Sinaloa had the highest rate corresponding to 12.9% of cases [49]. Although the routes of exposure to Salmonella in this agricultural region have not been well defined and are often epidemiologically unrelated to specific sources [18,20,27], several studies have suggested that surface waters in this agricultural region could potentially serve as the source of Salmonella contamination of produce. Therefore, tracing Salmonella in surface waters is of utmost importance, not only for human health, but also for the safety of the fresh produce industry in this region [18,27,50]. The Adolfo López Mateos and Sanalona dams supply water to three rivers in the Culiacán Valley region, where it is used for agricultural, recreation, and domestic purposes. Some surrounding communities do not have efficient wastewater treatment, leading to wastewater mixing with surrounding waters. Additionally, agriculture, livestock, and wildlife animals, which are prevalent in this region, can be potential sources of Salmonella, where the main source of water for this agricultural region are provided by the Humaya, Tamazula and Culiacan Rivers [51,52]. In the present study, the effect of environmental factors on the detected levels of Salmonella in this agricultural region were examined by using an ultrafiltration system to increases the probability and percentage of Salmonella recovery, which can have a great impact on the analysis of Salmonella levels [27]. Recent findings have indicated that the prevalence and survival of Salmonella in surface waters may be due, in part, to the temperature, rainfall, as well as the physicochemical conditions of the water, and even the type of hosts that coexist in this type of habitat [53,54]. In agreement with previous studies demonstrating that salmonellosis cases reaches its peak in the summer months worldwide [27,55–57], the present study reported that the highest levels of Salmonella were observed in the summer months (June, July, and August) with concentrations significantly higher than the other months. Among the environmental parameters examined in this study, only temperature was positively correlated with Salmonella levels, indicating warm temperatures may facilitate higher environmental persistence of the pathogen by promoting an increase in metabolic activity during the summer months [45,54,58,59]. By contrast, the levels of precipitation did not positively correlate with the detected levels of Salmonella, and the concentrations of Salmonella gradually decreased as a consequence of an increase in precipitation levels. Recent evidence has shown rainfall can affect Salmonella levels by causing the dilution of microbial populations [14,60–62]. During the survey in this study, precipitation levels in the Culiacan Valley reached 214.2 mm/m2 due to the presence of storms [46,51], and the high amounts of rainfall may have significantly affected the prevalence of Salmonella during the following months. Among the other physicochemical parameters for the river water, only Microorganisms 2022, 10, 1214 13 of 17 pH and salinity showed a significantly negative association with Salmonella levels. The sampling sites J and K, which were in close proximity to the Pacific Ocean, had the lowest Salmonella levels and the highest levels of salinity. In agreement with previous studies, high levels of salinity harbor a low prevalence of Salmonella [22,59,63,64]. In this survey, a total of 27 S. enterica serovars were recovered from the rivers in the Culiacan Valley. Further quantitative analysis determined a Simpson’s diversity index with combined values higher that 0.92, revealing that a high diversity of serovars was found in this survey. The diversity identified in this agricultural region in Mexico is significantly higher than that reported for other major agricultural regions for fresh produce in the United States, which have an overall index of 0.89 [65]. Although the clinically associated serovar Typhimurium was recovered at low levels during the sampling period in the Culiacan Valley, S. Oranienburg, S. Anatum, S. Saintpaul, and S. Pomona were commonly identified. Interestingly, S. Oranienburg was consistently recovered from most sampling sites and predominantly detected during the summer months. These findings are in agreement with previous studies reporting that S. Oranienburg was one of the most prevalent serovars in aquatic environments and animal hosts in northwest Mexico [21,27,50]. Recently, S. Oranienburg was identified as the causative agent in a salmonellosis outbreak in Mexico associated with the consumption of fresh produce [66], highlighting the need for close monitoring of this serovar in agricultural regions. According to the SNP phylogenetic analysis conducted in the present study, S. enterica isolates belonging to same serovar clustered together by MLST and had similar phylogenies according to the analysis of core genome SNPs and highly variable SNPs. Recently, a report on S. enterica isolates recovered from surface waters in central Mexico identified significantly higher intraserovar diversity based on the CSIPhylogeny analysis of highly variable SNPs [23]. Moreover, SNP phylogeny identified that some isolates belonging to the same serovar were grouped in different SNP clusters [23]. By contrast, the S. enterica isolates from river water in the Culiacan Valley in northwestern Mexico were closely related with lower intraserovar diversity and were found to have less than 400 SNP if they belonged to the same serovar. Importantly, isolates belonging to different serovars were found to differ by over 10,000 SNPs. To further characterize the differences among serovars, the SNP-based phylogeny for the recovered S. enterica isolates identified several nonsynonymous polymorphisms in the virulence gene stiA, a locus that is rpoS-dependent and required for stress situations, including thermotolerance, acid tolerance, and survival from carbon starvation [47]. Other SNPs were identified in sinH, also known as sivH, and ratA, which are both found in the 25-kb genetic island at centisome 54 (CS54 island) in S. Typhimurium and are required for persistent colonization of the animal host [48,67]. Future work will be aimed at characterizing relevant markers contributing to the fitness and survival of S. enterica in aquatic environments in the Culiacan Valley, an important agricultural region in northwestern Mexico. The recovery of S. enterica serovars from agricultural rivers warrants a more detailed gene expression analysis that will lead into the identification of fitness traits conferring increased survival in surface water habitats. Supplementary Materials: The following supporting information can be downloaded at: https: //www.mdpi.com/article/10.3390/microorganisms10061214/s1, Table S1: Genome characteristics of the Salmonella enterica isolates recovered from river waters in the Culiacan Valley, Mexico; Table S2: List of nonsynonymous single nucleotide polymorphisms identified in the sequenced Salmonella enterica isolates recovered from river water in the Culiacan Valley; Table S3: Summary of Salmonella levels, water and environment parameters for each sampling month; Table S4: Distribution of recovered Salmonella enterica serovars by sampling site in the Culiacan Valley, Northwestern Mexico. Author Contributions: Conceptualization, J.A.M.-F. and J.P.G.-G.; data curation, I.G.-L., B.G.L. and B.Q.; formal analysis, I.G.-L., J.B.V.-T., B.G.L. and B.Q.; funding acquisition, B.Q. and C.C.; investiga- tion, I.G.-L. and J.A.M.-F.; methodology, J.A.M.-F., J.P.G.-G., B.G.L. and B.Q.; project administration, C.C.; resources, B.Q. and C.C.; software, J.P.G.-G., J.B.V.-T., B.G.L. and B.Q.; supervision, C.C.; valida- tion, J.B.V.-T., B.G.L. and B.Q.; visualization, I.G.-L.; writing—original draft, I.G.-L.; writing—review & editing, J.A.M.-F., N.C.-d.C., O.L.-C., J.P.G.-G., J.R.A.-S., J.M.-U., B.G.-G., B.G.L., B.Q. and C.C. All Microorganisms 2022, 10, 1214 14 of 17 authors will be informed about each step of manuscript processing including submission, revision, revision reminder, etc. via emails from our system or assigned Assistant Editor. All authors have read and agreed to the published version of the manuscript. Funding: This material was based in part upon work supported by the Postgraduate Studies Scholar- ship from The National Council of Science and Technology in Mexico (CONACyT Grant No. 531137) to Irvin González-López, by the National Laboratory for Food Safety Research (LANIIA, CIAD), CONACyT Project 262200, and by the United States Department of Agriculture (USDA), Agricultural Research Service, CRIS Project 2030-42000-055-00D. Institutional Review Board Statement: Not applicable. Informed Consent Statement: Not applicable. Data Availability Statement: The complete genome sequences of the 73 S. enterica isolates have been deposited with links to BioProject accession number PRJNA832240 in the National Center for Biotechnology Information (NCBI) BioProject database (https://www.ncbi.nlm.nih.gov/bioproject/ (accessed on 1 May 2022). Accession numbers for each S. enterica isolate are found in the Supplemen- tary Table S1. Acknowledgments: The authors would like to thank Celida Isabel Martínez-Rodríguez (CIAD, Culiacan, Mexico) for providing excellent technical support. Mention of trade names or commercial products in this publication is solely for the purpose of providing specific information and does not imply recommendation or endorsement by the USDA. USDA is an equal opportunity provider and employer. Conflicts of Interest: The authors declare no conflict of interest. The funders had no role in the design of the study; in the collection, analyses, or interpretation of data; in the writing of the manuscript, or in the decision to publish the results. References 1. Besser, J.M. Salmonella epidemiology: A whirlwind of change. Food Microbiol. 2018, 71, 55–59. [CrossRef] [PubMed] 2. Majowicz, S.E.; Musto, J.; Scallan, E.; Angulo, F.J.; Kirk, M.; O’Brien, S.J.; Jones, T.F.; Fazil, A.; Hoekstra, R.M.; International Collaboration on Enteric Disease “Burden of Illness” Studies. The global burden of nontyphoidal Salmonella gastroenteritis. Clin. Infect. Dis. 2010, 50, 882–889. [CrossRef] [PubMed] 3. EFSA. European Food Safety Authority and European Centre for Disease Prevention Control. The European Union One Health 2019 Zoonoses Report. EFSA J. 2021, 19, e06406. [CrossRef] 4. Godínez-Oviedo, A.; Tamplin, M.L.; Bowman, J.P.; Hernández-Iturriaga, M. Salmonella enterica in Mexico 2000–2017: Epidemiology, Antimicrobial Resistance, and Prevalence in Food. Foodborne Pathog. Dis. 2020, 17, 98–118. [CrossRef] [PubMed] 5. Clemens, R.L. The Expanding U.S. Market for Fresh Produce; Center for Agricultural and Rural Development at Digital Repository Iowa State University: Ames, IA, USA, 2015; pp. 1–3. 6. Flores, D. Tomato Annual; MX8025; USDA Foreign Agricultural Service: Washington, DC, USA, 2018. 7. Huang, S.W. Imports Contribute to Year-Round Fresh Fruit Availability; FTS-356-01; U.S. Department of Agriculture, Economic Research Service: Washington, DC, USA, 2013. 8. Aguirre-Sanchez, J.R.; Ibarra-Rodriguez, J.R.; Vega-Lopez, I.F.; Martínez-Urtaza, J.; Chaidez-Quiroz, C. Genomic signatures of adaptation to natural settings in non-typhoidal Salmonella enterica Serovars Saintpaul, Thompson and Weltevreden. Infect. Genet. Evol. 2021, 90, 104771. [CrossRef] [PubMed] 9. Andino, A.; Hanning, I. Salmonella enterica: Survival, colonization, and virulence differences among serovars. Sci. World J. 2015, 2015, 520179. [CrossRef] 10. Cheng, R.A.; Eade, C.R.; Wiedmann, M. Embracing Diversity: Differences in Virulence Mechanisms, Disease Severity, and Host Adaptations Contribute to the Success of Nontyphoidal Salmonella as a Foodborne Pathogen. Front. Microbiol. 2019, 10, 1368. [CrossRef] 11. Liu, H.; Whitehouse, C.A.; Li, B. Presence and persistence of Salmonella in water: The impact on microbial quality of water and food safety. Front. Public Health 2018, 6, 159. [CrossRef] 12. Gonzalez-Lopez, I.; Medrano-Felix, J.A.; Castro-Del Campo, N.; Lopez-Cuevas, O.; Ibarra Rodriguez, J.R.; Martinez-Rodriguez, C.; Valdez-Torres, J.B.; Chaidez, C. Metabolic plasticity of Salmonella enterica as adaptation strategy in river water. Int. J. Environ. Health Res. 2021, 31, 1896682. [CrossRef] 13. Medrano-Félix, J.A.; Castro-del Campo, N.; Peraza Garay, F.d.J.; Martínez-Rodríguez, C.I.; Chaidez, C. Carbon source utilization- based metabolic activity of Salmonella Oranienburg and Salmonella Saintpaul in river water. Water Environ. J. 2018, 32, 118–124. [CrossRef] Microorganisms 2022, 10, 1214 15 of 17 14. Haley, B.J.; Cole, D.J.; Lipp, E.K. Distribution, diversity, and seasonality of waterborne salmonellae in a rural watershed. Appl. Environ. Microbiol. 2009, 75, 1248–1255. [CrossRef] [PubMed] 15. Harris, C.S.; Tertuliano, M.; Rajeev, S.; Vellidis, G.; Levy, K. Impact of storm runoff on Salmonella and Escherichia coli prevalence in irrigation ponds of fresh produce farms in southern Georgia. J. Appl. Microbiol. 2018, 124, 910–921. [CrossRef] [PubMed] 16. Ho, Y.N.; Tsai, H.C.; Hsu, B.M.; Chiou, C.S. The association of Salmonella enterica from aquatic environmental and clinical samples in Taiwan. Sci. Total Environ. 2018, 624, 106–113. [CrossRef] [PubMed] 17. Grimont, P.A.; Weill, F.-X. Antigenic Formulae of the Salmonella Serovars, 9th ed.; WHO Collaborating Centre for Reference and Research on Salmonella: Paris, France, 2007; pp. 1–166. 18. Contreras-Soto, M.B.; Medrano-Félix, J.A.; Ibarra-Rodríguez, J.R.; Martínez-Urtaza, J.; Chaidez, Q.C.; Castro-del Campo, N. The last 50 years of Salmonella in Mexico: Sources of isolation and factors that influence its prevalence and diversity. Rev. Bio Cienc. 2019, 6, e540. [CrossRef] 19. Estrada-Acosta, M.; Jiménez, M.; Chaidez, C.; León-Felix, J.; Castro-Del Campo, N. Irrigation water quality and the benefits of implementing good agricultural practices during tomato (Lycopersicum esculentum) production. Environ. Monit. Assess. 2014, 186, 4323–4330. [CrossRef] 20. Jiménez, M.; Martínez-Urtaza, J.; Chaidez, C. Geographical and temporal dissemination of salmonellae isolated from domestic animal hosts in the Culiacan Valley, Mexico. Microb. Ecol. 2011, 61, 811–820. [CrossRef] 21. López-Cuevas, O.; León-Félix, J.; Jiménez-Edeza, M.; Chaidez-Quiroz, C. Detección y resistencia a antibióticos de Escherichia coli y Salmonella en agua y suelo agrícola. Rev. Fitotec. Mex. 2009, 32, 119–126. [CrossRef] 22. Martínez-Urtaza, J.; Liebana, E.; García-Migura, L.; Pérez-Pineiro, P.; Saco, M. Characterization of Salmonella enterica serovar typhimurium from marine environments in coastal waters of Galicia (Spain). Appl. Microbiol. 2004, 70, 4030–4034. [CrossRef] 23. Ballesteros-Nova, N.E.; Sánchez, S.; Steffani, J.L.; Sierra, L.C.; Chen, Z.; Ruíz-López, F.A.; Bell, R.L.; Reed, E.A.; Balkey, M.; Rubio-Lozano, M.S.; et al. Genomic epidemiology of Salmonella enterica circulating in surface waters used in agriculture and aquaculture in Central Mexico. Appl. Environ. Microbiol. 2022, 88, e02149-21. [CrossRef] 24. Jiménez, M.; Chaidez, C. Improving Salmonella determination in Sinaloa rivers with ultrafiltration and most probable number methods. Environ. Monit. Assess. 2012, 184, 4271–4277. [CrossRef] 25. McEgan, R.; Rodrigues, C.A.; Sbodio, A.; Suslow, T.V.; Goodridge, L.D.; Danyluk, M.D. Detection of Salmonella spp. from large volumes of water by modified Moore swabs and tangential flow filtration. Lett. Appl. Microbiol. 2013, 56, 88–94. [CrossRef] 26. Rentería-Guevara, S.A.; Rangel-Peraza, J.G.; Rodríguez-Mata, A.E.; Amábilis-Sosa, L.E.; Sanhouse-García, A.J.; Uriarte-Aceves, P.M. Effect of Agricultural and Urban Infrastructure on River Basin Delineation and Surface Water Availability: Case of the Culiacan River Basin. Hydrology 2019, 6, 58. [CrossRef] 27. Jiménez, M.; Martínez-Urtaza, J.; Rodríguez-Álvarez, M.X.; León-Félix, J.; Chaidez, C. Prevalence and genetic diversity of Salmonella spp. in a river in a tropical environment in Mexico. J. Water Health 2014, 12, 874–884. [CrossRef] [PubMed] 28. EPA. Method 1200: Analytical Protocol for Non-Typhoidal Salmonella in Drinking Water and Surface Water. Available online: https://www.epa.gov/sites/default/files/2015-08/documents/epa817r12004.pdf (accessed on 1 May 2018). 29. Malorny, B.; Hoorfar, J.; Bunge, C.; Helmuth, R. Multicenter validation of the analytical accuracy of Salmonella PCR: Towards an international standard. Appl. Environ. Microbiol. 2003, 69, 290–296. [CrossRef] [PubMed] 30. Andrews, W.H.; Wang, H.; Jacobson, A.; Ge, B.; Zhang, G.; Hammack, T. Chapter 5: Salmonella. In Bacteriological Analytical Manual, 8th ed.; Food and Drug Administration: Silver Spring, MD, USA, 2014. Available online: https://www.fda.gov/food/ laboratory-methods-food/bam-chapter-5-salmonella (accessed on 1 May 2018). 31. Chen, S.; Zhou, Y.; Chen, Y.; Gu, J. fastp: An ultra-fast all-in-one FASTQ preprocessor. Bioinformatics 2018, 34, i884–i890. [CrossRef] 32. Coil, D.; Jospin, G.; Darling, A.E. A5-miseq: An updated pipeline to assemble microbial genomes from Illumina MiSeq data. Bioinformatics 2015, 31, 587–589. [CrossRef] [PubMed] 33. Bankevich, A.; Nurk, S.; Antipov, D.; Gurevich, A.A.; Dvorkin, M.; Kulikov, A.S.; Lesin, V.M.; Nikolenko, S.I.; Pham, S.; Prjibelski, A.D.; et al. SPAdes: A new genome assembly algorithm and its applications to single-cell sequencing. J. Comput. Biol. 2012, 19, 455–477. [CrossRef] 34. Yoshida, C.E.; Kruczkiewicz, P.; Laing, C.R.; Lingohr, E.J.; Gannon, V.P.; Nash, J.H.; Taboada, E.N. The Salmonella In Silico Typing Resource (SISTR): An Open Web-Accessible Tool for Rapidly Typing and Subtyping Draft Salmonella Genome Assemblies. PLoS ONE 2016, 11, e0147101. [CrossRef] 35. Hunter, P.R.; Gaston, M.A. Numerical index of the discriminatory ability of typing systems: An application of Simpson’s index of diversity. J. Clin. Microbiol. 1988, 26, 2465–2466. [CrossRef] 36. Simpson, E.H. Measurement of Diversity. Nature 1949, 163, 688. [CrossRef] 37. Treangen, T.J.; Ondov, B.D.; Koren, S.; Phillippy, A.M. The Harvest suite for rapid core-genome alignment and visualization of thousands of intraspecific microbial genomes. Genome Biol. 2014, 15, 524. [CrossRef] [PubMed] 38. Stamatakis, A. RAxML version 8: A tool for phylogenetic analysis and post-analysis of large phylogenies. Bioinformatics 2014, 30, 1312–1313. [CrossRef] [PubMed] 39. Tavare, S. Some probabilistic and statistical problems in the analysis of DNA sequences. In Some Mathematical Questions in Biology/DNA Sequence Analysis; Miura, R.M., Ed.; Lectures on Mathematics in the Life Sciences; American Mathematical Society: Providence, RI, USA, 1986; Volume 17, pp. 57–86. Microorganisms 2022, 10, 1214 16 of 17 40. Yang, Z. Maximum likelihood phylogenetic estimation from DNA sequences with variable rates over sites: Approximate methods. J. Mol. Evol. 1994, 39, 306–314. [CrossRef] [PubMed] 41. Letunic, I.; Bork, P. Interactive Tree of Life (iTOL) v5: An online tool for phylogenetic tree display and annotation. Nucleic Acids Res. 2021, 49, W293–W296. [CrossRef] 42. Kaas, R.S.; Leekitcharoenphon, P.; Aarestrup, F.M.; Lund, O. Solving the problem of comparing whole bacterial genomes across different sequencing platforms. PLoS ONE 2014, 9, e104984. [CrossRef] 43. Larsen, M.V.; Cosentino, S.; Rasmussen, S.; Friis, C.; Hasman, H.; Marvig, R.L.; Jelsbak, L.; Sicheritz-Pontén, T.; Ussery, D.W.; Aarestrup, F.M.; et al. Multilocus sequence typing of total-genome-sequenced bacteria. J. Clin. Microbiol. 2012, 50, 1355–1361. [CrossRef] 44. Zhou, Z.; Alikhan, N.-F.; Sergeant, M.J.; Luhmann, N.; Vaz, C.; Francisco, A.P.; Carriço, J.A.; Achtman, M. GrapeTree: Visualization of core genomic relationships among 100,000 bacterial pathogens. Genome Res. 2018, 28, 1395–1404. [CrossRef] 45. Zhang, Y.; Bi, P.; Hiller, J.E. Climate variations and Salmonella infection in Australian subtropical and tropical regions. Sci. Total Environ. 2010, 408, 524–530. [CrossRef] 46. Mancera González, O. The cyclone nobody saw: Risks that potentiated a disaster in the Northwest of Mexico. Rev. Estud. Latinoam. Sobre Reducc. Riesgo Desastre 2021, 5, 110–125. [CrossRef] 47. Spector, M.P.; Garcia del Portillo, F.; Bearson, S.M.D.; Mahmud, A.; Magut, M.; Finlay, B.B.; Dougan, G.; Foster, J.W.; Pallen, M.J. The rpoS-dependent starvation-stress response locus stiA encodes a nitrate reductase (narZYWV) required for carbon-starvation- inducible thermotolerance and acid tolerance in Salmonella typhimurium. Microbiology 1999, 145, 3035–3045. [CrossRef] 48. Kingsley, R.A.; Humphries, A.D.; Weening, E.H.; de Zoete, M.R.; Winter, S.; Papaconstantinopoulou, A.; Dougan, G.; Bäumler, A.J. Molecular and phenotypic analysis of the CS54 Island of Salmonella enterica serotype Typhimurium: Identification of intestinal colonization and persistence determinants. Infect. Immun. 2003, 71, 629–640. [CrossRef] [PubMed] 49. DGE. Dirección General de Epidemiología, Anuario de Morbilidad 1984–2020. Available online: https://epidemiologia.salud. gob.mx/anuario/html/index.html (accessed on 20 June 2021). 50. Estrada-García, T.; López-Saucedo, C.; Zamarripa-Ayala, B.; Thompson, M.R.; Gutierrez-Cogco, L.; Mancera-Martinez, A.; Escobar-Gutierrez, A. Prevalence of Escherichia coli and Salmonella spp. in street-vended food of open markets (tianguis) and general hygienic and trading practices in Mexico City. Epidemiol. Infect. 2004, 132, 1181–1184. [CrossRef] [PubMed] 51. Comisión Nacional del Agua (CONAGUA). Reporte del clima en México. Septiembre 2018. Coordinación General del Servicio Meteorológico Nacional Gerencia de Meteorología y Climatología Subgerencia de Pronóstico a Mediano y Largo Plazo. Avail- able online: https://smn.conagua.gob.mx/tools/DATA/Climatolog%C3%ADa/Diagn%C3%B3stico%20Atmosf%C3%A9rico/ Reporte%20del%20Clima%20en%20M%C3%A9xico/RC-Septiembre18.pdf (accessed on 13 September 2021). 52. Instituto Nacional de Estadística y Geografía (INEGI). Características Edafológicas, Fisiográficas, Climáticas e Hidrográficas de México. Available online: http://www.inegi.org.mx/inegi/spc/doc/internet/1geografiademexico/manual_carac_eda_fis_ vs_enero_29_2008.pdf (accessed on 20 June 2021). 53. Rhodes, M.W.; Kator, H. Survival of Escherichia coli and Salmonella spp. in estuarine environments. Appl. Environ. Microbiol. 1988, 54, 2902–2907. [CrossRef] [PubMed] 54. Yun, J.; Greiner, M.; Holler, C.; Messelhausser, U.; Rampp, A.; Klein, G. Association between the ambient temperature and the occurrence of human Salmonella and Campylobacter infections. Sci. Rep. 2016, 6, 28442. [CrossRef] 55. Akil, L.; Ahmad, H.A.; Reddy, R.S. Effects of climate change on Salmonella infections. Foodborne Pathog. Dis. 2014, 11, 974–980. [CrossRef] 56. Cherrie, M.P.C.; Nichols, G.; Iacono, G.L.; Sarran, C.; Hajat, S.; Fleming, L.E. Pathogen seasonality and links with weather in England and Wales: A big data time series analysis. BMC Public Health 2018, 18, 1067. [CrossRef] 57. Winfield, M.D.; Groisman, E.A. Role of nonhost environments in the lifestyles of Salmonella and Escherichia coli. Appl. Environ. Microbiol. 2003, 69, 3687–3694. [CrossRef] 58. Balamurugan, S.; Dugan, M.E. Growth temperature associated protein expression and membrane fatty acid composition profiles of Salmonella enterica serovar Typhimurium. J. Basic Microbiol. 2010, 50, 507–518. [CrossRef] 59. Simental, L.; Martinez-Urtaza, J. Climate patterns governing the presence and permanence of salmonellae in coastal areas of Bahia de Todos Santos, Mexico. Appl. Environ. Microbiol. 2008, 74, 5918–5924. [CrossRef] 60. Jacobsen, C.S.; Bech, T.B. Soil survival of Salmonella and transfer to freshwater and fresh produce. Food Res. Int. 2012, 45, 557–566. [CrossRef] 61. Levantesi, C.; Bonadonna, L.; Briancesco, R.; Grohmann, E.; Toze, S.; Tandoi, V. Salmonella in surface and drinking water: Occurrence and water-mediated transmission. Food Res. Int. 2012, 45, 587–602. [CrossRef] 62. McEgan, R.; Mootian, G.; Goodridge, L.D.; Schaffner, D.W.; Danyluk, M.D. Predicting Salmonella populations from biological, chemical, and physical indicators in Florida surface waters. Appl. Environ. Microbiol. 2013, 79, 4094–4105. [CrossRef] [PubMed] 63. Gu, L.; Chen, Q.; Guo, A.; Liu, W.; Ruan, Y.; Zhang, X.; Nou, X. Differential Effects of Growth Medium Salinity on Biofilm Formation of Two Salmonella enterica Strains. J. Food Prot. 2020, 83, 196–203. [CrossRef] [PubMed] 64. Martínez-Urtaza, J.; Liebana, E. Investigation of clonal distribution and persistence of Salmonella Senftenberg in the marine environment and identification of potential sources of contamination. FEMS Microbiol. Ecol. 2005, 52, 255–263. [CrossRef] [PubMed] Microorganisms 2022, 10, 1214 17 of 17 65. Gorski, L.; Liang, A.S.; Walker, S.; Carychao, D.; Aviles Noriega, A.; Mandrell, R.E.; Cooley, M.B.; Dudley, E.G. Salmonella enterica serovar diversity, distribution, and prevalence in public-access waters from a Central California coastal leafy green-growing region from 2011 to 2016. Appl. Environ. Microbiol. 2022, 88, e01834-21. [CrossRef] [PubMed] 66. CDC. Center for Disease Control and Prevention: Investigation Details. Salmonella Outbreak Linked to Onions. Available online: https://www.cdc.gov/salmonella/oranienburg-09-21/details.html (accessed on 31 March 2022). 67. Heinz, E.; Stubenrauch, C.J.; Grinter, R.; Croft, N.P.; Purcell, A.W.; Strugnell, R.A.; Dougan, G.; Lithgow, T. Conserved Features in the Structure, Mechanism, and Biogenesis of the Inverse Autotransporter Protein Family. Genome Biol. Evol. 2016, 8, 1690–1705. [CrossRef]