Taxonomic characterization of faecal microbiota
We got 192,389,261,700 base pair (bp) raw bases, with an average of 25,651,902 raw reads for each sample. After quality control, we obtained an average of 23,177,046 clean reads for each. Based on the clean reads, the composition of faecal microbiota communities was revealed at the phylum and genus level. Overall, Bacteroidetes, Firmicutes, Proteobacteria, Actinobacteria, Ascomycota, Synergistetes, Euryarchaeota, Verrucomicrobia, Spirochaetes, Mucorales, Basidiomycota and Fibrobacteres were detected in all samples. Among these phyla, Bacteroidetes, Firmicutes, Proteobacteria and Actinobacteria were the four most dominant phyla, accounting for more than 98% of microorganisms in all groups. In the UC group, Bacteroidetes increased and Firmicutes, Proteobacteria, Ascomycota, Actinobacteria, and Synergistetes decreased compared with NG rats. Compared with the UC group, the UC+MOX group had an increased abundance of Bacteroidetes, Actinobacteria, Ascomycota, Synergistetes and a decreased abundance of Firmicutes and Proteobacteria. In addition, Bacteroidetes, Actinobacteria, Ascomycota, Synergistetes were increased and Firmicutes, Proteobacteria and Deferribacteres were decreased in the UC+MES group compared with the UC group. Compared with the NG group, the NG+MOX group had a decreased abundance of Firmicutes, Proteobacteria and Actinobacteria and a increased abundance of Bacteroidetes (Fig.1A).
At the genus level, Bacteroides, Bacteroides_bacterium_M7, Bacteroides_bacterium_H3 followed by the Prevotella, Lactobacillus, Parabacteroides, Bacteroidales_bacterium_H2, Porphyromonas, Alistipes, Parasutterella, Barnesiella, Helicobacter, Clostridium, Bifidobacterium and Enterobacter were detected in all groups. Compared with the NG group, the relative abundance of Bacteroides, Bacteroides_bacterium_M7, Prevotella, Bacteroidales_bacterium_H2 and Barnesiella were increased and Bacteroides_bacterium_H3, Parabacteroides, Porphyromonas, Alistipes and Parasutterella reduced in the UC group. Interestingly, treatment of UC model rats with moxibustion led to a higher abundance of Bacteroides and Bacteroides_bacterium_H3 and a lower of Bacteroides_bacterium_M7, Prevotella, Bacteroidales_bacterium_H2 and Parasutterella. In addition, rats in the UC+MES group had a higher abundance of Bacteroides, Bacteroides_bacterium_M7, and Barnesiella and a lower abundance of Bacteroides_bacterium_H3, Prevotella, Parasutterella and Alistipes compared with the UC group (Fig.1B).
A total of 5984 species have been identified, of which 2069 are core species distributed in all groups (Fig.2). The average number of species detected in UC group was the highest (4926), while the lowest was found in NG group (2896). A total of 3611 species were detected in the UC+MOX rats, 3496 species in the UC+MES rats and 3830 species in the NG+MOX rats. Additionally, the UC group had the most unique species (800), even after treatment with moxibustion.
To gain further insight into the diversity of microbial community between groups, species-level differences was analysed (Fig.3). Achromobacter piechaudii, Achromobacter sp.LC458, Acidithiobacillus ferrivorans, Acinetobacter baumannii, Acinetobacter calcoaceticus, Acinetobacter johnsonii, Acinetobacter lwoffii, Acinetobacter rudis, Acinetobacter soli and Acinetobacter sp.BMW17 were the most affected species. Additionally, Achromobacter piechaudii and Achromobacter sp.LC458 were the absolute dominant species in the UC group, suggesting that it may be related to the state of UC disease. Compared with the NG group, Achromobacter piechaudii, Achromobacter sp.LC458 and Acinetobacter baumannii were increased and Acidithiobacillus ferrivorans, Acinetobacter johnsonii, Acinetobacter lwoffii, Acinetobacter rudis, Acinetobacter soli and Acinetobacter sp. BMW17 were decreased in the UC model group. In addition, we found significantly lower abundance of Acinetobacter baumannii in the moxibustion group than in the UC group. We also noted that Acidithiobacillus ferrivorans, Acinetobacter johnsonii, Acinetobacter lwoffii, Acinetobacter rudis, Acinetobacter soli and Acinetobacter sp. BMW17 were increased in moxibustion, mesalazine and normal rats with moxibustion groups than UC group.
A principal coordinates analysis (PCoA) was performed to compare the composition of bacterial community within 25 samples, showed a clear separation of the community composition between normal rats and UC rats (Fig.4). PCoA by distance of bray between samples explained 63.7% of variation in the first axes. The ordinations demonstrated that the disease state was the primary factor affecting the overall differences of bacterial communities between groups. The UC+MOX group was more similar to the UC+MES group on the first axis, whereas the NG+MOX group was more similar to the NG group. It was suggested that the diversity of gut microbiota was broken in UC rats. Furtherly, the flora community composition of UC rats were changed as well.
Correlation between the gut microbiota in species level and cytokines, LPS, SIgA in the UC model rats
Spearman correlation analysis between the gut microbiota and cytokines (IL-2, IL-10, IL-12, IL-17, IL-23, TGF-β, IFN-γ, TNF-α, TNFR1 and TNFR2), LPS, SIgA in serum in the UC model rats listed in Fig.5. Our previous studies showed that compared with the UC rats, the levels of IL-6, IL-12, IL-17, IL-23, IFN-γ, LPS, TNF-a, TNFR1, TNFR2 in the moxibustion group were decreased and IL-2, TGF-β were increased [17]. From the figure, we found that Bacteroides_massiliensis was negatively (P<0.05) correlated with IL-23, Bacteroides_eggerthii_CAG109 and Bacteroides_eggerthii were negatively (P<0.05) correlated with TGF-β. A significant negative (P<0.05) correlation also can be found between Helicobacter_rodentium and Bacteroides_intestinalis and IL-17. Helicobacter_rodentium was also negatively (P<0.05) correlated with TNFR1, Bacteroides_uniformis was negatively (P<0.05) correlated with IL-2. Additionally, IL-23 exhibited significant positive (P<0.05) correlation with the species Prevotella_sp_CAG1031 and Bacteroides_bacterium_H2.
Functional annotation analysis
We used metagenomics methodologies to understand the functional capacity potential of gut microbiome in the UC progression. KEGG orthologs (KOs) is a taxonomic system of proteins (enzymes) that are highly similar in sequence. In total, 30,3740 KEGG genes were found and assigned to 38 KEGG pathways. The five functional modules were metabolism (62.81%), genetic information processing (16.45%), environment information processing (10.78%), cellular processes (4.33%), human diseases (3.81%) and organismal systems (1.82%). In five metabolic pathways, the most predominant functional categories were carbohydrate metabolism (44,838 genes, 14.76%), amino acid metabolism (34,340 genes, 11.31%), nucleotide metabolism (24,483 genes, 8.06%), energy metabolism (22,533 genes, 7.42%), translation (22,389 genes, 7.37%), membrane transport (20,333 genes, 6.69%) and replication and repair (15,827 genes, 5.21%) (Fig.6A). Then we analyzed the differential KOs between NG and UC rats, and identified 33 KOs showing different abundance. Twenty out of these 33 KOs had the higher abundance in the NG rats, including those KOs associated with metabolism of cofactors and vitamins (K01633 and K03473), transport and catabolism (K04077), and carbohydrate metabolism (K15633). The other 13 KOs were enriched in UC rats, such as amino acid metabolism (K00215), translation (K02935), signaling molecules and interaction (K03699) (Fig.6B) . We also analyzed the differential KOs between UC and UC+MOX rats, and identified 12 KOs showing different abundance. Ten out of these 12 KOs had the higher abundance in the UC rats, including those KOs associated with nucleotide metabolism (K01933 and K01951) and membrane transport (K02471). The other 2 KOs were enriched in UC+MOX rats associated with Carbohydrate metabolism (K19265 and K00041) (Fig.6C). Similarly, we found that there was 25 KOs showing different abundance between UC and UC+MES rats. Nineteen out of these 25 KOs had the higher abundance in the UC rats, including those KOs associated with metabolism of cofactors and vitamins (K01906), nucleotide metabolism (K01937, K01933, K01493, K01951 and K01429) and membrane transport (K02471, K10544, K02282 and K02006). The other 6 KOs were enriched in UC+MES rats associated with Carbohydrate metabolism (K19265 and K00041) and Amino acid metabolism (K00657) (Fig.6D).
The EggNOG database is another common method for functional annotation using the Smith-Waterman alignment algorithm. We found that all genes were classified into 25 categories. Among these 25 functional classes, genes involved in certain metabolic pathways, such as energy production and conversion(25,339 and 27,102, respectively), amino acid transport and metabolism (36,275 and 38,286, respectively), nucleotide transport and metabolism (14,437 and 15,610, respectively), carbohydrate transport and metabolism (38,675 and 40,273, respectively), replication, recombination and repair (39,060 and 39,466, respectively), were under-represented in UC model rats compared with NG rats. It was suggested that the less vigorous microbial growth and metabolic stage in UC rats compared with that of normal rats. In addition, these changes in the metabolic pathways could be reversed by moxibustion treatment and mesalazine treatment (Fig.7).
KEGG annotation exposes metabolic differences between UC and moxibustion treatment rats
Markers genes involved in key processes of pyruvate metabolism that were significantly enriched in NG, UC or UC + MOX rats were annotated on KEGG pathway (Fig.8). Oxaloacetate, which is an important intermediate in the tricarboxylic acid cycle, is used to produce pyruvate by oxaloacetate decarboxylase (EC 4.1.1.3). Compared with the NG group, pyruvate-producing genes (ECs 4.1.1.3 and 1.1.1.28) were reduced in the UC group. Compared with the UC group, pyruvate-producing genes (EC 1.1.1.28) was increased in the UC + MOX group. In the UC group, EC 1.2.4.1 involved in the reduction of pyruvate to acetyl-CoA which is a central metabolite of energy metabolism was also reduced compared with the NG group. EC 6.4.1.2 (Acetyl-CoA carboxylase), reduced in the UC group, was responsible for Acetyl-CoA to become malonyl-CoA which plays an important role in the biological synthesis of fatty acids and the biosynthesis of polyketones. Compared with the UC group, EC 6.4.1.2 (Acetyl-CoA carboxylase), increased in the UC + MOX group. Fatty acid are metabolites of the intestinal microflora and have the function of anti-pathogenic microorganisms and inhibit inflammatory responses[25]. In summary, the metabolism of fatty acids in the UC group was decreased compared with the NG group, and increased in UC + MOX group to some extent compared with the UC group.