A total of 209 environment samples were collected from 2017–2019 from a large dairy farm (herd size = 25,000 animals) in Xinjiang province, China. The environmental samples included were fecal samples (n = 50), manure slurry from a storage tank (n = 36), raw milk (n = 90), water samples from the residential area (n = 9), soil samples (n = 12), and crop field soil (n = 12) based on random sampling technique (Figure 1). A 50 g of manure sample was collected from the animal living area and storage tank from five different sites using a five-point mixed sampling method (Sharp et al., 2012) and stored in sterile zipper bags. Raw milk samples (10 mL) were collected and transferred to sterile falcon tubes according to the guidelines of the National Mastitis Council (Hogan and Smith, 1992). The water samples (50 mL) were collected in sterile water bottles from the residential area by randomly selecting three different water outlets. The blank and crop field soil samples were collected from different sites on farm and fodder growing fields, respectively. All the collected samples were kept at 4°C and transferred to the laboratory within 24 h for further processing.

The 25 g of fecal, manure, and soil samples were first mixed in 225 mL of phosphate-buffered saline (PBS) to solubilize them. After mixing, 1 mL of the liquid was transferred to a 10 mL LB broth tube for bacterial enrichment by incubation at 37°C with continuous mixing at 160 rpm. From each tube, 100 μL of the enrichment culture was sub-cultured on MacConkey (MAC) agar under prior mentioned incubation conditions. However, the water and milk samples were swabbed directly onto MacConkey agar and incubated at 37°C overnight. Based on colony shape and color, large, smooth, and pink colonies were picked and further streaked onto Eosin Methylene Blue (EMB) agar and incubated overnight at 37°C. The appearance of a metallic green sheen with dark center colonies on EMB agar was indicative of E. coli growth (Peng et al., 2022). Further, presumptive identification was done by the VITEK2 system (BioMerieux, France) (Alfinete et al., 2022) and confirmation by 16S rRNA gene amplification by PCR using primers reported previously (Liu et al., 2021). The PCR-amplified product was visualized on 1% agarose gel under the GelDoc XR system (Supplementary Figure S1). The confirmed isolates were preserved in 20% glycerol at −80°C for further analysis.

The 16S rRNA-confirmed E. coli isolates were subjected to the identification of seven VAGs by the previously described method (Hu et al., 2022). The genomic DNA was extracted using a DNA extraction kit (Tiangen Biotech Beijing, Co., Ltd.) following the manufacturer’s guidelines. The virulence genes were identified by PCR amplification of target gene primers mentioned in Supplementary Table S1. The PCR reaction mixture (25 μL) consisted of 12.5 μL PreTaq Mix (Vazyme Biotech, China), 1 μL forward primer, 1 μL reverse primer, 1 μL genomic DNA, and 9.5 μL of deionized water under the following conditions; prior denaturation at 94°C for 5 min followed by 35 cycles of denaturation at 94°C for the 30s, annealing at varying temperatures mentioned in Supplementary Table S1 for 30s, initial extension at 72°C for 30s followed by a final extension at 72°C for 5 min. After amplification, the PCR product was run on 1% agarose gel electrophoresis at 180 V/200 mA followed by ethidium bromide staining for visualization, and images were taken under the GelDoc XR system (Supplementary Figure S2).

The AST was done by broth micro-dilution assay following the EUCAST guidelines.¹ Briefly, the preserved isolates were thawed at room temperature and re-suspended in LH broth by vigorous mixing (120 rpm) at 37°C for 12 h. The loopful enriched broth was streaked on MacConkey agar following the overnight incubation. The bacterial inoculum was prepared by adjusting the cell density at 5 × 10⁵ CFU/mL. The 96-well round bottom plate was used for broth dilution assay and 100 μL of Mueller Hinton (MH) broth was added from the 1st well to the 12th well with a micropipette. Next, 50 μL of prepared bacterial inoculum was added from the 1st to 11th well by keeping the 12th well as a negative control. The antibiotics were selected based on medical and veterinary use which includes trimethoprim-sulfamethoxazole (SXT), ampicillin (AMP), cefotaxime (CTX), tetracycline (TET), ciprofloxacin (CIP), gentamicin (GEN), amikacin (AMK), colistin sulfate (CS), florfenicol (FFC), meropenem (MEM), and tigecycline (TIG) were added from 1st to 10th well by keeping 11th well as a positive control. The reference strain E. coli ATCC 25922 was used as a quality control. The MIC (minimum concentration that inhibits visible growth of bacteria) of fosfomycin was calculated by the agar dilution method, recommended by EUCAST. The MIC of all antibiotics was evaluated by visualizing the growth in the bottom of the plate well as tinny buttons/turbidity. The MIC values were compared with standard EUCAST MIC breakpoints (Supplementary Table S2). The strains showing resistance to at least one antibiotic from ≥3 classes were classified as MDR.

