Role of Fusobacteria in the serrated pathway of colorectal carcinogenesis

Fusobacteria are associated with colorectal cancer (CRC) and are amplified during colorectal carcinogenesis. Compared to the adenoma-carcinoma sequence of carcinogenesis, serrated neoplasm has distinct clinical features and a different molecular background. We aimed to compare the gut microbiome between tubular adenoma (TA) and sessile serrated adenoma/polyp (SSA/P). Patients with TA, SSA/P, or CRC were recruited. Three pieces of colorectal mucosal tissue were obtained from each patient by endoscopic biopsy. 16S rRNA gene pyrosequencing and phylogenetic investigation of communities by reconstruction of unobserved states (PICRUSt) were performed. Among 26 enrolled patients, 8, 10, and 8 had TA, SSA/P, and CRC, respectively. The relative abundance of Fusobacteria did not differ significantly between the TA and SSA/P groups (4.3% and 1.9%, P = 0.739) but was higher in the CRC group (33.8%) than in the TA or SSA/P group, respectively (TA vs. CRC, P = 0.002, false discovery rate [FDR] = 0.023; SSA/P vs. CRC, P < 0.001, FDR = 0.001). PICRUSt revealed that most functions in the TA metagenome were similar to those in the SSA/P metagenome. The gut microbiome, including relative abundance of Fusobacteria, did not differ between TA and SSA/P, suggesting that Fusobacteria may contribute to both the serrated pathway and the adenoma-carcinoma sequence.


Results
Baseline characteristics and microbiome reads. Table 1 shows baseline patient characteristics and microbiome reads. Among 26 enrolled patients, 8 had CRC, 8 had TA, and 10 had SSA/P. Median age was 68, 54, and 63 years in the TA, SSA/P, and CRC group, respectively (P = 0.119). The proportion of males did not differ significantly among the groups (P = 0.244). Five (62.5%) TAs and 7 (70.0%) SSA/Ps were in the proximal colon, whereas all lesions in the CRC group were in the distal colon.
A total of 136,414 read were obtained from 26 enrolled patients through 16S rRNA gene pyrosequencing. The median number of reads was 5,540, 6,162 and 4,969 in TA, SSA/P, and CRC groups, respectively (P = 0.247). Median operational taxonomic units (OTUs), Chao1 estimator, and Shannon's and Simpson's diversity index did not differ among the groups.
Relative bacterial abundance. Figure 1A shows relative bacterial abundance at the phylum level for each group. The five most abundant phyla were Fusobacteria, Firmicutes, Proteobacteria, Bacteroidetes, and Actinobacteria in all three groups. Relative abundance of Fusobacteria did not differ between the TA and SSA/P group (4.3% and 1.9%, P = 0.739), but was significantly higher in the CRC group (33.8%) than in the TA or SSA/P group, respectively (Table 2; TA vs. CRC, P = 0.002, FDR = 0.023; SSA/P vs. CRC, P < 0.001, FDR = 0.001).
Because all lesions of the CRC group were in the distal colon, we performed a subgroup analysis for lesions located in the distal colon. As shown in Fig. 1B, relative abundance of Fusobacteria in the CRC group was higher  than that in the TA or SSA/P group. Bacterial abundance at the phylum level was not associated with tumor location.
Relative abundance of Fusobacteria in each tissue. Figure 2 shows relative abundance of Fusobacteria in each tissue. Although the relative abundance of Fusobacteria varied, 37.5% (3 of 8) and 50.0% (5 of 10) of lesions contained Fusobacteria in the TA and SSA/P group, respectively. In the CRC group, Fusobacteria was identified in all lesions. Interestingly, Fusobacteria was identified not only in the lesions but also in adjacent normal tissues. Four (50.0%), three (30.0%), and eight (100.0%) adjacent normal tissues in the TA, SSA/P, and CRC group, respectively, showed the presence of Fusobacteria.
The relative abundance of functional category according to the level 3 KEGG module is illustrated in Fig. 4. Overall, distribution of each predictive function at the level 3 KEGG module did not differ among the three groups. However, several functions of the cellular process category including bacterial chemotaxis, bacterial motility proteins, and flagella assembly tended to be highly expressed in the TA and SSA/P metagenomes compared to the CRC metagenomes. On the contrary, several functions of the genetic information processing category including RNA degradation, DNA repair and recombination proteins, and mismatch repair tended to be more highly expressed in the CRC metagenome than in the TA and SSA/P metagenomes. Detailed relative abundance of predicted functions for specific level 3 KEGG modules is shown in supplementary Dataset 1.

