Forty Holstein cows were enrolled at the start of the dry-off period in June-November of 2019 at the Michigan State University (MSU) Dairy Cattle Teaching and Research Center as described (Vasco et al., 2023). After the last milking, twenty cows (cases) received a single IMM infusion containing 500 mg of ceftiofur hydrochloride (CHCL; SpectramastDC^®; Zoetis Animal Health) along with a non-antibiotic teat sealant with bismuth subnitrate (Orbeseal^®; Zoetis Animal Health) in each teat (n = 4; total of 2 g of CHCL). Cows in the control group (n = 20) received only the IMM teat sealant. All cows had a somatic cell count (SCC) of <150,000 cells/mL at the most recent Dairy Herd Improvement Association test and none received antibiotics in the prior 90 days of lactation. The study protocol was approved by the Institutional Animal Care and Use Committee at MSU (IACUC number ROTO201800166) prior to sampling.

To avoid a parity effect, controls were matched to the ceftiofur-treated cows based on parity as well as monthly milk production. Diet regimens were formulated using Spartan Dairy 3™ software per guidelines outlined in the Nutrient Requirements of Dairy Cattle report (National Research Council, 2001). Based on the dietary information extracted from farm records, different diets were given to the cows in accordance with their production demands across the sampling period (Table 1); the matched treated and control cows were given the same diets within each lactation stage. Near the end of lactation, which corresponded to a day prior to the IMM treatment (day -1), cows received the maintenance diet containing 14% more metabolizable energy and 2.5 times more metabolizable protein (g) than is provided in the dry-off diet (weeks 1 through 5) (Supplementary Table S1). At week 9, animals were given a diet for fresh cows consisting of 64% of dry matter intake when compared to lactating cows, but with transitioning levels of energy and protein that were 15% and 64% higher than during dry-off, respectively.

Animals in all phases received corn silage, soybean meal with 47.5% crude protein, CFE MSU dairy base, and haylage. It was only during the lactation and fresh periods that the ration included corn (ground fine and fed dry) as well as MSU fresh high supplement to increase the energy density and provide essential nutrients such as calcium, magnesium, potassium, and niacin to prevent metabolic disorders that can occur during the transition into lactation. Comparatively, soybean hull pellets, cottonseed, and long bulk brown midrib (BMR) corn silage (CS) were exclusively given to cows in late lactation, while alfalfa hay was only provided to fresh cows. Although grass silage, MSU straw, grass pasture, and ignite supplement (Quality Liquid Feeds^®, WI, Unites States) were given to all animals, they were only provided during the dry period in small quantities. The ignite supplement contains 25% protein, fat, trace minerals, and vitamins A, D, and E (https://onealsfarmandgarden.com/products/tub-qlf-ignite-30).

Fecal samples were collected from all cows at four time points corresponding to the different stages of lactation (Figure 1). These stages included late lactation (day -1), dry-off (weeks 1 and 5), and the periparturient period (week 9). Along with physiological differences, the dietary needs of the animals and feed formulations differ across the three lactation stages. Samples were collected before ceftiofur administration at the end of lactation when the cows received a maintenance diet (day -1) as well as after dry-off (week 1; week 5) when they were given an early dry diet. The final sample was collected during or just prior to calving (week 9) when the cows were given the fresh formulation. One cow could not be sampled at week 9 due to birthing by cesarean section that required antibiotics treatment, leaving 79 samples for analysis from the 9-week sampling point and 159 samples in all. As indicated previously (Vasco et al., 2023), the fecal samples were collected via the rectum using clean obstetric sleeves and transported in a cooler to MSU in sterilized sampling bags for processing.

Each sample was homogenized by hand and aliquots containing 0.25 g of feces and 0.25 g of feces in 750 µL of 190 Proof ethanol were stored for metabolite and DNA extractions, respectively. All fecal aliquots were flash frozen with liquid nitrogen for 1 min and were stored at −80°C until further processing. The 119 samples collected at the day -1, week 1, and week 9 samplings were used for both untargeted metabolomics and metagenomic sequencing, while metagenomic sequencing data were also available for the 40 samples collected at week 5 through our prior study (Vasco et al., 2023). Since we previously observed similar microbiota diversity and composition in the samples collected at weeks 1 and 5, metabolomics was not applied to the week 5 (dry-off) samples.

Metabolite extractions were performed on all 119 fecal samples collected at day -1 (n = 40), week 1 (n = 40), and week 9 (n = 39). Internal standard solutions were prepared for quality control and normalization including: 1) labeled short-chain fatty acids (SCFAs) (10 µM each of [¹³C]sodium formate, [¹³C₂]sodium acetate, [¹³C₃]sodium propionate, and [¹³C₄]sodium butyrate) in 50:50 (v/v) methanol/water; 2)) [¹³C₁₆]palmitic acid (10 µM in 100% isopropanol); 3) phenylalanine-d₇ (10 µM in 50:50 methanol/water); 4) succinic acid-d₄ (10 µM in 50:50 methanol/water); and 5) labeled bile acids (10 µM each of glycocholic acid-d ₄ and glycoursodeoxycholic acid-d ₄ in 50:50 methanol/water). A total of 20 mg of feces was weighed under sterile conditions and 350 µL of ice-cold methanol containing 0.1% butylated hydroxytoluene (BHT) was added. The sample was homogenized and incubated on ice for 10 min. For feces sedimentation, 10 µL of each standard was mixed into the samples, agitated for 30 s, and centrifuged at 10,000 × rpm at 4°C for 10 min. The supernatant was pipette-transferred to a sterile microcentrifuge tube on ice, while ice-cold HPLC-grade isopropanol (200 µL) was added to the pellet, homogenized for 30 s, and centrifuged at 10,000 × rpm and 4°C for 10 min. Finally, the isopropanol supernatant was combined with the initial extract and 100 µL aliquots of the mixed extracts were stored into glass vials inserted in 2-mL amber glass autosampler vials sealed with 9 mm screw septum caps. Metabolite extracts were preserved at −80°C until analyzed.

