Human Colon Mucosal Biofilms and Murine Host Communicate via Altered mRNA and microRNA Expression during Cancer

Bacteria and bacterial biofilms have been implicated in colorectal cancer (CRC), but it is still unclear what genes these microbial communities express and how they influence the host. MicroRNAs regulate host gene expression and have been explored as potential biomarkers for CRC. An emerging area of research is the ability of microRNAs to impact growth and gene expression of members of the intestinal microbiota. This study examined the bacteria and bacterial transcriptome associated with microbes derived from biofilm-positive human cancers that promoted tumorigenesis in a murine model of CRC. The murine response to different microbial communities (derived from CRC patients or healthy people) was evaluated through RNA and microRNA sequencing. We identified a complex interplay between biofilm-associated bacteria and the host during CRC in mice. These findings may lead to the development of new biomarkers and therapeutics for identifying and treating biofilm-associated CRCs.

KEYWORDS colorectal cancer, transcriptomics, microRNAs, microbiota, biofilm, germfree N umerous 16S rRNA and shotgun metagenomic studies have demonstrated that colorectal cancer (CRC) patients have an altered intestinal microbiota compared to healthy controls (1,2). Colibactin-producing Escherichia coli, enterotoxigenic Bacteroides fragilis, and Fusobacterium nucleatum among others, are implicated in CRC pathogenesis due to their abilities to produce genotoxins and adhesins which promote proliferation and modulate immune responses in preclinical models (3,4). How these bacteria interact with the rest of the complex microbiota to influence CRC initiation and/or progression is still unclear. In fact, recent studies with F. nucleatum suggest these bacteria may be associated with later stages of disease and have less of an influence on CRC initiation (5)(6)(7)(8). Testing the functional role of human CRC-associated bacterial communities in chemically induced mouse models of CRC have led to mixed results (9,10). One group demonstrated stool communities from either individual CRC or healthy patients promoted polyp formation depending on the composition of the microbiota that established in mice (9). Another recent report revealed an increased tumorigenic phenotype in mice that received stools pooled from multiple CRC patients compared to stools from controls (10).
The lack of a consensus carcinogenic CRC-associated microbiota from patient stools suggests that other factors, including how the bacteria are organized/located or the genes they express, may be just as important to CRC pathogenesis. Polymicrobial bacterial biofilms, spanning Ͼ200 m of epithelial surface were recently identified in ϳ52% of human CRC patients (11) and also found in ϳ13% of the healthy patients who were screened (12). We previously showed that human biofilm-forming bacterial communities from either CRC or healthy patients play a functional role in CRC development in multiple preclinical mouse models, emphasizing the contribution of bacterial organization to CRC (13).
CRC is an evolving disease, characterized by a series of molecular and microbial changes (14)(15)(16), suggesting a dynamic interplay between the host and intestinal microbiota as the disease progresses. MicroRNAs (miRNAs) have emerged as potential mediators of these host-microbe interactions with their ability to modulate both host (17) and bacterial genes, which can result in shifts in microbiota composition (18,19). In turn, the microbiota is able to modulate host miRNA expression (18,20,21), with F. nucleatum targeting several miRNAs related to CRC pathogenesis (5,22). However, it is uncertain how human CRC-associated microbial communities as a whole impact fecal miRNA expression and whether host miRNAs affect bacterial composition/gene expression during CRC.
To examine the bacterial activities associated with biofilm-positive microbes from CRC patients, we examined mouse and bacterial gene expression from colon tissues and mouse small RNA sequencing from stools collected from biofilm-positive associated Apc MinΔ850/ϩ ;Il10 Ϫ/Ϫ mice. We found that a number of bacterial virulence genes were increased in biofilm-positive communities and identified a conserved core group of transmissible biofilm-positive associated bacteria. Additionally, we demonstrate that biofilm status and CRC development alter miRNA expression and specific miRNAs correlate with biofilm-positive associated taxa.

