key: cord-0970328-40tgpl01 authors: Zha, Hua; Li, Qian; Chang, Kevin; Xia, Jiafeng; Li, Shengjie; Tang, Ruiqi; Li, Lanjuan title: Characterising the Intestinal Bacterial and Fungal Microbiome Associated With Different Cytokine Profiles in Two Bifidobacterium strains Pre-Treated Rats With D-Galactosamine-Induced Liver Injury date: 2022-03-24 journal: Front Immunol DOI: 10.3389/fimmu.2022.791152 sha: 9c38104bfd5d2ddeca9aa72abc827b5c469b48ac doc_id: 970328 cord_uid: 40tgpl01 Multiple probiotics have protective effects against different types of liver injury. Different intestinal microbes could be beneficial to the protective effects of the probiotics on the treated cohorts in different aspects. The current study was designed to determine the intestinal bacterial and fungal microbiome associated with different cytokine profiles in the Bifidobacterium pseudocatenulatum LI09 and Bifidobacterium catenulatum LI10 pretreated rats with D-galactosamine-induced liver injury. In this study, partition around medoids clustering analysis determined two distinct cytokine profiles (i.e., CP1 and CP2) comprising the same 11 cytokines but with different levels among the LI09, LI10, positive control (PC), and negative control (NC) cohorts. All rats in PC and NC cohorts were determined with CP1 and CP2, respectively, while the rats with CP1 in LI09 and LI10 cohorts had more severe liver injury than those with CP2, suggesting that CP2 represented better immune status and was the “better cytokine profile” in this study. PERMANOVA analyses showed that the compositions of both bacterial and fungal microbiome were different in the LI10 cohorts with different cytokine profiles, while the same compositions were similar between LI09 cohorts with different cytokine profiles. The phylotype abundances of both bacteria and fungi were different in the rats with different cytokine profiles in LI09 or LI10 cohorts according to similarity percentage (SIMPER) analyses results. At the composition level, multiple microbes were associated with different cytokine profiles in LI09 or LI10 cohorts, among which Flavonifractor and Penicillium were the bacterium and fungus most associated with LI09 cohort with CP2, while Parabacteroides and Aspergillus were the bacterium and fungus most associated with LI10 cohort with CP2. These microbes were determined to influence the cytokine profiles of the corresponding cohorts. At the structure level, Corynebacterium and Cephalotrichiella were determined as the two most powerful gatekeepers in the microbiome networks of LI09 cohort CP2, while Pseudoflavonifractor was the most powerful gatekeeper in LI10 cohort with CP2. These identified intestinal microbes were likely to be beneficial to the effect of probiotic Bifidobacterium on the immunity improvement of the treated cohorts, and they could be potential microbial biomarkers assisting with the evaluation of immune status of probiotics-treated cohorts. Liver injury has caused great illness in human beings (1) . It could be induced by multiple factors, e.g., drug, virus, alcohol, food additives, and dietary supplements (2) (3) (4) (5) . A variety of products and materials have been determined to have protective effects against different types of liver injury (6) (7) (8) . The protective effects of probiotics against liver injury have been widely reported. For example, Lactobacillus plantarum C88 was capable of preventing ethanol-induced mice liver injury (9) . Bacillus spores could protect rats from acetaminophen-induced acute liver injury (10) . L. plantarum C88 was found to protect mice from aflatoxin B1-induced liver injury (11) . Different intestinal microbes were likely to work in concert with probiotics to promote health. Bifidobacterium pseudocatenulatum LI09 and Bifidobacterium catenulatum LI10 were found to alleviate the severity of D-galactosamine (D-GalN)-induced rat liver injury (12) . However, the intestinal microbes that can enhance the protective effects of the two Bifidobacterium on the improvement of cytokine profiles have not been determined. The current study was designed to (1) characterize the intestinal bacterial and fungal microbiome associated with different cytokine profiles of LI09 and LI10 pretreated rats with D-GalN-induced liver injury and (2) investigate the microbes with the biomarker potentials to assist with the evaluation of better immune status in the probiotics-treated cohorts. The animal experiments were performed as previously described (12) , with a few modifications. Briefly, B. pseudocatenulatum LI09 and B. catenulatum LI10 were streaked on the trypticasephytone-yeast agar and revived anaerobically at 37°C. The two bacterial strains were then prepared in physiological saline (PBS) at a final concentration of 3 × 10 9 CFU/ml. Sprague-Dawley male pathogen-free rats weighting 250-350 g were fed with a standard laboratory rat chow and raised under the 12:12 light/ dark cycle at 22°C for 7 days to adapt to the environment. The 122 rats were randomly allocated into four cohorts, i.e., LI09 cohort (n = 40), LI10 cohort (n = 40), positive control (PC) cohort (n = 22), and negative control (NC) cohort (n = 20). The rats in LI09 and LI10 cohorts were orally administrated with a 1ml aliquot of LI09 or LI10 (3 × 10 9 CFU) for a week, while the rats in PC and NC cohorts were orally administrated with 1-ml aliquot of PBS for the same period. Afterwards, an intraperitoneal injection of D-GalN was given to each of the rats in LI09, LI10, and PC cohorts at a dose of 700 mg/kg body weight. Twenty-four hours after the induction of liver injury, all the living rats were anesthetized by an intraperitoneal injection of 10 mg/kg xylazine and 80 mg/kg ketamine, before being subjected to laparotomy for collection of the blood, liver, and cecal content. The study was approved by Animal Care and Use Committee of the First Affiliated Hospital, Zhejiang University School of Medicine. Serum was extracted from blood samples by centrifugation and stored at −80°C. Concentrations of liver function variables in serum, i.e., gamma glutamyl transferase (GGT), total bilirubin (TB), total bile acid (TBA), albumin (ALB), aspartate aminotransferase (AST), alanine aminotransferase (ALT), and alkaline phosphatase (ALP), were measured by an automatic biochemical analyzer (Roche Diagnostics, Ottweiler, Germany) according to the manufacturer's instructions. The concentrations of 23 cytokines in the serum samples were measured using a Bio-Plex Pro ™ Rat Cytokine 23-Plex Assay kit (Bio-Rad Ltd., Hercules, CA, USA) as per the protocol of the manufacturer. These cytokines included macrophage inflammatory protein (MIP)-1a, MIP-3a, macrophage colonystimulating factor (M-CSF), granulocyte colony-stimulating factor (G-CSF), granulocyte-macrophage colony-stimulating factor (GM-CSF), interferon gamma (IFN-g), tumor necrosis factor alpha (TNF-a), interleukin (IL)-1a, IL-1b, IL-2, IL-4, IL-5, IL-6, IL-7, IL-10, IL-12p70, IL-13, IL-17A, IL-18, monocyte chemoattractant protein-1 (MCP-1) , growth-related oncogene (GRO)/keratinocyte chemoattractant (KC), regulated upon activation, normal T cell expressed, and secreted (RANTES), and vascular endothelial growth factor (VEGF). The liver tissue from the left liver lobe of each rat was dissected and fixed in 10% formalin solution, before being dehydrated and processed in paraffin using standard histological methods. The liver samples were stained and mounted on microscope slides. The liver injury severity was evaluated by a professional pathologist based on the Ishak scoring system (13) . primers (i.e., 341F/785R) and fungal primers (i.e., ITS3F/ITS4R) (14, 15) . The PCR products were purified by using a DNA Clean and Concentrator Kit (Zymo Research, Irvine, CA, USA), and their concentrations were measured by using a Qubit ™ dsDNA HS Assay Kit (Thermo Fisher Scientific Inc., Waltham, MA, USA). The purified PCR products were submitted for sequencing on Illumina NovaSeq 6000 platform (Illumina Inc. USA). The sequencing data were processed with standard bioinformatic procedures, e.g., merge, quality filter, singleton removal, and chimera check. The sequences were clustered into groups of amplicon sequence variants (ASVs). Bacterial phylotypes were classified using QIIME2 against the Silva 138 database, while fungal phylotypes were classified against the UNITE fungal database. One rat in the LI10 cohort was not recruited in the current study, as it did not provide enough sequencing reads. All the other rats were recruited for the subsequent analyses. One-way ANOVA tests and Mann-Whitney tests were used to determine the differences in 23 cytokines in LI09, LI10, PC, and NC cohorts. The cytokines with significant differences among the four cohorts were selected and transformed in log 2 (X+1) for the cytokine profile analyses. Partition around medoids (PAM) clustering analysis was performed to cluster the cytokine profiles of all the four cohorts, after an average silhouette analysis being conducted to determine the optimal number of clusters. The four cohorts were determined with two distinct cytokine profiles, i.e., CP1 and CP2. The rats in LI09 and LI10 cohorts with two different cytokine profiles were defined as CP1_LI09, CP2_LI09, CP1_LI10, and CP2_LI10 cohorts. T-tests and Mann-Whitney tests were used to compare the liver function variables in CP1_LI09 and CP2_LI09 cohorts. The same tests were performed to compare the liver function variables in CP1_LI10 and CP2_LI10 cohorts. Mann-Whitney tests were carried out to compare the Ishak scores of CP1_LI09 and CP2_LI09 cohorts and those of CP1_LI10 and CP2_LI10 cohorts. The alpha diversity indices (i.e., observed species and Shannon and Pielou indices) of the bacterial and fungal microbiome in the LI09 cohorts with different cytokine profiles, and their control cohorts (i.e., PC and NC cohorts), were all calculated. One-way ANOVA was used to compare the alpha diversity indices of the four cohorts, and t-tests were performed for the comparisons of LI09 and control cohorts. The same strategy was carried out for the same microbiome composition comparisons of LI10 cohorts with different cytokine profiles and their control cohorts. Permutational analysis of variance (PERMANOVA) was conducted in R 4.1.0 to compare the CP1_LI09 and CP2_LI09 cohorts for their bacterial and fungal microbiome compositions and compare them with their control cohorts (i.e., PC and NC cohorts). The same strategy was performed for the same microbiome composition comparisons of CP1_LI10 and CP2_LI10 cohorts and their control cohorts. Similarity percentage (SIMPER) analysis was performed to determine the bacterial and fungal microbiome similarities within the LI09 and LI10 cohorts with different cytokine profiles. The same analysis was used to determine the bacterial and fungal microbiome dissimilarities between CP1_LI09 and CP2_LI09 cohorts. The same strategy was used for the comparisons of the bacterial and fungal microbiome dissimilarities between CP1_LI10 and CP2_LI10 cohorts. Linear discriminant analysis (LDA) effect size (LEfSe) was carried out to compare the bacterial and fungi microbiome of CP1_LI09 and CP2_LI09 cohorts, respectively, to determine the bacteria and fungi associated with each of the two cohorts. The same analysis was performed to determine the bacteria and fungi associated with CP1_LI10 and CP2_LI10 cohorts. Spearman test was used to determine the individual correlations between the cytokines in the cytokine profiles and the microbes associated with CP1_LI09, CP2_LI09, CP1_LI10, and CP2_LI10 cohorts. Distance-based redundancy analysis (db-RDA) was performed to determine the effect of bacteria and fungi associated with CP1_LI09, CP2_LI09, CP1_LI10, and CP2_LI10 cohorts on the corresponding cytokine profiles. Spearman test was carried out to investigate the correlations between the bacteria and fungi associated with CP1_LI09 cohort. The same strategy was performed for the bacteria and fungi correlations in CP2_LI09, CP1_LI10, and CP2_LI10 cohorts. The correlations between the bacteria and fungi in the intestinal microbiome networks of LI09 and LI10 cohorts with different cytokine profiles were determined by co-occurrence network inference (CoNet) analysis. Five metrics, i.e., Bray-Curtis, Spearman, Pearson, mutual information, and Kullback-Leibler dissimilarities, were used to calculate the ensemble inference in CP1_LI09, CP2_LI09, CP1_LI10, and CP2_LI10 cohorts. The top 10 microbes (i.e., bacteria and fungi) with the largest numbers of correlations in the microbiome networks of the four cohorts were determined. Network fragmentation calculations and generation of a null distribution were carried out in R as previously described (16) to explore the gatekeeper(s) in the microbiome networks of LI09 and LI10 cohorts with different cytokine profiles. Statistical significance was defined as the number of times a fragmentation score over that resulting from the removal of the bacteria or fungi observed within the null distribution. Eleven serum cytokines were determined with significant differences among the LI09, LI10, PC, and NC cohorts, i.e., IL-1a, IL-2, IL-4, IL-5, IL-6, IL-12p70, IL-17A, M-CSF, MCP-1, MIP-3a, and RANTES (all p < 0.03), and they were used for the subsequent clustering analysis of cytokine profiles. Higher levels of IL-1a and M-CSF were determined in both LI09 and LI10 cohorts than in PC cohort, and IL-2 and IL-17A were greater in the LI10 cohort than in the PC cohort (Supplementary Figure S1) , suggesting IL-1a, M-CSF, IL-2, and IL-17A were enhanced by LI09 and/or LI10. In contrast, MCP-1, IL-5, and MIP-3a were determined with lower levels in both LI09 and LI10 cohorts than in the PC cohort (Supplementary Figure S1 ), suggesting that they were suppressed by LI09 and LI10. Two was determined as the most optimal number for clustering the cytokine profiles ( Figure 1A) , and two distinct cytokine profiles (i.e., CP1 and CP2) were determined in the LI09, LI10, PC, and NC cohorts ( Figure 1B ). All rats in the PC cohort were determined with CP1, and all rats in the NC cohort had CP2, suggesting that CP2 represented better immune status and was the "better cytokine profile." Twenty-one rats in the LI09 cohort and 18 rats in the LI10 cohort were determined with CP1, while 19 rats in the LI09 cohort and 21 rats in the LI10 cohort had CP2 ( Figure 1B ). Six out of the seven measured liver function variables, i.e., ALT, AST, ALP, TBA, TB, and GGT, were greater in the CP1_LI09 cohort than in the CP2_LI09 cohort ( Figure 2 ). By contrast, ALB was at similar level between the two cohorts. Similarly, ALT, AST, ALP, TBA, TB, and GGT were determined at higher levels in the CP1_LI10 cohort than in the CP2_LI10 cohort (Figure 3) , while ALB was lower in the CP1_LI10 cohort than in the CP2_LI10 cohort ( Figure 3 ). Ishak scoring system was used to help evaluate the liver histopathology, and in this study, a higher Ishak score represented greater liver injury severity. The Ishak score was greater in the CP1_LI09 cohort than in the CP2_LI09 cohort ( Figure 4A ), and the same score was greater in the CP1_LI10 cohort than in the CP2_LI10 cohort ( Figure 4B ). The rats with CP1 in LI09 and LI10 cohorts had more severe liver injury than those with CP2, further suggesting that CP2 was the "better cytokine profile" compared with CP1. Firmicutes, Bacteroidota, and Verrucomicrobiota were determined as the three most abundant bacterial phyla in the LI09 and LI10 cohorts with different cytokine profiles. At the family level, Lachnospiraceae and Bacteroidaceae were determined with the largest abundances in the bacterial microbiome of all the LI09 and LI10 cohorts with different cytokine profiles ( Figure 5) . Akkermansiaceae was the third most abundant bacterial family in CP1_LI09, CP1_LI10, and CP2_LI09 cohorts, while Tannerellaceae was determined with the third most abundance in CP2_LI10 cohort ( Figure 5 ). Shannon and Pielou indices were both similar in bacterial microbiome of the CP1_LI09, CP2_LI09, PC, and NC cohorts (both p > 0.09), while a significant difference was determined in observed species of the four cohorts (p < 0.001). CP1_LI09 and CP2_LI09 cohorts were found with similar observed species, while they were both greater than in the NC cohort (Supplementary Figure S2A) . Likewise, Shannon and Pielou indices were both similar between CP1_LI10, CP2_LI10, PC, and NC cohorts (both p > 0.05), while observed species was different between the four cohorts (p < 0.001). No difference was found in the observed species of CP1_LI10 and CP2_LI10 cohorts, but they were both greater than in the PC and NC cohorts (Supplementary Figure S2B) . PERMANOVA analysis suggested that the bacterial composition was similar between CP1_LI09 and CP2_LI09 cohorts (R 2 = 0.028, p > 0.26), and they were both different from the PC and NC cohorts (p < 0.002) The same analysis determined a significant difference in the bacterial composition between CP1_LI10 and CP2_LI10 (R 2 = 0.063, p< 0.001), and they were both different from the PC and NC cohorts (p < 0.001). SIMPER analysis determined that the similarity of bacterial phylotype abundances within CP1_LI09 was higher than that of CP2_LI09, i.e., 51% versus 46%. The same analysis determined a dissimilarity of 52% between CP1_LI09 and CP2_LI09 cohorts. Likewise, SIMPER analysis determined that the similarity of bacterial phylotype abundances was greater within CP1_LI10 (i.e., 51%) than within CP2_LI10 (i.e., 44%). The same analysis determined a dissimilarity of 55% between CP1_LI09 and CP2_LI09 cohorts. LEfSe analysis determined five bacteria associated with CP1_LI09 cohort and one bacterium (i.e., Flavonifractor) associated with CP2_LI09 cohort (Figure 6A ), among which Lachnospiraceae_UCG_006 was most associated with CP1_LI09 cohort. The same analysis revealed that 30 bacteria were associated with CP1_LI10 and CP2_LI10 cohorts ( Figure 6B ), among which Lachnospiraceae_NK4A136_group and Parabacteroides were most associated with CP1_LI10 and CP2_LI10 cohorts, respectively. Ascomycota and Basidiomycota were determined as the two most abundant fungal phyla in all the LI09 and LI10 cohorts with different cytokine profiles. At the family level, Aspergillaceae and Trichocomaceae were determined with most abundances in the LI09 and LI10 cohorts with different cytokine profiles (Figure 7) . Debaryomycetaceae, Hypocreales_Incertae_sedis, Trichosphaeriaceae, and Mrakiaceae were determined as the third largest abundances in the mycobiome of CP1_LI09, CP1_LI10 and CP2_LI09 and CP2_LI10 cohorts, respectively (Figure 7) . Shannon and Pielou indices were both similar in the mycobiome between the CP1_LI09, CP2_LI09, PC, and NC cohorts (both p > 0.68), while observed fungal species was different in the four cohorts (p < 0.009). No difference was determined in the fungal observed species of CP1_LI09 and CP2_LI09, but they were both less than that in the NC cohort (Supplementary Figure S3) . Likewise, the three alpha diversity indices were all similar in the mycobiome of CP1_LI10, CP2_LI10, PC, and NC cohorts (all p > 0.53). PERMANOVA analysis revealed that the mycobiome composition was similar between CP1_LI09 and CP2_LI09 cohorts (R 2 = 0.03, p > 0.22), but they were both different from PC and NC cohorts (p < 0.004). The same analysis determined a significant difference in the mycobiome composition between CP1_LI10 and CP2_LI10 cohorts (R 2 = 0.049, p < 0.007), and they were both different from PC and NC cohorts (p < 0.025). SIMPER analysis determined that the similarity of fungal phylotype abundances within the CP1_LI09 cohort (i.e., 29%) was lower than that within the CP2_LI09 cohort (i.e., 34%). The same analysis determined a dissimilarity of 69% between CP1_LI09 and CP2_LI09 cohorts. Likewise, SIMPER analysis determined that the similarity of fungal phylotype abundances was lower within the CP1_LI10 cohort than that within the CP2_LI10 cohort, i.e., 34% versus 40%. The same analysis determined a dissimilarity of 65% between CP1_LI09 and CP2_LI09 cohorts. LEfSe analysis determined seven fungi associated with CP1_LI09 and CP2_LI09 cohorts ( Figure 8A ), among which Meyerozyma and Penicillium were most associated with CP1_LI09 and CP2_LI09 cohorts, respectively. The same analysis determined that 12 fungi were associated with CP1_LI10 and CP2_LI10 cohorts ( Figure 8B ), among which Talaromyces and Aspergillus were associated with CP1_LI10 and CP2_LI10 cohorts, respectively. Multiple correlations were determined between the cytokines in cytokine profile and the microbes associated with each of the CP1_LI09, CP1_LI10, and CP2_LI10 cohorts, except the CP2_LI09 cohort (Figure 9 ). ASF356 and Meyerozyma were determined with more correlations in the CP1_LI09 cohort ( Figure 9A) . Lachnospiraceae_NK4A136_group, Meyerozyma, Stemphylium, and Talaromyces were determined with more correlations with the cytokines in the CP1_LI10 cohort, and Talaromyces seemed to have an opposite effect on IL-a and M_CSF when comparing with Meyerozyma and Stemphylium ( Figure 9C) . GCA_900066575 was negatively correlated with most cytokines in the cytokine profile of CP2_LI10 cohort ( Figure 9D) . Some of the bacteria and fungi closely associated with CP1_LI09, CP2_LI09, CP1_LI10, and CP2_LI10 cohorts were determined to influence the cytokine profiles in the corresponding cohorts ( Figure 10 ). Different correlations were determined between the bacteria and fungi associated with CP1_LI09, CP2_LI09, CP1_LI10, or CP2_LI10 cohorts. Three correlations were determined in the bacterial and fungi associated with the CP1_LI09 cohort, i.e., a positive correlation between Prevotellaceae_UCG_001 and Meyerozyma, negative correlations between Lachnospiraceae_ UCG_006 and Alternaria and between ASF356 and Issatchenkia. By contrast, no correlation was determined between the bacteria and fungi associated with the CP2_LI09 cohort. All the nine correlations were determined to be positive between six bacterial and six fungi associated with the CP1_LI10 cohort ( Figure 11) . Defluviitaleaceae_UCG_011, Candidatus_Christensenellaceae, and Pygmaiobacter were positively correlated with Oidiodendron in the CP2_LI10 cohort. CoNet results revealed 10 microbes with most correlations in the intestinal microbiome networks of LI09 and LI10 cohorts with different cytokine profiles (Tables 1, 2) . Seven out of the top 10 microbes with most correlations in the CP2_LI09 cohort were not determined in the top 10 microbes in the CP1_LI09 cohort ( Table 1) . Likewise, the top eight microbes with most correlations in the CP1_LI10 cohort were all not determined in Table 2) . Eubacterium and Clostridium were both determined with many correlations in the microbiome networks of both CP2_LI09 and CP2_LI10 cohorts but not in the probioticstreated cohorts with CP1. Among these bacteria and fungi, Hungatella, Papillibacter, and Myceliophthora were determined as gatekeepers in the network of the CP1_LI09 cohort (fragmentation analyses, all p < 0.05). In the CP2_LI09 cohort, Corynebacterium, Eubacterium, Papillibacter, and Cephalotrichiella were identified as the microbiome network gatekeepers (fragmentation analyses, all p < 0.05). By contrast, Udeniomyces and Pseudoflavonifractor were determined as the only gatekeepers in the microbiome networks of CP1_LI10 and CP2_LI10 cohorts, respectively (fragmentation analyses, both p < 0.03). Fragmentation analysis showed that the fragmentation level of the microbiome network of CP1_LI09 cohort was lower than that of the CP2_LI09 cohort, i.e., 0.497 versus 0.575. Similarly, the fragmentation level was lower in the CP1_LI10 cohort (i.e., 0.416) than that in the CP2_LI10 cohort (i.e., 0.527). Multiple probiotics have been used to alleviate different types of liver injury (17) (18) (19) . Liver injury severity, gut microbiome alterations, changes in liver function, and cytokine variables have been used to evaluate the effects of probiotics on the treated cohort (20, 21) . Different intestinal microbes could have different abilities of improving these different aspects of the probiotics-treated cohorts. The present study demonstrated the intestinal microbiome associated with distinct cytokine profiles of the LI09 and LI10 pre-treated rats with liver injury and explored the microbes with the biomarker potentials to assist with the evaluation of better immune status in the two probiotics-treated cohorts. Cytokine profiles have been investigated in the cohorts treated with probiotics (22, 23) . In this study, the two cytokine profiles (i.e., CP1 and CP2) were determined by PAM clustering analysis based on the overall pattern of the same 11 cytokines, which were determined with significant different levels among LI09, LI10, PC, and NC cohorts. All rats in PC and NC cohorts were determined with CP1and CP2, respectively, while the rats with CP1 in LI09 and LI10 cohorts were determined with more severe liver injury than those with CP2, suggesting that CP2 represented better immune status and was the "better immune profile" compared with CP1. PAM analysis has been widely used in the microbiome studies to cluster the microbiome in the same or different cohorts (24, 25) but was seldom used to cluster the immune variables. The cytokines in the cytokine profile in this study, i.e., IL-1a, IL-2, IL-4, IL-5, IL-6, IL-12p70, IL-17A, M-CSF, MCP-1, MIP-3a, and RANTES, were all important to the immune system of mammals (26, 27) . Their alterations were associated with the cohorts with different types of liver injury. T-cell activation could be promoted by IL-1a in mice with carbon tetrachloride-induced liver injury (28) . IL-2 and IL-17A were independent risk factors that could lead to liver injury in coronavirus disease 2019 (COVID-19) patients at admission (29) . Elevated level of plasma IL-4 was determined in the mice with dicloxacillininduced liver injury (30) . Serum IL-5 was increased in the mice with Schistosoma mansoni-induced granulomatous liver injury (31) . The increased level of serum IL-6 was found in the rats with thioacetamide-induced liver injury (32) . DEREG mice with liver fibrosis were found with higher level of IL-12p70, MCP-1, and RANTES (33) . M-CSF and MIP-3a were increased in the rats with D-GalN-induced acute liver injury (12) . The majority of the measured liver function variables were determined at higher levels in the LI09 and LI10 cohorts with CP1 (i.e., CP1_LI09 and CP1_LI10 cohorts) than the corresponding cohorts with CP2 (i.e., CP2_LI09 and CP2_LI10 cohorts). As previous study has determined that the majority of A B FIGURE 4 | Comparisons of Ishak scores in the (A) LI09 cohorts with different cytokine profiles (i.e., CP1_LI09 and CP2_LI09 cohorts) and control cohorts and (B) LI10 cohorts with different cytokine profiles (i.e., CP1_LI10 and CP2_LI10 cohorts) and control cohorts. PC represents positive control; NC represents negative control. * represented 0.01 < P < 0.05; ** represented P < 0.01. liver function variables were lower in the NC cohort than in the PC cohort (12) , it suggested that the cohorts with CP2 had better liver function compared with those with CP1. Similarly, Ishak score was greater in CP1_LI09 and CP1_LI10 cohorts than the corresponding cohorts with CP2. As lower Ishak score represented lower liver injury severity (34, 35) , it indicated that the cohorts with CP2 profile were with lower liver injury severity. Bacterial microbiome alterations in probiotics-treated cohorts have been well studied (36, 37) . In the current study, PERMANOVA results showed that the compositions of both bacterial and fungal mycobiome were similar between the CP1_L09 and CP2_LI09 cohorts, while the same microbial compositions were different between the CP1_LI10 and CP2_LI10 cohorts. These suggest that the better cytokine profile in the LI10 cohort were associated with the altered bacterial and fungal microbiome compositions. SIMPER analysis has been widely used in the microbiome studies for different objectives (38, 39) . SIMPER analyses revealed that the phylotype abundances of bacterial and fungal microbiome were both different between CP1_LI09 and CP2_LI09 cohorts, and between CP1_LI10 and CP2_LI10 cohorts, suggesting that the different cytokine profiles in LI09 and LI10 cohorts were associated with the altered phylotype abundances of the intestinal microbiome. LEfSe has been carried out in multiple microbiome studies to determine the phylotypes associated with the different cohorts (40, 41) . LEfSe analysis determined that multiple bacteria were associated with LI09 cohorts with different cytokine profiles, among which Lachnospiraceae_UCG_006 and Flavonifractor were most associated with CP1_LI09 and CP2_LI09 cohorts, respectively. The enriched Lachnospiraceae_UCG_006 was associated with the intervention of Nostoc commune Vaucher by polysaccharides (42) . Flavonifractor plautii was capable of attenuating inflammatory responses in the obese adipose tissue (43) . As for fungi, Meyerozyma and Penicillium were most associated with CP1_LI09 and CP2_LI09 cohorts, respectively. Meyerozyma has been found in the patients with vulvovaginal candidiasis infection, while some Penicillium species have been used for the beneficial products and cheese-making (44, 45) . Some alternative fungi were also closely associated with CP2_LI09 cohort, e.g., Sporobolomyces and Fusicolla. Sporobolomyces could accumulate beneficial metabolites (46) , while Fusicolla was found as a fungus of soil origin (47) . Likewise, multiple bacteria and fungi were associated with the L I 1 0 c o h o r t s w i t h d i ff e r e n t c y t o k i n e p r o fi l e s . Lachnospiraceae_NK4A136_group and Talaromyces were the bacterium and fungus most associated with the CP1_LI10 cohort, while Parabacteroides and Aspergillus were the bacterium and fungus most associated with the CP2_LI10 cohort. Lachnospiraceae_NK4A136_group was determined with low abundance in the obese mice (48) , while Talaromyces has been found as pathogenic fungus in human beings (49) . Parabacteroides was a commensal gut bacterium and has been used to alleviate 2,4,6-trinitrobenzene sulfonic acid (TNBS)induced colitis in mice (50) . Aspergillus has been found to produce beneficial protease to improve the colonic luminal environment in rats with high-fat diet (51) . Some alternative bacteria and fungi were also closely associated with CP2_LI10 cohort, e.g., Pygmaiobacter, GCA_900066575, and Oidiodendron. More abundant Pygmaiobacter was found in the intestines of type 2 diabetes mice treated with debranched corn starch (52) . Increased gut GCA_900066575 has been determined in the highfat diet mice (53) , while Oidiodendron could improve the root biomass of Vaccinium corymbosum (54) . The effect of clinical variables or environmental factors on the microbiome has been well reported (24, 55) . However, the potential effect of microbe on the immune profile was seldom studied. In the current study, although the bacteria and fungi associated with each of the four cohort (i.e., CP1_LI09, CP2_LI09, CP1_LI10, and CP2_LI10) were determined to have varied correlations with the individual cytokines, there was no obvious difference in the correlation patterns between the CP1 and CP2 in LI09 and LI10 cohorts. However, the db-RDA results revealed that the microbes associated with each of the four cohorts seemed to influence the corresponding cytokine profiles and were likely to be associated with the formation distinct cytokine profiles. We acknowledge that the detailed mechanisms of the effects of microbes on the cytokine profiles need further investigation. The correlations between bacteria and fungi were determined in multiple studies for different objectives (56, 57) . In this study, different correlations were determined in the LI09 and LI10 cohorts with different cytokine profiles, but it seemed that no obvious difference was found in the correlation patterns or types (i.e., positive and negative) between the different cytokine profiles in the same cohorts. CoNet and fragmentation analyses have been used to investigate the microbiome networks in multiple studies (16) . Corynebacterium, Eubacterium, Papillibacter, and Cephalotrichiella were identified as the microbiome network gatekeepers in the CP2_LI09 cohort, while Pseudoflavonifractor was the only gatekeeper in the CP2_LI10 cohort. Some Corynebacterium and Eubacterium species have been determined to have beneficial potentials (58, 59) . Papillibacter was determined to have the potential to assist Enterococcus faecium in enhancing the absorption and utilization of phosphorus (60) . Pseudoflavonifractor was associated with the regulation of inflammation response in the aged mice with Listeria monocytogenes infection (61) . Lower fragmentation levels in the microbiome indicate greater co-occurrence patterns and more biotic interactions (16) . The network fragmentation levels of LI09 and LI10 cohorts with CP1 were lower than the corresponding cohorts with CP2, suggesting that more biotic interactions were in the LI09 and LI10 cohorts with CP1 than those with CP2. This could be partly supported by the finding that there are more correlations between bacteria and fungi in the LI09 and LI10 cohorts with CP1 than in the corresponding cohorts with CP2. In conclusion, the intestinal microbiome associated with distinct cytokine profiles in LI09 and LI10 pre-treated rats with liver injury was characterized. Multiple bacteria and fungi were associated with the better cytokine profiles in LI09 and LI10 "B_" and "F_" represent the microbes belonging to bacteria and fungi, respectively. Rank represents the rank of correlation number. *Represents the microbes with most correlations found in CP2_LI09 cohort but not CP1_LI09 cohort. cohorts, some of which were determined to influence the cytokine profiles in the corresponding cohorts. Their biomarker potentials in assisting with the evaluation of better cytokine profiles in the probiotics-treated cohorts deserve further investigation. The datasets presented in this study can be found in online repositories. The names of the repository and accession number(s) can be found below: https://www.ncbi.nlm.nih.gov/, PRJNA755955 and PRJNA767956. F_Udeniomyces B_Bifidobacterium "B_" and "F_" represented the microbes belonging to bacteria and fungi, respectively. Rank represented the rank of correlation number. *Represented the microbes with most correlations found in CP2_LI10 cohort but not CP1_LI10 cohort. Liver Injury and Failure in Critical Illness Hepatitis E Virus Detection in Liver Tissue From Patients With Suspected Drug-Induced Liver Injury Fecal Microbiota Manipulation Prevents Dysbiosis and Alcohol-Induced Liver Injury in Mice The Development of a Database for Herbal and Dietary Supplement Induced Liver Toxicity Cholestatic Liver Injury Induced by Food Additives, Dietary Supplements and Parenteral Nutrition Characteristics of Intestinal Microecology During Mesenchymal Stem Cell-Based Therapy for Mouse Acute Liver Injury Protective Effect of Flavonoids From Cyclocarya Paliurus Leaves Against Carbon Tetrachloride-Induced Acute Liver Injury in Mice A Crucial Strategy for Targeted Therapy of Liver Diseases Protective Effects of Lactobacillus Plantarum C88 on Chronic Ethanol-Induced Liver Injury in Mice Probiotic Bacillus Spores Protect Against Acetaminophen Induced Acute Liver Injury in Rats Lactobacillus Plantarum C88 Protects Against Aflatoxin B 1-Induced Liver Injury in Mice via Inhibition of NF-kb-Mediated Inflammatory Responses and Excessive Apoptosis Bifidobacterium Pseudocatenulatum LI09 and Bifidobacterium Catenulatum LI10 Attenuate D-Galactosamine-Induced Liver Injury by Modifying the Gut Microbiota. Sci Rep Histological Grading and Staging of Chronic Hepatitis Gut Mycobiota Alterations in Patients With COVID-19 and H1N1 Infections and Their Associations With Clinical Features New Strain of Pediococcus Pentosaceus Alleviates Ethanol-Induced Liver Injury by Modulating the Gut Microbiota and Short-Chain Fatty Acid Metabolism Bacterial Community Collapse: A Meta-Analysis of the Sinonasal Microbiota in Chronic Rhinosinusitis Protective Effects of Probiotics on Acute Alcohol-Induced Liver Injury in Mice Through Alcohol Metabolizing Enzymes Activation and Hepatic TNF-a Response Reduction Probiotic Species in the Modulation of Gut Microbiota: An Overview Gut-Liver Axis, Gut Microbiota, and its Modulation in the Management of Liver Diseases: A Review of the Literature Saccharomyces Boulardii Administration Changes Gut Microbiota and Attenuates D-Galactosamine-Induced Liver Injury Probiotics may Delay the Progression of Nonalcoholic Fatty Liver Disease by Restoring the Gut Microbiota Structure and Improving Intestinal Endotoxemia Immunomodulatory Effects of Probiotics on Cytokine Profiles Lactic Acid Bacteria and Probiotic Organisms Induce Different Cytokine Profile and Regulatory T Cells Mechanisms Multiple Bacteria Associated With the More Dysbiotic Genitourinary Microbiomes in Patients With Type 2 Ascitic Bacterial Composition is Associated With Clinical Outcomes in Cirrhotic Patients With Culture-Negative and Non-Neutrocytic Ascites Persistent Systemic Inflammation in Patients With Severe Burn Injury is Accompanied by Influx of Immature Neutrophils and Shifts in T Cell Subsets and Cytokine Profiles Cytokine Expression in Placenta-Derived Mesenchymal Stem Cells in Patients With Pre-Eclampsia and Normal Pregnancies Secreted IL-1a Promotes T-Cell Activation and Expansion of CD11b+ Gr1+ Cells in Carbon Tetrachloride-Induced Liver Injury in Mice Inflammatory Cytokines, T Lymphocyte Subsets, and Ritonavir Involved in Liver Injury IL-4 Mediates Dicloxacillin-Induced Liver Injury in Mice Impact of Paracoccidioides Brasiliensis Coinfection on the Evolution of Schistosoma Mansoni-Induced Granulomatous Liver Injury in Mice Tranilast Reduces Serum IL-6 and IL-13 and Protects Against Thioacetamide-Induced Acute Liver Injury and Hepatic Encephalopathy Depletion of Foxp3+ Regulatory T Cells Promotes Profibrogenic Milieu of Cholestasis-Induced Liver Injury A Novel Nomogram to Predict Evident Histological Liver Injury in Patients With HBeAg-Positive Chronic Hepatitis B Virus Infection Arctigenin Protects Against Liver Injury From Acute Hepatitis by Suppressing Immune Cells in Mice Bacillus Amyloliquefaciens SC06 Protects Mice Against High-Fat Diet-Induced Obesity and Liver Injury via Regulating Host Metabolism and Gut Microbiota Bifidobacterium Adolescentis CGMCC 15058 Alleviates Liver Injury, Enhances the Intestinal Barrier and Modifies the Gut Microbiota in D-Galactosamine-Treated Rats Vital Members in the More Dysbiotic Oropharyngeal Microbiotas in H7N9-Infected Patients The Microbiome of Crohn's Disease Aphthous Ulcers Characteristics of Three Microbial Colonization States in the Duodenum of the Cirrhotic Patients Urinary Microbiome and Psychological Factors in Women With Overactive Bladder Polysaccharides Isolated From Nostoc Commune Vaucher Inhibit Colitis-Associated Colon Tumorigenesis in Mice and Modulate Gut Microbiota Oral Administration of Flavonifractor Plautii Attenuates Inflammatory Responses in Obese Adipose Tissue Domestication of the Emblematic White Cheese-Making Fungus Penicillium Camemberti and its Diversification Into Two Varieties Endangerment of Ostrya Rehderiana Chun and its Relationship With Rhizosphere Soil Microflora Streptomyces Lydicus M01 Regulates Soil Microbial Community and Alleviates Foliar Disease Caused by Spermidine Improves Gut Barrier Integrity and Gut Microbiota Function in Diet-Induced Obese Mice Talaromyces (Penicillium) Marneffei Infection in Non-HIV-Infected Patients In Vitro Characterization of Gut Microbiota-Derived Commensal Strains: Selection of Parabacteroides Distasonis Strains Alleviating TNBS-Induced Colitis in Mice Beneficial Effects of Protease Preparations Derived From Aspergillus on the Colonic Luminal Environment in Rats Consuming a High-Fat Diet The Protective Mechanism of a Debranched Corn Starch/Konjac Glucomannan Composite Against Dyslipidemia and Gut Microbiota in High-Fat-Diet Induced Type 2 Diabetes Probiotic Mixture of Lactobacillus Plantarum Strains Improves Lipid Metabolism and Gut Microbiota Structure in High Fat Diet-Fed Mice Ericoid Fungal Inoculation of Blueberry Under Commercial Production in South Africa Using Soil Bacterial Communities to Predict Physico-Chemical Variables and Soil Quality Longitudinal Study of the Bacterial and Fungal Microbiota in the Human Sinuses Reveals Seasonal and Annual Changes in Fungi Participate in the Dysbiosis of Gut Microbiota in Patients With Primary Sclerosing Cholangitis Beneficial and Safety Properties of a Corynebacterium Vitaeruminis Strain Isolated From the Cow Rumen High Intake of Dietary Fructose in Overweight/Obese Teenagers Associated With Depletion of Eubacterium and Streptococcus in Gut Microbiome Enterococcus Faecium Modulates the Gut Microbiota of Broilers and Enhances Phosphorus Absorption and Utilization Aging-Induced Dysbiosis of Gut Microbiota as a Risk Factor for Increased Listeria Monocytogenes Infection The animal study was reviewed and approved by Animal Care and Use Committee of the First Affiliated Hospital, Zhejiang University School of Medicine. HZ and LL designed the study. HZ, QL, JX, SL, and RT conducted the experimental work and prepared the dataset. HZ and KC performed data analyses. HZ and LL wrote the manuscript. All research was conducted under supervision of LL. All authors reviewed the manuscript and approved the submission. Thanks to Jinling Cheng and Pei Zhang for their help in the laboratory. The Supplementary Material for this article can be found online at: https://www.frontiersin.org/articles/10.3389/fimmu.2022. Conflict of Interest: The authors declare that the research was conducted in the absence of any commercial or financial relationships that could be construed as a potential conflict of interest.Publisher's Note: All claims expressed in this article are solely those of the authors and do not necessarily represent those of their affiliated organizations, or those of the publisher, the editors and the reviewers. Any product that may be evaluated in this article, or claim that may be made by its manufacturer, is not guaranteed or endorsed by the publisher.Copyright © 2022 Zha, Li, Chang, Xia, Li, Tang and Li. This is an open-access article distributed under the terms of the Creative Commons Attribution License (CC BY). The use, distribution or reproduction in other forums is permitted, provided the original author(s) and the copyright owner(s) are credited and that the original publication in this journal is cited, in accordance with accepted academic practice. No use, distribution or reproduction is permitted which does not comply with these terms.