Analysis of mass-spectra identified twice the amount of nonpolar (n = 11,007) metabolite clusters than polar (n = 5,390) clusters. Likewise, molecular networks aggregating metabolites based on their MS² spectral similarity resulted in 1,122 nonpolar and 658 polar components (Supplementary Figure S1A). Network components connect clusters (nodes) that are structurally related via edges representing a modified cosine score that is calculated based on ions that differ by the mass difference (Supplementary Figure S1B). Only a small fraction of clusters had annotations based on library matches, corresponding to 2.48% of the total metabolites (polar, n = 135; nonpolar, n = 270), of which 68 were found with both polar and nonpolar modes.

In the metabolite clusters with library matches, various classes of metabolites were observed in all samples (Figure 3). These included amino acids, lactones, carboxylic acids, cyclic anhydrides, phospholipids, glycerides, ketones, sugars, fatty acids, nucleosides, chromones, vitamins, and butanoate derivatives. Additionally, several highly abundant metabolites including guanosine, benzofuran-2-one, phthalic anhydride, phosphocholine, monoelaidin, anzacyclotridecan-2-one, cytidine, glycan lacto-N-biose, glycan lacto-N-biose, hexanedioic acid, propanoic acid, octadecenoic acid, and others, were found in all samples. A comprehensive list of the relative abundance of each metabolite and class detected is shown in Supplementary Table S2. Importantly, known metabolites categorized at the class level clustered together based on the stage of lactation but not the ceftiofur treatment status. The samples collected during the fresh period were most similar to each other.

Among the 159 samples, 262 metabolic pathways were identified that were assigned to bacterial taxa representing 797 pathways with different bacterial contributions. Only nine bacterial genera, however, were assigned to the pathways and included: Bifidobacterium (n = 75 pathways), Clostridium (n = 4), Escherichia (n = 25), Methanobrevibacter (n = 18), Olsenella (n = 2), Ruminococcaceae unclassified (n = 26), Sarcina (n = 12), Turicibacter (n = 12), other (n = 262), and unclassified (n = 200). On average, 93% of the reads were classified as unmapped and 6% as unintegrated for microbial metabolic pathways.