RESULTS
Bacterial activities associated with biofilm status. In order to elucidate microbial activities associated with biofilm-forming bacteria derived from human CRC patients that promote tumorigenesis in Apc MinΔ850/ϩ ;Il10 Ϫ/Ϫ mice (13), we characterized mouse and microbial gene expression from colon tissue snips using RNA sequencing (see   Table S1 at https://figshare.com/s/bd06b560c635de3ac830]) compared to BF-bx mice. Pathways related to protein export, bacterial secretion systems, carbon fixation, flagellar assembly, and biosynthesis of amino acids were increased in BFϩT mice compared to BF-bx mice (Table 1 and Fig. S2B and C). Additional genes related to virulence and biofilm formation, including stress response, toxins, iron acquisition, mucin cleavage/transport, outer membrane polysaccharide importers, and adhesins were also significantly increased in BFϩT mice (see Table S1 at https://figshare.com/s/ bd06b560c635de3ac830). Increased toxin genes included Clostridium difficile toxins A and B, Clostridium perfringens Mu toxin, and E. coli colibactin (clbG and clbI). Weighted gene coexpression network analysis (23) identified 34 hub genes (see Table S1 at https://figshare.com/s/bd06b560c635de3ac830) from modules detected in BFϩT mice and included outer membrane proteins involved in protein export and heat shock proteins involved in the stress response. Despite the high number of microbial genes with increased expression in BFϩT microbial communities, no separation by PCA analysis was found at the host gene expression level, and only 62 significant DE genes ( Fig. 2A; see also Table S2 at https://figshare.com/ s/4b593b780f756a4acf69) between BFϩT and BF-bx mice were detected. Instead, the host was more responsive to the microbiota in general as opposed to the type of microbiota, since the host transcriptomes of either BFϩT-or BF-bx-associated mice clustered separately from germfree (GF) mice ( Fig. S1A and B). There were Ͼ3,000 significant DE host genes (ϳ2,000 upregulated, ϳ1,300 downregulated) in the BFϩT and BF-bx groups compared to GF mice (see Tables S3 at https://figshare.com/s/ 163676e591c87b0d6c35, S4 at https://figshare.com/s/652055bdb15e48b866ef, and S5 at https://figshare.com/s/c9adfcd56af278666c55) and pathway analysis revealed that the majority of upregulated genes in colonized mice belonged to immune-related pathways ( Fig. S1C and D). Only the peroxisome proliferator-activated receptor (PPAR) signaling pathway was significantly upregulated (P ϭ 1.88eϪ05, P FDR ϭ 0.004) in BFϩT and reassociation (III) experiments, along with the analyses done on the stool and tissue samples (II and IV) at the end of the 12-week experiments. Twelve-week stool and/or DC tissue samples were used for RNA, miRNA, and 16S rRNA sequencing analyses (II). Tissue was collected from 12-week-associated BF-bx mice and 16-to 20-week-associated BFϩT mice to make the reassociation inoculums (III). (B) PCA of bacterial transcriptomes from BFϩT-and BF-bx-associated Apc MinΔ850/ϩ ;Il10 Ϫ/Ϫ mice generated from Trinity de novo assembly (N ϭ 5 for BFϩT and BF-bx). mice compared to BF-bx mice. The majority of the significant DE genes were upregulated (56 versus 6 downregulated) in BFϩT mice and included genes related to lipid metabolism, iron scavenging, and the extracellular matrix (see Table S2 at https:// figshare.com/s/4b593b780f756a4acf69). Cross-referencing the list of upregulated genes with The Cancer Genome Atlas (TCGA) colorectal microarray data using Oncomine (24) revealed 10 of the upregulated genes were also significantly overexpressed in colon adenocarcinomas compared to control tissues with SCD (stearoyl-coenzyme A [CoA] desaturase, PPAR pathway member), MMP10, and SLC22A3 increased more than 2-fold (see Table S6 at https://figshare.com/s/7248a5727e5972879086). Biofilm status alters the fecal miRNA profile. Given the differential bacterial and host gene expression observed in BFϩT mice that developed colon tumors and the potential of miRNAs to modulate interkingdom interactions (18), we next profiled stool miRNA expression. To examine the interplay between miRNAs and bacterial communities derived from human CRC or healthy patients, we sequenced the fecal small RNAs from GF, BF-bx, and BFϩT Apc MinΔ850/ϩ ;Il10 Ϫ/Ϫ mice. The PCA plot shows clear separation between GF mice and BF-bx or BFϩT-associated mice, demonstrating that microbial colonization alters the fecal miRNA profile in Apc MinΔ850/ϩ ;Il10 Ϫ/Ϫ mice ( Fig. 3A; see also Table S7 at https://figshare.com/s/b7d7ff960f787eb31954). Biofilm/ cancer status of the initial human-derived microbes also modulates miRNA expression since PCA analysis demonstrates separation of BFϩT and BF-bx miRNAs ( Fig. 3A; see also Table S7 at https://figshare.com/s/b7d7ff960f787eb31954). Pairwise comparisons between the three groups of mice revealed that 25 unique miRNAs were significantly DE (out of 142 total detected) ( Fig. 3B; see also Tables S7 at https://figshare.com/s/b7d7ff960f787eb31954, S8 at https://figshare.com/s/2eaad72ee a8a5890f7f2, and S9 at https://figshare.com/s/20a8d9c108e990242bc7). Next, we compared the significantly different miRNAs in BFϩT or BF-bx versus GF mice and found that nine significant DE miRNAs overlapped (mmu-miR-6538, -146b-5p, -215-5p, -194-5p, -192-5p, -2137, and -5126 and mmu-let-7b-5p and mmu-let-7i-5p), suggesting host miRNAs targeting the microbiota ( Fig. 3B; see also Tables S8 at https://figshare.com/s/ 2eaad72eea8a5890f7f2 and S9 at https://figshare.com/s/20a8d9c108e990242bc7). We were also able to identify eight miRNAs (mmu-miR-709, -690, -21a-5p, -142a-5p, -6240, -6239, and -148a-3p and mmu-let-7a-5p) that were significantly DE according to biofilm status (BFϩT versus BF-bx [ Fig. 3B; see also Table S7 at https://figshare.com/s/ b7d7ff960f787eb31954]). These findings suggest that the microbiota modulates host miRNA expression. Together, both microbiota composition and disease status drive miRNA expression (BFϩT versus BF-bx DE miRNAs).
Fecal miRNAs correlate with specific bacterial taxa and target bacterial and mouse genes. We next examined whether the fecal miRNAs of human microbiotaassociated mice correlated with bacterial taxa previously identified in the mice (13). We identified 11 miRNAs that were significantly DE between GF, BF-bx, and BFϩT mice and correlated with the relative abundances of bacterial taxa in the stool, distal colon tissue, or both compartments (Fig. 3B [colored green, orange, or underlined]; see also Tables S10 at https://figshare.com/s/55e156ec52d9d3e53394 and S11 at https://figshare.com/s/1c47b438d806062915c9). Additionally, we identified miRNAs that although not significantly DE between groups, correlated with five and eight bacterial genera in the stool and distal colon tissues, respectively ( Fig. 4A and B; see also Tables S10 at https://figshare.com/s/55e156ec52d9d3e53394 and S11 at https:// figshare.com/s/1c47b438d806062915c9). Five of these genera (Bacteroides, Lachnospiracea incertae sedis, Anaerostipes, Clostridium XVIII, and Roseburia) were significantly increased in BFϩT mice ( Fig. 4A and B; see also Tables S10 at https://figshare.com/s/ 55e156ec52d9d3e53394 and S11 at https://figshare.com/s/1c47b438d806062915c9). Furthermore, mmu-miR-140-3p correlated with Lachnospiraceae incertae sedis abundance and tumor number, both of which were increased in BFϩT mice ( Fig. 4A and Fig. S3A). Thus, CRC-associated microbial communities elicit a specific host miRNA profile and maintain a subset of miRNAs that correlate with specific microbial taxa. Interestingly, computationally predicting the bacterial genes targeted by the significant miRNAs revealed several miRNAs (mmu-miR-2137, mmu-miR-5126, and mmu-miR-6538) that primarily target bacterial genes (see Table S12 at https://figshare.com/ s/33fb65da16763a5e3347). In contrast, other miRNAs (mmu-let-7s, mmu-miR-21a-5p, mmu-miR-142a-5p, mmu-miR-148a-3p, mmu-miR-194-5p, mmu-miR-690, and mmu-miR-709) are predicted to primarily target mouse genes (see Table S13 at https:// figshare.com/s/be3f38607b2e9875e7f6), and there is a significant negative correlation between the number of predicted bacterial versus mouse gene targets for the significant miRNAs we identified ( Fig. 4C and D and Fig. S3B and C). Taken together, these predictions suggest that specific miRNAs have differential roles in mediating host- BFϩT microbiota elicited more host miRNAs than the GF and BF-bx microbiota with six miRNAs uniquely DE in response to it. The arrows next to miRNA names indicate whether the miRNA is increased (1) or decreased (2) in the first group listed relative to the second for each comparison. Blue arrows indicate direction of miRNA expression in BFϩT mice compared to BF-bx mice. Pink arrows indicate direction of miRNA expression in BFϩT mice compared to GF mice. Purple arrows indicate direction of miRNA expression in BF-bx mice compared to GF mice. Black arrows indicate direction of miRNAs that are shared between comparison groups and go in the same direction. The miRNAs that are significantly correlated with taxa identified from 16S rRNA sequencing of the stool are highlighted in green, those that significantly correlated with taxa identified from the DC tissue are highlighted in orange, and those correlating with both are in black font and underlined (total of 11 miRNAs).

DC tissue
Genus Relative Abundance (Continued on next page) Tomkovich et al. microbe interactions, with some mainly modulating host gene expression and others primarily impacting bacterial gene expression and/or abundance. A core BF؉T microbiota is transmissible. We previously showed that carcinogenic properties are retained over time as microbial inoculums from homogenized mouse colon tissues collected from the initial BFϩT but not BF-bx tissue-associated microbes promoted tumors in new cohorts of GF Apc MinΔ850/ϩ ;Il10 Ϫ/Ϫ mice (13) (Fig. 1A.III).

Lachnospiracea incertae sedis
We next examined microbial compositional differences between BFϩT and BF-bx reassociated mice. 16S rRNA sequencing revealed separation of stool and distal colon (DC) tissue microbiota from reassociated mice according to the biofilm status of the initial association ( Fig. 5A and B). Thirteen genera were significantly different (12 genera enriched and one genus depleted) in the stool and/or DC tissue of both the BFϩTassociated and reassociated Apc MinΔ850/ϩ ;Il10 Ϫ/Ϫ mice (Fig. 5C), three of which also correlated with specific miRNAs (Fig. 4A and B). Six of these genera (Clostridium XVIII, Erysipelotrichaceae incertae sedis, Escherichia/Shigella, Eubacterium, Parabacteroides, and Robinsoniella) were increased in both the stool and DC tissue of BFϩT-associated and reassociated mice compared to BF-bx mice. To determine how much the microbiota composition shifted after reassociation, we compared the microbiotas from the reassociated mice to the initial associated mice, whose tissues were used to generate the inoculums ( Fig. S4 and S5). The BF-bx communities were highly transmissible, with no significant differences based on principal-coordinate analysis (PCoA) of the stool and DC tissue communities (Fig. S4). The BFϩT communities shifted more after reassociation, with distinct clustering seen in the PCoA ( Fig. S5A and B). Depending on the colon region the BFϩT inoculum was derived from, there were 5 and 13 significantly different genera within the DC tissue or stool compartment; however, only 1 and 5 of these genera were significantly different based on biofilm status in the initial associations ( Fig. S5C and D). We also examined how the location of the colon tissues used to make the BFϩT reassociation inoculums impacts community composition by comparing the stool and DC tissue communities from mice reassociated with PC or DC tissue inoculums (Fig. S6A). We found that only 2 and 5 genera were significantly different (Fig. S6B) and only 2 of these genera (Coprobacillus and Holdemania) differed according to biofilm status in the initial associations. Out of the 12 genera that were increased in BFϩTassociated and reassociated mice, only 1 (Coprobacillus) was not maintained in both groups of BFϩT reassociated mice. Thus, regardless of the murine colon region of origin (proximal or distal), the majority of BFϩT microbes (11 genera total; Anaerostipes, Clostridium XI, Clostridium XIVa, Clostridium XVIII, Erysipelotrichaceae incertae sedis, Escherichia/Shigella, Eubacterium, Flavonifractor, Lachnospiraceae incertae sedis, Parabacteroides, and Robinsoniella) are able to reestablish and promote cancer when transmitted to a new cohort of GF mice. Taken together, the data suggest there is a core set of bacteria and bacterial gene expression associated with biofilm-positive cancers.

DISCUSSION
In contrast to previous studies that tested the carcinogenicity of CRC-associated microbiotas by gavaging mice with stools from either CRC or healthy control patients (9, 10), mucosa-associated BFϩT bacteria retain their carcinogenicity when transplanted into a new set of GF mice, regardless of whether they were extracted from the proximal or distal colon (13). The number of taxa that overlap between the initial The name of each miRNA is shown below the genus it correlates with. Genera in blue font are significantly different between BF-bx and BFϩT Apc Min⌬850/ϩ ;Il10 Ϫ/Ϫ mice. The red underlined genera were significantly different based on biofilm status in both the initial association and reassociation experiments. The direction of correlation is shown within parentheses. Relative abundance data are from the subset of BF-bx and BFϩT mice that were used for miRNA sequencing (n ϭ 7 for BF-bx; n ϭ 10 for BFϩT). See Tables S10 at https://figshare.com/s/55e156ec52d9d3e53394 and S11 at https:// figshare.com/s/1c47b438d806062915c9 for the full list of miRNAs that correlate with bacterial taxa and the corresponding correlation coefficients. (C) Heatmap comparing the log 2 -transformed number of predicted bacterial versus mouse gene targets for the miRNAs that were significantly DE between the GF, BF-bx, or BFϩT group. There is a significant negative correlation (Pearson) between the number of predicted bacterial versus mouse gene targets for the set of significant DE miRNAs. (D) Scatter plot demonstrating the significant negative Pearson correlation between the log 2 -transformed number of mouse versus bacterial gene targets where each circle represents a unique mouse miRNA that was significantly DE between GF, BF-bx, or BFϩT group. association and the reassociation experiments, shown herein (Fig. 5), suggest that there is a core, transmissible, cancer-promoting microbiota associated with biofilm-positive cancers.

Significant genera in reassociated Apc
Four of the BFϩT core transmissible cancer-promoting genera overlap with human stool-derived taxa that established and were associated with cancer in mice. Parabacteroides correlated with high tumor numbers in a gnotobiotic azoxymethane (AOM)/ dextran sulfate sodium (DSS) C57BL/6 model (9), and Erysipelotrichaceae, Lachnospiraceae, and a Clostridium XIVa derived from CRC patients were associated with polyp formation in an antibiotic-treated AOM C57BL/6 model (10). Additionally, increased Erysipelotrichaceae, E. coli, Lachnospiraceae, and Parabacteroides and decreased Bifidobacterium have previously been associated with human CRC patient samples by 16S rRNA gene sequencing studies (1,25).
Metagenomic predictions generated from 16S rRNA data identified bacterial secretion systems and motility genes associated with human CRC stool communities (26) and host glycan utilization genes correlated with tumor numbers in human stoolassociated AOM/DSS mice (9). Recent meta-analysis studies of metagenomes from CRC patient fecal samples identified gluconeogenesis, mucin degradation, and colibactin as associated with the CRC microbiome (27,28). These bacterial pathways and genes were also increased in our BFϩT metatranscriptome. PICRUSt analysis of biofilm-positive versus biofilm-negative human CRC tissues also demonstrated an increased sporulation capacity associated with biofilm-positive CRCs contributed by several taxa, among them the Lachnospiraceae family (11). Similarly, herein, we saw genera from the Lachnospiraceae family (Anaerostipes, Clostridium XIVa, Lachnospiracea incertae sedis, and Robinsoniella) and upregulation of multiple sporulation genes in mice transplanted with the BFϩT community. There is also overlap between BFϩT community gene expression and human periodontitis polymicrobial metatranscriptomes, particularly for genes related to iron acquisition, flagellar synthesis, and the stress response (29). These findings are interesting since biofilms have also been associated with periodontal disease and F. nucleatum and Porphyromonas spp. are associated with both oral biofilms and colon cancer (4,11,25,30).
Multiple genes related to nutrient, envelope, DNA damage, and environmental stress responses were increased in the BFϩT community that could be indicative of host immune pressures, but could also be associated with a competitive polymicrobial environment, a feature of biofilms (31). Host iron metabolism changes during inflammation and cancer can promote competition for iron within the intestinal microbiota (32). Multiple iron acquisition genes, including siderophores and transport receptors, were increased in BFϩT mice (32). The expression of these metabolic and iron acquisition genes could be indicative of a low-nutrient environment, fostering interbacterial competition.
Bacterial adhesion genes are a critical colonization determinant and may also contribute to biofilm formation, in which attachment to host cells represents a key initiating step (33). There are a number of adhesins expressed in the BFϩT community, including type I and IV pili, capsule genes, and proteins that bind to host extracellular matrix (ECM) components such as fibronectin and laminin. On the host side, BFϩT mice exhibit upregulation of a laminin subunit and the ECM-degrading matrix metalloproteinase MMP10 (34), suggesting alterations to the host ECM. BFϩT communities also expressed numerous moonlighting adhesins (such as flagellin, GroEL, DnaK, and elongation factor Tu), putative multifunctional proteins which have been demonstrated in some bacterial strains to bind host cells, mucin, or ECM components (35). Bacterial (C) Heatmap depicting the mean log 10 -normalized relative abundances of genera that were significantly different in the stool and/or DC tissues of reassociated mice inoculated with murine colon tissue homogenates derived from human BFϩT or BF-bx tissue-associated mice. The underlined genera were significantly different based on biofilm status in both the initial association and the reassociation, and the red underlines represent the three genera that correlated with specific miRNAs (Fig. 4A and B).

Host-Microbiota Genomic Interactions during CRC in Mice
January/February 2020 Volume 5 Issue 1 e00451-19 msystems.asm.org 11 adherence has also been identified as an important feature of CRC-associated bacteria, including F. nucleatum (36) and Streptococcus gallolyticus subsp. gallolyticus (37), the latter of which is capable of forming biofilms on collagen IV, an ECM component (38). Furthermore, bacteria expressing adhesins that bind to ECM and host glycoproteins may have an additional colonization advantage, as host glycosylation is disrupted during inflammation and cancer with increases in sialylation and fucosylation that can result in decreased host cell adhesion to ECM components (39). Additionally, some of the effects of CRC-associated bacteria may be contact dependent; for example, colibactin-induced DNA damage requires direct contact between the bacteria and epithelial cells (40). Coupled with its role in facilitating colonization and attachment, mucin also represents a source of nutrition for intestinal bacteria (41). Antibiotic treatment was shown to increase sialic acid levels in the lumen and promoted the expansion of pathogenic bacteria such as C. difficile and Salmonella enterica serovar Typhimurium (42). Interestingly, sialic acid and other mucin sugar cleavage and transport expression were increased in BFϩT mice along with an increased abundance of Clostridium XI and Salmonella. The SusC and SusD outer membrane proteins, involved in oligosaccharide binding and transport (43), were also increased in the BFϩT community (from Bacteroides and Parabacteroides spp.). The upregulation of stress responses, mucin, and other nutrient acquisition genes indicate environmental conditions that could in turn promote virulence expression. Nutrient-and iron-responsive global transcriptional factors such as cyclic AMP receptor protein and ferric uptake regulator (Fur), which were increased in BFϩT mice, have also been shown to regulate bacterial virulence expression (44).
Shotgun metagenomic sequencing of patient stools has revealed that host glycan utilization and virulence factor genes are associated with the CRC microbiome (16) and that genes in these categories were also overexpressed in the BFϩT microbial community. Host inflammation, bacterial iron (Fur), and stress response (Hsp90 chaperone) genes have all been implicated in colibactin regulation, and all of these genes were increased in BFϩT mice (45)(46)(47). Additionally, iron acquisition genes have previously been associated with E. coli mucus colonization (48), and mucins have the capacity to induce E. coli virulence gene expression (41). Mucin and its components may also serve as a cue for virulence regulation of other BFϩT community members, as they have also been linked to virulence regulation in S. enterica and Campylobacter jejuni (41). Although the expression of multiple B. fragilis genes were increased in the BFϩT community, B. fragilis toxin (bft) was not detected. One possible explanation is that expression of the RprXY two-component system, recently implicated in bft suppression (49), was significantly increased in BFϩT mice. However, even intermittent enterotoxigenic B. fragilis colonization as short as 2 weeks appears to be sufficient to induce tumor formation (50), so it is possible that bft expression occurred at an earlier time point in the Apc MinΔ850/ϩ ; Il10 Ϫ/Ϫ model. Taken together, the metatranscriptomic data suggest that the BFϩT community expresses more pathogen-related virulence factors and metabolism genes that provide competitive advantages over commensals but may have detrimental side effects to the host.
Members of the Erysipelotrichaceae family, part of the core transmissible bacteria in BFϩT mice, were also increased in the microbiota of Western or high-fat diet-fed mice (Clostridium innocuum, Eubacterium dolichum, and Clostridium ramosum) and were associated with increased fat deposition (51,52). Diets high in fat and obesity are established risk factors for CRC (53). Conceivably, the activation of the PPAR signaling pathway within BFϩT mice could be a response to colonization with these obesityassociated taxa. PPARs are nuclear hormone receptors that regulate key aspects of lipid and carbohydrate metabolism, including fatty acid synthesis, uptake, and storage (54), and have both suppressive and promotional effects in CRC (55).
Branched-chain amino acid (BCAA) biosynthesis was increased in the BFϩT community, and serum BCAA have also been associated with metabolic disorders and correlated with intestinal microbiota members such as Bacteroides vulgatus (56). BCAAs, which include valine, leucine, and isoleucine, may also contribute to fatty acid synthesis, which is also a feature of cancer metabolism (57). Valine and/or leucine secretion were previously associated with E. coli (58) and polymicrobial environmental biofilms (59). Polyamines were previously shown through metabolomics to be increased in CRC patients, with a rare polyamine, N1, N12-diacetylspermine detected in biofilm-positive CRCs (60), and we found that multiple microbial polyamine-related genes were increased in BFϩT mice (61). Polyamines increase in proliferating cells and can promote tumor growth and invasion (62) and are also important to bacterial biofilm formation (61). Several transporter genes (Slc22a3, Abcb1a, and Abcb1b) that have been previously implicated in polyamine uptake (63) were upregulated in BFϩT mice. Thus, biofilmassociated communities and their associated metabolism pathways have the potential to modulate host metabolism, which may promote cancer.
In addition to bacterial components and metabolites, miRNAs represent another mode of host-microbe interplay during cancer. Profiling the fecal miRNAs of Apc MinΔ850/ϩ ;Il10 Ϫ/Ϫ mice under different microbial conditions allowed us to identify specific miRNAs that were associated with biofilm/CRC status. A few of the CRC-associated miRNAs we identified have conserved sequences with human miRNAs (hsa-miR-21-5p, hsa-miR-142-5p, and hsa-miR-146a-5p) that are increased in CRC patients (17,22,64,65). Mmu-miR-21a-5p was significantly increased in the BFϩT mice, and F. nucleatum has previously been shown to increase miR-21 (22), suggesting that miR-21 may be targeted by multiple CRC-associated bacterial genera.
Although the depleted miRNAs in BFϩT mice, miR-690 and miR-709, are not found in humans, they do share several CRC-related gene targets (such as Ctnnb1, Il6ra, Stat3, Src, and Zeb1) with other miRNAs depleted in CRC (17). Though none of these genes were significantly DE according to biofilm status in the mouse colon tissue RNA-seq data set, it is possible the luminal miRNAs might target other regions of the colon or specific intestinal epithelial cell types. Additionally, even though miRNAs have been shown to primarily control gene expression through mRNA degradation, translational repression is also possible (66). Alternatively, another mechanism of miRNA regulation could relate to targeting the RNA-induced silencing complex genes like Argo1, Argo2, Argo3, Argo4, Cnot6, and Dcp2 (66), which are predicted targets of multiple significant DE miRNAs, the majority of which are increased in BFϩT mice. miRNAs also have the capacity to target bacterial genes and impact microbial composition (18), and our computational predictions indicate that newly discovered miRNAs (those with higher numbers in their names) preferentially target bacterial genes. These miRNAs include miR-2137, miR-5126, miR-6239, miR-6240, and miR-6538, which were also increased in DSS-treated mice (67), where microbiota composition also contributes to disease susceptibility (68). Many of these miRNAs were predicted to have redundant bacterial targets (including genes regulating motility, secretion, outer membrane proteins, stress response, iron acquisition, and carbohydrate utilization/transport) that overlap with genes that were increased in the BFϩT microbial community. miR-6239 and -6240 were decreased in BFϩT mice, but miR-2137, -5126, and -6538 were increased in mice, regardless of biofilm status.
Conclusions. The mechanisms by which miRNAs enter and regulate bacterial gene expression and how bacterial organization affects bacterial community gene expression warrant further investigation. A recent publication showed that plant-derived exosomelike nanoparticles contain miRNAs that altered intestinal microbial composition and gene function (69). Whether intestinal epithelial cell-derived exosomes carrying miRNAs could similarly target intestinal microbiota is under investigation. Our findings suggest a complex interplay between BFϩT-associated bacteria, their gene expression, the host transcriptome, and miRNAs that may contribute to CRC pathogenesis (Fig. 6). Deciphering this complex interplay will likely identify new regulatory pathways and molecules with potential therapeutic implications.

MATERIALS AND METHODS
Animals. Germfree (GF) 129/SvEv Apc MinΔ850/ϩ ;Il10 Ϫ/Ϫ mice were transferred to separate gnotobiotic experimental isolators based on inoculum type for the duration of the association.
Initial associations with human tissue-associated microbes. GF 129/SvEv Apc MinΔ850/ϩ ;Il10 Ϫ/Ϫ mice were inoculated with pooled tissue-derived microbes from biofilm-negative tissues collected from healthy patients via colonoscopy biopsy (BF-bx) or biofilm-positive tumor tissues collected from CRC patients (BFϩT) via surgical resection. Patient tissues were collected and screened for biofilm status via fluorescence in situ hybridization (FISH) (biofilm-positive criteria, polymicrobial, within the mucus layer, spanning 200 m, and Ͼ10 9 bacteria/ml) as described previously (12). Tissues were analyzed with the universal bacterial probe (EUB338), and a nonsense probe (non338) was used as a negative control (12,13). Additional FISH analysis was conducted on biofilm-positive tumor tissues with probes to detect Bacteroidetes, Lachnospiraceae, Fusobacteria, and Proteobacteria (12,13). The probe sequences are listed in Table S4 in reference 12. Each inoculum was prepared anaerobically by homogenizing tissue (tissue pooled from five different patients) in phosphate-buffered saline (PBS), and FISH images for these tissues can be found in Fig. S1 in reference 13. Each mouse received 100 to 200 l of inoculum, and associations were carried out for 12 to 20 weeks (Fig. 1AI). Tissues and/or stools from mice collected 12 weeks after association were used for transcriptome sequencing (RNA-seq), microRNA sequencing (miRNA-seq), and 16S rRNA gene sequencing analyses (Fig. 1AII).
Mouse reassociation inoculums. Mouse reassociation inoculums (Fig. 1AIII) were made from colon tissues from 12-week BF-bx-associated (cohort 2) and 16-to 20-week BFϩT-associated (cohort 3) Apc MinΔ850/ϩ ;Il10 Ϫ/Ϫ mice (13). After the colon was flushed 1ϫ with PBS, tissue snips were taken from both the distal colon (DC) and proximal colon (PC) and stored at -80°C until time of inoculum preparation. Each inoculum was prepared from colon tissue snips pooled from four mice. All inoculums were prepared anaerobically by mincing and homogenizing tissue snips in prereduced PBS. The BF-bx reassociation inoculum was prepared from inflamed (average inflammation score of 2.5/4) distal colon tissues (BF-bx DC). The two BFϩT reassociation inoculums were prepared from mostly normal (average PC inflammation score of 0.9/4) proximal colon tissues (BFϩT PC) or distal colon tissues (BFϩT DC) from the same four mice with colitis and tumors (average DC inflammation score of 3.6/4; average number of tumors ϭ 5.5, range ϭ 3 to 10 tumors).
Reassociation with mouse tissue-associated microbes. Six-to 13-week GF 129/SvEv Apc MinΔ850/ϩ ; Il10 Ϫ/Ϫ (males and females) were transferred to gnotobiotic isolators (separate isolator for each exper-  Lipid metabolism, Laminin, Matrix metalloproteinase, Iron scavenging FIG 6 Schematic depiction of the major findings. The core, transmissible bacteria found in biofilm-positive tumor (BFϩT) associated and reassociated mice are listed under core bacteriome. Some of the bacterial and mouse genes that were differentially expressed in the BFϩT associated mice compared to biofilm-negative (BF-bx) associated mice are listed. Fecal miRNAs were also differentially expressed according to biofilm status, correlated with the relative abundances of some bacterial taxa, and were predicted to target mouse and bacterial genes.
imental group) and gavaged with 100 to 200 l of inoculum (Fig. 1AIII). Mice were euthanized after 12 weeks, and the colon was flushed once with PBS and then cut and splayed longitudinally for macroscopic tumor counts. About 2 ϫ 0.5 cm tissue snips from the proximal and distal colon were collected by flash freezing in liquid nitrogen and stored at -80°C until analysis.
Stool and distal colon tissue DNA extraction. Stool and DC tissue DNA ( Fig. 1AII and IV) was extracted via phenol-chloroform separation by lysing the cells with phenol-chloroform-isoamyl alcohol (25:24:1) and 0.1-mm zirconium glass beads on a bead beater (Precellys), followed by phase separation with chloroform-isoamyl alcohol (24:1), DNA precipitation with ethanol, and subsequent purification with the DNeasy Blood & Tissue kit (Qiagen catalog no. 69506).
16S rRNA sequencing analysis. Reads were preprocessed using Quantitative Insights into Microbial Ecology (QIIME) (70) version 1.9.1 including trimming and filtering at Q20. The final set of reads was fed to the RDP (Ribosomal Database Project) classifier (71) version 2.12 with the confidence set at 80%. Reads were grouped by genera, and samples with less than 1,000 total reads and genera with less than 5 reads were removed. The resulting counts were normalized and log 10 transformed (72) using the following formula: where RC is the read count for a particular taxon in a particular sample, n is the total number of reads in that sample, the sum of x is the total number of reads in all samples, and N is the total number of samples. The principal-coordinate analysis (PCoA) was generated from the Bray-Curtis distance of the normalized and log 10 -transformed counts using the phyloseq (73) R package (74). Genera significant for biofilm group (BF-bx, BFϩT, BF-bx DC, BFϩT PC, and BFϩT DC) were detected using the lme function in the R nlme package, with the REML method (75) to fit a mixed linear model of the form: genus ϳ variable ϩ 1|cage ϩ where genus indicates the log 10 normalized abundance of a particular genera, variable indicates either the biofilm group or PCoA axis, and 1|cage indicates that we used the cage as a random effect. We then ran an analysis of variance (ANOVA) on the above model to generate P values for biofilm group or PCoA axis. We checked for possible cage effect by comparing the above model and a model with the cage removed (genus ϳ variable ϩ ) using ANOVA. The P values were adjusted for multiple hypothesis testing in R using the p.adjust function employing the method of Benjamini and Hochberg (76). The heatmaps were generated using the R function ggplot2 (77).
We performed two additional analyses on the 16S rRNA data, the first utilizing QIIME (70) v. 1.9.1 closed-reference at 97% similarity level using the Greengenes reference data set release 13_8 and the second employing Deblur (78) workflow v. 1.0.3 with the default parameters (using Deblur's default positive and negative reference filtering) and trim length set to 100 bases. Both pipelines showed no significant separation between the BF-bx and BF-bx DC samples (see Fig. S4C and D in the supplemental material).
RNA extraction, rRNA depletion, and RNA sequencing. Total RNA was extracted from frozen proximal colon tissue snips (Fig. 1AII) using the mirVana miRNA isolation kit, with phenol (ThermoFisher Scientific catalog no. AM1560) according to the manufacturer's instructions, with the addition of an ϳ1:1 mix of 1-mm acid-washed glass beads and 0.1-mm zirconia beads and a Precellys24 (Bertin Instruments catalog no. EQ03119-200-RD000.0) bead beater for tissue disruption and lysis. Extracted RNA was treated with the Turbo DNA-free kit (ThermoFisher Scientific catalog no. AM1907) to remove DNA. Quality control, rRNA depletion, and cDNA library preparation were performed by the University of Florida's Interdisciplinary Center for Biotechnology Research (ICBR) Gene Expression and Genotyping core using the Agilent 2100 bioanalyzer (Agilent Genomics catalog no. G2939BA), Ribo-Zero Gold rRNA removal kit (Epidemiology) (Illumina catalog no. MRZE724) and ScriptSeq v2 RNASeq library preparation kit (Illumina catalog no. SSV21124) starting with 1 g total RNA. Samples were sequenced by the University of Florida ICBR NextGen DNA Sequencing core on the Illumina HiSeq 3000 (2 ϫ 100 run), multiplexing each sample into three lanes to avoid lane effect.
Mouse RNA-seq analysis. Reads were quality filtered at Q20 and trimmed to remove remaining adapters using Trimmomatic (79) (87) pathways, and genes were mapped to KEGG pathways using Pathview (88). We considered a pathway significant if its GAGE false-discovery rate (q value) was less than 0.05. We tested the effect of sequencing lane on the clustering of the samples ( Fig. 2A) and found it to be insignificant (P value Ͼ 0.05) (data not shown). Metatranscriptome analysis. Quality-filtered and trimmed reads from above were aligned to iGenome Mus musculus Ensembl GRCm38 reference genome using BWA (89) version 0.7.16a, and reads with alignments were excluded from further analysis. The remaining reads were then filtered from rRNA and tRNA sequences by aligning (using BWA) to a collection of NCBI rRNA and tRNA sequences and SILVA database sequences, resulting in an average of 1,208,429 reads per sample, which were then submitted for de novo assembly using Trinity (90) version 2.4.0. The resulting assembly was annotated using Trinotate (91) version 3.0.1 (http:// trinotate.github.io) with the following databases: uniprot_sprot (92), Pfam (93), and Virulence Factor Database (VFDB) (94). The resulting annotations were examined, and sequences annotated as nonbacterial were removed. Transcript abundance was determined using RNA-seq by expectation maximization (RSEM) (95) through Trinity's align_and_estimate_abundance.pl script, and the counts were imported to edgeR version 3.16.5 for differential expression analysis. A gene was considered for the differential expression test if it was present in at least 50% of the samples. We considered a transcript DE if its edgeR FDR-adjusted P value was Ͻ0.05. To account for normalization artifacts, we also examined the ratios of DE genes between the BFϩT and BF-bx groups generated from a rarefied data set that was based on 140,000 reads per sample ( Table 2). The similar ratios of DE genes from the complete and rarefied data sets suggest that our findings are not an artifact of normalization. We conducted a second analysis (reference-based analysis) by aligning the reads submitted for de novo assembly to the human gut microbiome integrated gene catalog (IGC) (96) using Bowtie2 (81) (v.2.3.4.2) followed by quantification using featureCounts (85) from the subread package (v.1.5.3) and obtained similar results (Fig. S2A).
Pathway analysis was conducted through GAGE version 2.24 using Kyoto Encyclopedia of Genes and Genomes (KEGG) reference pathways on the assembled transcript and The Human Microbiome Project (HMP) Unified Metabolic Analysis Network (HUMAnN) (97) on the unassembled reads. Genes were mapped to KEGG pathways using Pathview (88). We considered a pathway significant if its q value was Ͻ0.05. Weighted gene coexpression network analysis (WGCNA) version 1.68 (23) was utilized to detect modules in each biofilm status samples using the blockwiseConsensusModules function which performs the network construction and consensus module detection. The hub gene in each detected module was identified using the WGCNA function chooseTopHubInEachModule. The sequencing lane had no effect on the clustering of the samples in Fig. 1A (P value Ͼ 0.05, data not shown). miRNA extraction and sequencing. Small RNAs were extracted from snap-frozen stool samples (Fig. 1AII) using the mirVana miRNA isolation kit. Because of the low amount of small RNA, GF stools were pooled from 13 Apc MinΔ850/ϩ ;Il10 Ϫ/Ϫ mice (20-to 44-week age range) total, or stools from 2 to 5 mice per sample (n ϭ 4). BF-bx and BFϩT small RNAs were extracted from the stools of 12-week-associated BF-bx (n ϭ 7) and BFϩT (n ϭ 10) Apc MinΔ850/ϩ ;Il10 Ϫ/Ϫ mice. cDNA libraries were synthesized with the NEBNext Multiplex Small RNA Library Prep Set for Illumina kit (New England Biolabs catalog no. E7300) and small RNAs for each library (21-to 30-nucleotide size range) were gel purified. For the GF, BF-bx, and BFϩT comparisons, a pool of 21 libraries (equivalent molar concentrations; 4 GF, 7 BF-bx, and 10 BFϩT) were multiplexed and sequenced using the Illumina Miseq. miRNA analysis. CAP-miRSeq (98) was used to process the miRNA sequences. We used the databases and reference sequences that ship with CAP-miRSeq for all the analyses. Briefly, sequences were filtered and trimmed using cutadapt (99). Quantification of miRNA was done using miRDeep2 (100), and DE miRNAs were detected using edgeR version 3.16.5. We considered a miRNA DE if its edgeR FDR-adjusted P value was Ͻ0.05. Principal-component analysis (PCA) was created using R's prcomp function from the normalized and log 10 -transformed miRNA counts according to the equation above.
Correlations with microbiota taxon abundance were done using R lm function, P values were determined using R's ANOVA function, and FDR correction was done using R's p.adjust function  (76) and only those with a FDR-adjusted P value of Ͻ0.05 were considered. Mouse miRNA targets were predicted using miRDB (101), and bacterial target prediction for mice miRNA was done using PITA (102) on the assembled bacterial transcripts from the RNA-seq data described above. We considered a bacterial transcript a potential target for a particular mouse miRNA if its ΔΔG score was less than or equal to Ϫ15 kcal/mol. The bacterial and mouse gene targets of miRNAs significantly different between GF, BF-bx, and BFϩT groups are listed in Tables S12 at https://figshare .com/s/33fb65da16763a5e3347 and S13 at https://figshare.com/s/be3f38607b2e9875e7f6.
For miRNA expression correlation with tumor numbers, two BF-bx samples were excluded because they were fixed without splaying the colon so tumor counts could not be generated. Correlation was done using Spearman's rank correlation through R's cor.test function.
Statistical analyses. For all sequencing analyses, a taxon or miRNA was considered only if it was present in at least 30% of the comparison samples, and statistics are described in above 16S rRNA sequencing, mouse RNA-seq, metatranscriptome and miRNA analysis sections. P values of Ͻ0.05 were considered statistically significant.
Ethics. All animal experiments were approved by the University of Florida Institutional Animal Care and Use Committee (protocol 201308038). All patient tissues were collected as previously described (11,12,103).
Data availability. The data supporting the results of this article have been deposited in the National Center for Biotechnology Information Gene Expression Omnibus (NCBI GEO) under accession number GSE108165 for the RNA-seq and miRNA data sets and the National Center for Biotechnology Information Sequence Read Archive (NCBI SRA) under BioProject identifier (ID) PRJNA422588 for 16S rRNA Run02 sequences. Run01 sequences are in the article by Tomkovich et al. (13).

SUPPLEMENTAL MATERIAL
Supplemental material is available online only.