Discussion
Many studies have revealed that Fusobacteria are prevalent in CRCs compared to adjacent normal tissues [6][7][8] . It has been shown that Fusobacteria are present in colorectal premalignant lesions including TA and SSA/P although the relative abundance is lower than in CRC tissues 13 . Moreover, it was shown that abundance of F. nucleatum in CRC tissue was associated with a lower density of CD3+ T cell 14 . It has been suggested that Gut bacteria including Fusobacteria may have an influence on the development of CRC through interaction with the innate immune system or host factors 3 . These evidences support that Fusobacteria may have a distinct role in the carcinogenesis of the colorectum.
Recently, Kostic et al. demonstrated that F. nucleatum can induce intestinal tumorigenesis in APC Min/+ mice, and showed that Fusobacteria can affect the early stage of the adenoma-carcinoma sequence such as APC-mutated adenoma. However, the influence of Fusobacteria on CRC development through the serrated pathway remained relatively unknown. Compared to the adenoma-carcinoma sequence, sessile serrated adenoma/polyp has an obviously different pathway to the development of CRC 1,9,15,16 . For example, abnormal hypermethylation of promoter CpG island is common in the serrated pathway to CRC 9 . In addition, BRAF mutations are strongly associated with CpG island methylator phenotype (CIMP)-high CRCs, acting as an alternative to the KRAS mutations that commonly occur in chromosomal instability cancers 9 . Moreover, APC mutation is relatively rare in SSA/Ps 15 . We hypothesized that if Fusobacteria are not associated with SSA/P, the gut microbiome, including relative abundance of Fusobacteria, might differ between TA and SSA/P tissues. In our study, however, the relative abundance  of Fusobacteria did not differ between the TA and SSA/P groups. PCoA also showed that the gut microbiome in the TAs was more similar to that of the SSA/Ps than the CRCs. These results suggested that Fusobacteria is a possible contributor to both the adenoma-carcinoma sequence and the serrated polyp-carcinoma sequence.
In this study, the relative abundance of Fusobacteria in the TA or SSA/P groups was lower than in the CRC group. Considering that Fusobacteria are extremely rare in healthy individuals 17 , their abundance in 37.5-50% of lesions may be higher than the relative abundance of Fusobacteria in normal colorectal mucosal of healthy individuals. These findings were consistent with previous studies. One recent study by Ito et al. using quantitative polymerase chain reaction for F. nucleatum in formalin-fixed paraffin-embedded tumor tissues showed that only 30-40% of serrated adenomas expressed positivity for F. nucleatum 13 . Additionally, F. nucleatum expression did not differ between serrated and non-serrated adenomas. Moreover, F. nucleatum was more expressed in CRC tissues than in both serrated and non-serrated adenomas. Based on their results, Ito et al. concluded that expression of F. nucleatum may contribute to the progression of colorectal neoplasia.
Differences in the gut microbiome, including the relative abundance of Fusobacteria, between TA or SSA/P and CRC can be explained by the physiologic and metabolic changes that occur during colon carcinogenesis, including changes in colonic barrier function and rupture and bleeding of the colonic epithelium 18 . For example, Bacteroides fragilis (ETBF) strains produce B. fragilis metalloprotease toxin 19,20 . This metalloprotease facilitates cleavage of E-cadherin, that is tumor suppressor protein, in intestinal epitheliums, causing cell proliferation and permeabilization of the intestinal barrier 21 . As a result of changes in colonic barrier permeability and cellular metabolism, Fusobacteria might be more prevalent in CRC tissues than in adenomas 18 .
One interesting finding of our study is that Fusobacteria was identified in the adjacent normal tissues as well as colorectal neoplasms. The relative abundance of Fusobacteria between lesions and adjacent normal tissues was similar in each group. This finding can be explained by the field cancerization effect, i.e., a field of cellular and molecular alterations that predisposes to the development of neoplasms within that territory 22,23 . Other studies on gut microbiomes in patients with CRC also showed that Fusobacteria are abundant in adjacent normal tissues as well as CRC tissues 24 .   An additional merit of our study is functional analysis for the predicted metagenome using PICRUSt in TAs, SSA/Ps, and CRCs. Although statistical power was not achieved because of the small sample size, a tendency for functional differences among the three metagenomes could be identified. Overall, the relative abundance of the cellular processes category including bacterial chemotaxis, bacterial motility proteins, and flagella assembly tended to be higher in the TA and SSA/P groups, whereas that of the genetic information processing category including RNA degradation, DNA repair and recombination proteins, and mismatch repair tended to be higher in the CRC group. These results, however, were metagenomic inferences without statistical significance. Therefore, more experimental studies on these bacterial functions should be required to clarify the role of the gut microbiome in CRC development.
Although this study was the first to compare the gut microbiome between two types of CRC precursor, TA and SSA/P, it has several limitations. The small sample size is the first limitation. Although some significant differences in the gut microbiome were identified at the phylum level, statistical power could not be achieved at a deeper level because of the small sample size. Further large-scaled studies are needed for identifying detailed differences of gut microbiome among disease status. Second, all CRCs in our study were located in the distal colon. Unfortunately, patients with proximal colon cancer could not be recruited during the study period therefore we performed a subgroup analysis for lesions located in the distal colon. Subgroup analysis revealed that the relative abundance of bacteria including Fusobacteria did not differ according to the tumor location in our study. Finally, although our data showed the comparison of gut microbiome between patients with TA and those with SSA/P, we could not make a definitive conclusion about the role of Fusobacteria in serrated pathways. Experimental studies using a SSA/P animal model may provide data supporting our current study.
Despite these limitations, our data provide a better understanding of the role of Fusobacteria in the development of CRC. The relative abundance of Fusobacteria in CRCs was higher than that in the TAs or SSA/Ps. In addition, the gut microbiome, including relative abundance of Fusobacteria, did not differ between the TAs and SSA/Ps. Therefore, Fusobacteria may contribute to the serrated pathway as well as the adenoma-carcinoma sequence.