At the class level the most abundant microbial pathways were associated with the biosynthesis of amino acids, nucleoside and nucleotides, carbohydrates, and vitamins (Supplementary Figure S2). Aminoacyl-tRNA charging, fermentation and glycolysis pathways were also highly abundant in all samples. Pathways associated with the biosynthesis of fatty acid/lipids, cell structures, aromatic compounds, pentose phosphate and secondary metabolites were also common but were found in lower abundance. Similarly, pathways linked to the degradation of nucleoside/nucleotides as well as carbohydrates, carboxylates, and amine polyamines were commonly found (Supplementary Table S3). The abundance of microbial metabolic pathways varied between phases, with late lactation predominantly displaying elevated levels for most pathways. Nonetheless, no distinct clustering at the class level was evident among the lactation stages or by treatment status (Supplementary Figure S2).

The within-sample diversity was measured by comparing the number of observed features and the Shannon index for the metabolomes and metagenomes among samples collected from the same time point, between two time points, and over the entire sampling period while accounting for repeated measures. For the number of nonpolar metabolites, considerable fluctuations in alpha diversity were detected over time (ANOVA, nonpolar: p = 0.001). The number of polar metabolites, however, displayed a more consistent number of clusters (Observed, ANOVA, p = 0.288) but with varying evenness (Shannon, ANOVA, p = 0.005) (Figures 4A,D). A significantly greater number of nonpolar metabolites was detected during lactation (day -1) relative to the dry-off period (week 1) (t-test: p < 0.0001), though no difference was observed relative to the fresh period (t-test, p = 0.086). Fresh cows exhibited a similar number of metabolites as when they were dry (t-test, polar: p = 0.71, nonpolar: 0.98) but with lower evenness for polar (Shannon, t-test, p = 0.0028) and nonpolar (Shannon, t-test, nonpolar: p = 0.0025) metabolites (Figures 4B,E). Although the metabolite richness was similar between fresh and lactating cows, the Shannon index showed lower diversity in the fresh phase denoting a transition like the one detected for diet composition. No differences were observed between polar and nonpolar metabolites when comparing the ceftiofur-treated and control cattle at any of the time points.

For the metabolic pathway predications, fluctuations in alpha diversity were also observed across samplings (ANOVA, p < 0.0001) (Figures 4C,F). Interestingly, the alpha diversity was more similar between the dry and fresh cows than between fresh and lactating cows despite the similarity in diet. Compared to cows in the dry and fresh periods, cows in late lactation had a significantly greater quantity and diversity of microbial metabolic pathways (t-test, p < 0.0001). Stratifying by treatment status did not result in significant differences in alpha diversity between the groups, except that the number of metabolic pathways at week 5 was significantly higher in control cows as compared to the cows given IMM ceftiofur at dry off (t-test, p < 0.001). Because samples from week 5 were not analyzed with untargeted metabolomics, however, we could not compare between treatment groups at this time point. Nonetheless, no long-term effects in the number of metabolites were observed at week 9 (fresh cows).

Bray-Curtis dissimilarity distances showed significant differences between samplings for polar and nonpolar metabolome composition comprising all metabolite clusters (PERMANOVA, p < 0.001) (Figure 5A). Although the microbial pathways had overlapping composition between dry and fresh cows (PERMANOVA, p > 0.3), the lactating cows showed a significantly higher dispersion in the PCoA compared to the samples from dry and fresh cows (PERMDISP, F = 53.32, p = 1.34e-11) as well as a different average composition (PERMANOVA, F = 63.69, p = 0.001) (Figure 5B). Despite the metabolome composition differences associated with the sampling period, the microbial metabolic pathways were similar in the dry and fresh phases. Furthermore, cows treated with IMM ceftiofur had an identical mean metabolite composition as the controls (PERMANOVA, p > 0.38) even though differences in the composition of the microbial pathways were observed in weeks 5 (PERMANOVA, F = 4.25, p = 0.007) and 9 (PERMANOVA, F = 2.67, p = 0.045).

After comparing the abundance of a total of 16,589 metabolite clusters, 3,753 metabolite components, and 797 microbial-metabolic pathways, only one metabolite cluster was significantly different between controls and cows treated with IMM ceftiofur 1 week after the treatment. This cluster corresponds to the nonpolar-positive metabolite cluster #6574 with a parent mass of 245.07 m/z, and a consensus retention time of 1.08 min (Supplementary Figure S3). This cluster is not identifiable and was not part of a network component, thereby limiting our understanding of its occurrence in the ceftiofur-treated animals. Similarly, RF could not correctly classify the metabolomic composition by treatment group at any time point; the OOB estimate of error rate was >55% when groups were compared based on metabolite or microbial-pathway composition at each time point.