Phenotypically resistant E. coli strains were subjected to the detection of 20 ARGs from eight antibiotic classes (Supplementary Table S3) according to the method described previously (Yu et al., 2020). The bacterial DNA was extracted using a DNA extraction kit (Tiangen Biotech Beijing, Co., Ltd.) and used as a template for PCR amplification of 20 ARGs (listed in Supplementary Table S3). The PCR reaction mixture (25 μL) consisted of 12.5 μL PreTaq Mix (Vazyme Biotech, China), 1 μL of forward and reverse primer each, 1 μL of bacterial DNA, and 9.5 μL of deionized water. The reactions were performed under the following conditions: denaturation at 94°C for 5 min followed by 35 cycles of denaturation at 94°C for 30s, annealing at varying temperatures (see Supplementary Table S3) for 30s, and extension at 72°C for 30s, followed by a final extension at 72°C for 5 min. After amplification, the PCR product was separated on 1% agarose gel at 180 V/200 mA and stained with ethidium bromide for visualization using the GelDoc XR system.

The phylogenetic grouping of E. coli strains was carried out by 2 sets of PCR using primers listed in Supplementary Table S4. The quadruple PCR reaction mixture (25 μL) consists of Premix Taq TM 12.5 μL, 1 μL of each forward and reverse primers (chuA, yjaA, tspE4C2), 2 μL of each arpA forward and reverse primers, 1.5 μL DNA template, and 1 μL dd H₂O. PCR was carried out under the following conditions; pre-denaturation at 94°C for 4 min, 30 cycles of denaturation at 94°C for 5 s followed by annealing at 59°C for 20 s, and extension at 72°C for 5 min. The PCR reaction for group E and C identification consisted of Premix Taq TM 12.5 μL, 0.6 μL of trpBA forward and reverse primers each, 1 μL of each group-specific primer (Supplementary Table S4), 1.5 μL DNA template, and 2.8 μL dd H₂O. In the PCR reaction solution, trpBA primers were added as an internal control. The PCR amplifications conditions were pre-denaturation at 94°C for 4 min, 30 cycles of denaturation at 94°C for 5 s, annealing at 57°C (group E) or 59°C (group C) for 20 s, and final extension at 72°C for 5 min. After the PCR amplification, the PCR product was run on 1% agarose gel and visualized under the GelDoc XR system (Supplementary Figure S3), and the phylogenetic group was identified by comparing the results with Supplementary Table S4.

The prevalence was calculated using the formula described by Thrusfield (2018).

The antimicrobial susceptibility data were analyzed by descriptive statistics using Microsoft Excel. Moreover, the data for various factors such as sampling source and sampling year affecting the prevalence, AMR and virulence rates were analyzed using the Pearson’s Chi-Squared test keeping the level of significance, α = 5% (Zhao et al., 2021; Ma et al., 2022). p-value < 0.05 was considered statistically significant and vice versa. The graphical representation of data was done by GraphPad Prism version 8.2.1 and Microsoft Excel.

A total of 209 samples were collected from different sites of dairy environment including fecal samples (n = 50), manure slurry from the storage tank (n = 36), raw milk (n = 90), water samples (n = 9), soil samples (n = 12), and crop field soil (n = 12) samples. In total, 534 suspected E. coli strains were isolated from 141/209 (67.5%) samples based on colony characteristics. Subsequently, 338 E. coli strains were confirmed by 16S rRNA gene amplification. The isolation rates were comparable over the years, with 30.8% (104/338) in 2017, 34.9% (118/338) in 2018, and 34.3% (116/338) in 2019 (Figure 2A). Overall, most of the E. coli strains were isolated from manure slurry (39.3%, 133/338), followed by fecal samples (34.9%, 118/338), raw milk (24.8%, 84/338), crop field soil (0.59%, 2/338), and least from blank soil (0.29%, 1/338). However, none of the E. coli strains was isolated from water samples (Figure 2B). In 2019, a higher number of E. coli strains were isolated from fecal and milk samples compared to other sampling years while more E. coli strains were isolated from slurry samples in 2018. Moreover, only 1 and 2 strains were isolated in 2017 from blank and crop field soil samples, respectively, while none in 2018 and 2019 (Figure 2B).

The AST of 338 E. coli strains showed that 284/338 (84.0%) were resistant to at least one antibiotic and 54/338 (26.0%) were susceptible strains (Figure 3A). Most of the strains were resistant to trimethoprim/ sulfamethoxazole (62.43%, 211/338), followed by cefotaxime (44.08%, 149/338), ampicillin (33.73%, 114/338), ciprofloxacin (31.36%, 106/338), tetracycline (28.99%, 98/338), and less to florfenicol (7.99, 27/338), gentamicin (4.44%, 15/338), amikacin (1.77%, 6/338), and fosfomycin (1.18%, 4/338). All of the E. coli strains were susceptible to meropenem, tigecycline, and colistin sulfate (Figure 3B). All E. coli strains from 2017–2019 were found 100% susceptible to meropenem, tigecycline, and colistin sulfate. Additionally, the AMR rate of AMP (44.23%), CIP (40.38%), TET (44.23%), and GEN (4.81%) was noted higher with a significant difference (p < 0.05) in 2017 than in other sampling years. Moreover, E. coli strains showed higher AMR to CTX (55.93%), SXT (68.64%), and AMK (3.39%) in 2018 with a significant difference (p < 0.05). However, none of the E. coli strain from 2017 and 2019 exhibited resistance to AMK and FOS (Figure 3C).