Methods
Study population. Patients that were pathologically diagnosed with TA, SSA/P, or CRC were recruited for this study. Written informed consent was obtained from all patients. Patients who were diagnosed with inflammatory bowel disease and human immunodeficiency virus enteropathy were excluded. Additionally, patients who had been treated with antibiotics or antidiarrheal agents within 3 months prior to enrollment were excluded. The study protocol was approved by the Institutional Review Board of Hanyang University Guri Hospital (GURI 2014-10-010). All experiments were performed in accordance with relevant guidelines and regulations.
Tissue sampling and bacterial DNA extraction. Three pieces of mucosal tissue were obtained from patients with TA, SSA/P, or CRC by endoscopic biopsy. Additionally, four biopsy samples were taken from adjacent normal tissues of each patient.
Extraction of bacterial DNA was performed from mucosal biopsy samples as we previously reported 25 . Briefly, 100 mg of frozen gastric mucosal tissues was suspended in 750 μ L of sterile bacterial lysis buffer (200 mmol/L NaCl, 100 mmol/L EDTA [pH 8.0], 20 mmol/L Tris base, 20 mg/mL lysozyme) and incubated at 37 °C for 30 minutes. Then, we added 20 μ L of proteinase K and 80 μ L of 10% SDS to the mixture and incubated it at 65 °C for 30 minutes. Finally, bead beating was performed for 90 seconds at 5,300 rpm (PRECELLYS 24; Bertin Technologies, Le Bretonneux, France) after adding a 300 mg of 0.1-mm zirconium beads (BioSpec Products, Bartlesville, OK, USA) for finishing homogenization. The homogenized mixture was cooled on ice and then centrifuged at 14,000 rpm for 5 minutes. Bacterial DNA was extracted from the supernatant by phenol/chloroform/iso-amyl alcohol (25:24:1) followed by chloroform/iso-amyl alcohol (24:1), and precipitated by absolute ethanol at − 20 °C for 1 hour. The precipitated DNA was suspended in DNase-free H 2 O and cleaned up using a DNA clean-up kit according to the manufacturer's instructions. Isolated DNA was stored at − 80 °C until use in microbial characterization Amplification of 16S rRNA gene and sequencing. Extracted gDNA was amplified using primers targeting the V1 to V3 regions of the 16S rRNA gene. For bacterial gDNA amplification, barcoded primers of 27F 5′ -CCATCTCATCCCTGCGTGTCTCCGACTCAG-MID-GAGTTTGATCMTGGCTCAG-3′ , which was consisted of the GS FLX+ adapter sequence A (underlined sequence), linker nucleotides (bolded sequence), multiplex identifiers sequence (MID), and the universal 16S rRNA-specific sequence (italic sequence). The reverse primer 518R 5′ -CCTATCCCCTGTGTGCCTTGGCAGTCTCAG-MID-WTTACCGCGGCTGCTGG-3′ consisted of the GS FLX+ adapter sequence B (underlined sequence), linker nucleotides (bolded sequence), MID, and the universal 16S rRNA-specific sequence (italic sequence). The amplifications were performed at the following conditions: initial denaturation at 95 °C for 5 minutes, followed by 35 cycles of denaturation at 95 °C for 40 seconds, primer annealing at 57 °C for 40 seconds, and extension at 72 °C for 60 seconds; followed by a final elongation at 72 °C for 60 seconds. 16S rRNA PCR products were quantified, pooled, and purified for the sequencing reaction. Sequencing was performed using a 454 Life Sciences Genome Sequencer FLX+ machine (Roche, Florence, SC, USA) according to the manufacturer's instructions.
Analysis of 16S rRNA Sequences. Initially generated partial reads of 16S rRNA were trimmed for quality using the standard SFF software tools from Roche/454. Pre-processing and clustering cluding raw read filtering and trimming, reads picking, and OTU clustering using the CD-HIT-OTU 26 . In addition, bacterial 16S rRNA sequence data of the gastric mucosal microbiome were processed through the QIIME pipeline 27 . Reads which showed low quality, contained incorrect primer sequences, or contained more than one ambiguous base were excluded. The remaining reads were classified into CRC, tubular adenoma, and serrated adenoma groups based on their unique nucleotide barcodes. Chimeric sequences were filtered out using the UCHIME algorithm after nearest alignment space termination (NAST) based on the SILVA database 28,29 . Taxonomic composition from phylum to species levels and bacterial diversity for each sample were evaluated based on 97% similarity using the QIIME 30 . A UniFrac distance was calculated from representative reads of each OTU using the UniFrac program 31 . Additionally, Unifrac and weighted Unifrac distances were fed into a principal coordinate analysis (PCoA) 31 .
DNA sequences obtained from this metagenomic project have been deposited in the NCBI Short Read Archive under the Accession No. SRA306790.
Metagenome prediction. Functional composition of the microbiome was predicted using phylogenetic investigation of communities by reconstruction of unobserved states (PICRUSt) 12 . PICRUSt uses evolutionary modeling to predict metagenomes from 16S data compared with a reference genome database 12 . Metagenome prediction in this study was conducted using Galaxy, which is an open, web-based platform for computational biomedical research (http://huttenhower.sph.harvard.edu/galaxy) 12 . To evaluate dissimilarity between reference genomes and the metagenome, the nearest sequenced taxon index (NSTI), which is the sum of phylogenetic distances for each organism in the OTU table to its nearest relative with a sequenced reference genome, was calculated 12 . NSTI has a negative correlation with accuracy of PICRUSt. Functions of predicted metagenomes were classified according to the Kyoto Encyclopedia of Genes and Genomes (KEGG).
Statistical analysis. Baseline patient characteristics and results of 16S rRNA gene pyrosequencing were described as median with interquartile range (IQR) or number with proportion. Fisher's exact test or Kruskal Wallis test was used for group comparisons. Relative abundance of different bacteria phyla in each group was presented as mean with standard deviation. Assigned OTUs were compared individually between study groups using Mann-Whitney U test with false discovery rate (FDR) control to correct for multiple comparisons 32 . All statistical procedures were conducted using R (version 2.15.3; R Foundation for Statistical Computing, Vienna, Austria).