Comparing groups with RF when considering only metabolite clusters with a library ID, which excluded unknown metabolites, the classifier error between treatment groups was reduced. Specifically, the OOB was 30% at day -1, 25% at week 1, and 18% at week 9. The primary identifiable metabolites that influenced the classification with RF varied across the time points. A week after drying off, for instance, inosine and palmitoylcarnitine were higher and lower, respectively, in cows treated with IMM ceftiofur (Mean decrease accuracy (MDA): 10.43, 7.83) than controls (Figure 6A). At week 9, the most important classifier was Leu-Val, followed by (E)-5-(4-methoxy-5-methyl-6-oxopyran-2-yl)-3-methylhex-4-enoic acid, which was higher in the controls relative to the antibiotic-treated cows (MDA: 14.31) (Figure 6B). By contrast, aleuretic acid (MDA: 7.49) was higher in the antibiotic-treated cows at week 9 compared to the controls.

In the microbial metabolic pathway analysis, we identified nuanced distinctions among treatment groups. During the first week, for example, only seven microbial metabolic pathways exhibited disparities between cows treated with ceftiofur and the control group. These differences included an elevated presence of thiamine phosphate pathways by yeasts, as well as reduced levels of biosynthesis of cis-vaccenate and the degradation of punine ribonucleosides, D-fructuronate, and 4-deoxy-L-threo-hex-4enopyranuronate in the ceftiofur-treated cows.

At week 5, 17 microbial pathways were significantly less abundant in the cows treated with ceftiofur than the control cows; these included pathways involved in the biosynthesis of amino acids (i.e., L-ornithine, L-isoleucine, L-lysin, L-threonine, L-methionine), peptidoglycan, glycogen, isoprene, preQ0, chorismate, and coenzyme A, as well as in the degradation of L-arginine (Figure 7A). The unintegrated pathways of Ruminococcaceae bacterium P7 and methylerythritol phosphate pathway I, however, were significantly lower in ceftiofur-treated cows. During week 5, several pathways were primarily identified in control cows but were absent in those given ceftiofur (Figure 7B). These pathways include L-arginine degradation XIII (controls, n = 5), the incomplete reductive TCA cycle (controls, n = 5), and L-ornithine biosynthesis II (controls, n = 6; ceftiofur, n = 1).

At week 9 (fresh phase), eleven microbial pathways explained differences between treatments using RF (Figure 7C). The most important differentially abundant pathways among the ceftiofur-treated cows included those related to peptidoglycan maturation, degradation of D-galacturonate, glycogen, purines, D-galactose, and the biosynthesis of 6-hydroxymethy-dihydropterin diphosphate and L-cysteine. Mirroring the observations from week 5, three microbial pathways were exclusive to controls (n = 4) at week 9, while another pathway was detected in just one ceftiofur-treated cow and seven controls. Curiously, a higher abundance of pathways related to anaerobic energy metabolism, and the biosynthesis of L-methionine and L-alanine were identified in cows treated with ceftiofur compared to the controls.

The OOB error was 1.68% across the samplings, as the metabolomes of two fresh cows were misclassified as lactating cows. Hierarchical clustering of the 50 most important metabolite components identified with RF showed a transitional composition in fresh cows between the lactation stages (Supplementary Figure S4). Only four of these 50 components had clusters with library identification representing long-chain fatty acids, which were more abundant in the dry-off period (week 1), and amino acids that were increased during lactation (Supplementary Figure S5).

Differential abundance tests also identified 9,850 features that differed between samplings, corresponding to 46.59% of the total features (metabolites and microbial pathways) (Supplementary Figure S6). Dry and fresh cows showed a lower number of different features (53.88% of the total) than lactating and dry and cows (64.64%) or lactating and fresh cows (66.9%). In particular, the dry and fresh cows only differed in a few microbial pathways, whereas approximately one-third of the microbial metabolic pathways had varying abundances during lactation as compared to the other stages.