Most of the E. coli strains from all samples were resistant to trimethoprim/sulfamethoxazole (SXT) and 100% susceptible to meropenem (MEM), tigecycline (TIG), and colistin sulfate (CS). Moreover, E. coli strains from fecal samples exhibited higher resistance to ampicillin (AMP), ciprofloxacin (CIP), and tetracycline (TET) in 2017 than other sampling years with a significant difference (Figure 4A). A similar trend was observed for E. coli strains isolated from milk and manure slurry (Figures 4B,C). Furthermore, E. coli strains isolated from blank and crop field soil in 2017 exhibited 100% resistant to CTX, CIP, TET, and SXT, while none of the strains isolated from 2018 and 2019 was resistant (Figures 4D,E). In addition, E. coli strains from crop field showed 50% resistance to AMP and florfenicol (FFC).

Among the resistant strains, 126/284 (44.4%) were identified as multi-drug resistant (MDR) and 158/284 (55.6%) were recognized as non-MDR (Figure 5A). Most of the strains showed resistance to 2 antibiotics (63.38%), followed by 3 (50.0%), 1 (36.61%), 4 (30.28%), 5 and 6 (7.75% each), 7 (3.17%), and 8 (0.70%) antibiotics (Figure 5B). Moreover, diverse AMR patterns were recognized such as CTX + AMP, CTX + AMP + SXT, CTX + AMP + SXT + CIP, CTX + TET + SXT + CIP, CTX + AMP + SXT + CIP + GEN, CTX +AMP + SXT + CIP + TET + FFC, AMP + CTX + GEN + TET + SXT + FFC + FOS, and AMP + CTX + CIP + GEN + TET + SXT + FFC. Only 1 strain showed resistance to 8 antibiotics, AMP + CTX + CIP + AMK + GEN + TET + SXT + FFC (Table 1).

The genotypic analysis was done by targeting 20 ARGs among 8 classes of antibiotics (mentioned in Table 1). Eighteen out of 20 ARGs were identified, and the prevalent genotypes included sul2 (67.3%, sulfonamides), bla_TEM (56.3%, beta-lactams), gyrA (73.6%, quinolones), tet(B) (70.4%, tetracycline’s), aph(3)-I (85.7%, aminoglycosides), floR (44.4%, amphenicol), and fosA3 (100%, phosphonic). The percentage distribution of other ARGs identified was as follows: sulfonamides (sul1, 27.9%; sul3, 18.1%), β-lactams (bla_OXA, 25.8%; bla_CTX-M, 22.4%), aminoglycosides (aac(3)-IV, 14.3%; aac(3)-II, 33.3%; aadA, 0.00%; rmtB, 4.76%), quinolones (qnrB, 0.94%; qnrS, 9.43%), polymyxin (pmrB, 0.35%) and tetracycline’s (tet(A), 11.2%; tet(D), 0.00%) (Figure 6).

The correlation between phenotypic resistance and genotypic detection of ARGs was noted variable. For example, no strains showed resistance to colistin sulfate phenotypically but one strain was carrying the ARG upon genotypic analysis. Moreover, the number of strains carrying ARGs was noted higher as compared to phenotypic resistance among sulfonamides, beta-lactams, and aminoglycosides-resistant strains while the inverse was noted among quinolone, tetracycline, and amphenicol-resistant strains. However, the phenotypic and genotypic expression was observed 100% correlated for fosfomycin-resistant strains (Table 2).

Among the VAGs, ompA was most prevalent (86.69%), followed by ibeB (85.0%), traT (84.91%), ompT (73.96%), fyuA (23.1%), iroN (23.1%), and irp2 (21.9%) (Figure 7A). All of the E. coli strains carrying the fyaA gene were also carrying the iroN gene. Moreover, 93.59% of E. coli strains carrying irp2 were also harboring the fyuA gene. VAGs such as ompT, traT, iroN, ibeB, and ompA were detected in E. coli strains from all sources while irp2 and fyuA genes were not observed from manure slurry and fecal samples, respectively. However, both irp2 and fyuA genes (36.9%, 31/84) were identified in strains of milk origin. Collectively, a higher percentage of VAGs was identified in strains of milk origin as compared to feces and slurry (Figure 7B).

The phylogenetic analysis of 338 E. coli strains showed most of the strains belong to the B1 group (75.45%, 255/338), followed by A (18.34%, 62/338), C (2.96%, 10/338), D (1.18%, 4/338), E (1.18%, 4/338), and F (0.30%, 1/338). However, the phylogenetic group for 2 of the strains was not identified. The most prevalent VAGs among the various phylogenetic groups were as follows; B1 (ompA, 87.4%), A (ibeB and ompA, 88.7%), C (traT and ompA, 90.0%), D (traT, ibeB, and ompA, 100%), E (traT, 100%), and F (traT, 100%). Moreover, the percentage distribution of other VAGs among the phylogenetic groups is presented in Table 3.