Among the top 25 most important metabolites classified by RF with library matches (Figure 8), cows in late lactation had a greater abundance of oxoheptadecanedioic acid, tri(propylene glycol) butyl ether, apigenin tetramethoxyflavone, phosphocholine, and phenylalanine compared to dry and fresh cows. Conversely, dry cows had higher concentrations of coenzyme Q10, myristoyl ethanolamide, riboflavin, harmol, oxobutanoic acid, guanidineacetic acid, citrulline, and phosphoethanolamide. Fresh cows demonstrated increased abundance of glutamate, histidine, and lysine, alongside histamine, which was elevated in both fresh and lactating cows compared to dry cows. Similarly, an evaluation of the top 25 microbial metabolic pathways showed that multiple pathways were significantly higher in lactating cows compared to those cows in the other stages (Figure 9). Pathways involved in the biosynthesis of nucleotides, amino acids, cell wall and glycolysis are among the most abundant. Intriguingly, during dry-off (week 5), higher levels of 2-oxobutanoate degradation I pathways were detected.

Notably, the top 8 most important pathways that predicted the sampling period through RF corresponded to three categories: 1) cell division, 2) amino acid biosynthesis, and 3) carbohydrate biosynthesis. These pathways were significantly higher during lactation, particularly those related to cell division that were mostly absent during the dry-off period. Cell-division pathways included inosine-5′-phosphate biosynthesis III, pyrimidine deoxyribonucleotides de novo biosynthesis IV, UDP-N-acetyl-D-glucosamine biosynthesis I, and O-antigen building blocks biosynthesis (E. coli), which were assigned to Bifidobacterium, Turicibacter, Olsenella, and Escherichia coli (Supplementary Figure S7). The three main pathways related to amino acid biosynthesis involved the superpathway of L-lysine, L-lysine biosynthesis VI, L-valine biosynthesis, and L-threonine and L-methionine biosynthesis II (Supplementary Figure S8). During lactation, these amino acid biosynthesis pathways were mainly assigned to Bifidobacterium spp. No taxa could be assigned to these pathways in the samples collected during the dry and fresh periods. Finally, carbohydrate biosynthesis pathways included glycogen biosynthesis I (from ADP-D-Glucose), sucrose biosynthesis II, and horismite biosynthesis I (Supplementary Figure S9). Although Sarcina was the main taxa assigned to glycogen and sucrose biosynthesis during the dry-off and fresh periods, Bifidobacterium was mainly associated with carbohydrate biosynthesis during lactation.

Positive Spearman’s correlations among metabolites, microbial pathways, microbial species, viruses, and antimicrobial resistance genes were analyzed. Potential functional relationships between metabolites and microbial species included uncultured Firmicutes and Bacteroidetes species with triacylglycerol (TAG); Clostridium with N,N,N-trimethyllysine, N,N-dimethyldodecylamine N-oxide and myristamidopropyl betaine; Proteobacteria and Blautia with histamine; monolinolenin, beauvericin [+NH4+], and neobavaisoflavone, with Campylobacter, Ruminococcaceae uncultured bacteria (GGB3236), and bacteriophages from enterobacteria (Vectrevirus Vec3); 2′-deoxyadenosine with Gammaproteobacterial; 1,2-Diheptadecanoyl-sn-glycero-3-phosphocholine with uncultured Bacteroidetes. Fruchterman Reingold networks show other correlations between the microbiome and polar metabolites (Supplementary Figure S10) as well as nonpolar metabolites (Supplementary Figure S11).

To better explore the functional associations of differentially abundant features, hierarchical clustering was performed with only those known metabolic clusters that differed significantly between samplings using LefSe, MaAsLin2 and ANCOM-BC. A tree was constructed based on a distance matrix highlighting the fold-change of metabolites and metagenomic features that differed significantly (n = 684). One relevant group showed concomitant higher levels of Actinobacteria, Proteobacteria, and histamine during lactation (Figure 10). This finding suggests a role of these taxa in the production of pro-inflammatory compounds such as histamine, which was also observed in the correlation networks.

Fresh cows showed higher levels of urate and nonpolar plant-derived compounds, which was likely due to a diet rich in alfalfa hay that was provided to the animals during this phase. Other clusters also showed patterns with lower contrast across the stages of lactation, with the most relevant HC groups displayed in Supplementary Figure S12. For instance, higher quantities of amino acids and dipeptides were detected with the polar mode on samples from lactation, which was related to a higher abundance of bacterial amino acid synthesis pathways mentioned priorly (Supplementary Figure S8). Not surprisingly, higher levels of androstane were also detected in fresh cows, which are expected to have peak levels of estrogens at this time (Supplementary Figure S13).