Tumor microbiome diversity influences papillary thyroid cancer invasion | Panda Anku

Tumor microbiome communities are significantly associated with tumor invasion in patients with resected PTC

The demographic and clinical characteristics of the patients with PTC are shown in Supplementary Tables S1 and S2. The tumor microbial diversity was compared among surgically resected patients in different clinical stages. The tumor microbial diversity was measured using different methodologies (Sobs, Shannon, and Simpson indices). The alpha-diversity (α-diversity) of the tumor microbiome was significantly lower in patients with T1/T2 PTC than in those with T3/T4 PTC, as shown by the Shannon and Simpson indices (p = 0.0309 and p = 0.0088; Wilcoxon rank-sum test was also performed Fig. 1a). This indicated lower microbiome diversity (Shannon and Simpson indices) in patients with T1/T2 PTC. Microbiome richness was measured by the number of observed OTUs (Sobs index), and no significant differences were found among different clinical stages (Fig. 1b). To gain a better understanding of the role of microbiome diversity, beta-diversity (β-diversity) was used to carry out a principal coordinate analysis (PCoA) using Bray–Curtis metric distances (Fig. 1c; p = 0.001 [T1/T2/T3/T4]; p = 0.001 [T1_2/T3_4]). Significant differences in β-diversity were observed, further suggesting that the tumor microbial communities varied during tumor progression.

Fig. 1: Changes in the intratumor microbiomes of patients with papillary thyroid carcinoma (PTC) at different clinical stages.
figure 1

a Alpha-diversity (α-diversity) at the operational taxonomic unit (OTU) level (estimated using Sobs, Shannon, and Simpson indices) in patients with PTC in groups T1_2 and T3_4. b α-Diversity at the OTU level (estimated using Sobs, Shannon, and Simpson indices) in patients with PTC of clinical stages T1, T2, T3, and T4 (*p < 0.05 and **p < 0.01). Note: The Sobs values indicate the number of observed OTUs, Shannon diversity index, and Simpson diversity index for an OTU definition, respectively. The lower the value of the Simpson index, the higher the diversity (opposite of the Shannon index). Error bars indicate SD. c Principal coordinate analysis (PCoA) score plots based on the Bray–Curtis distance at the OTU level in patients with PTC in groups T1_2 and T3_4. d PCoA score plots based on the Bray–Curtis distance at the OTU level in patients with PTC at clinical stages T1, T2, T3, and T4. Wilcoxon rank-sum test and analysis of similarities (ANOSIM) were used for the analysis of intragroup differences in α-diversity and PCoA analysis, respectively.

The general landscape of the tumor microbiome composition was assessed by considering the relationship between PTC intratumoral bacterial diversity and clinical stages (Fig. 2). Enrichment of specific bacterial communities was evaluated at different taxonomic levels (Supplementary Fig. S1). At the genus level, Pseudomonas was the dominant bacterium, followed by Rhodococcus, Ralstonia, Acinetobacter, and Sphingomonas. Next, we sought to determine whether tumor microbiome composition differed among patients with PTC at stages T1, T2, T3, and T4. We found similar communities among patients with PTC at various stages (Fig. 2b, c, Supplementary Data 1–3, 8). However, considerable differences in the percentage of community abundance were also observed via Circos analysis (Fig. 2d, e), showing the percentage of the dominant genera in distinct groups. For example, the percentage of Pseudomonas was higher in T1 or T2 than in T3 or T4.

Fig. 2: Tumor microbiome communities significantly differ among clinical stages.
figure 2

ac Bar plots of the genus level in patients with PTC. a Relative abundance for all patients with PTC. The genus level shows a more refined taxonomic level than the phylum level in bacteria and can help find particular species. b Relative abundance for patients with PTC in groups T1_2 and T3_4. c Relative abundance for patients with PTC in clinical stages T1–T4. d and e Distribution of microbial communities in each group at the genus level. d Distribution of microbial communities for patients with PTC in groups T1_2 and T3_4. e Distribution of microbial communities for patients with PTC at clinical stages T1–T4. The data are visualized using Circos. The width of the bars for each genus indicates the relative abundance of that genus in the group. f Microbiome alterations at the genus level in patients with PTC of clinical stages T1–T4. Intragroup differences were analyzed using the Kruskal–Wallis test, and differences between groups were analyzed based on the post-hoc test using Welch’s uncorrected test, adjusted by false discovery rate. The top three differential bacteria (genus) identified were tested individually. g Taxonomic cladogram from LEfSe, depicting taxonomic associations between microbiome communities from patients with T1_2 and T3_4 PTC. Each node represents a specific taxonomic type. Yellow nodes denote taxonomic features that are not significantly differentiated between T1_2 and T3_4. Red nodes denote taxonomic types more abundant in T1_2 than in T3_4, whereas blue nodes represent taxonomic types more abundant in T3_4. h LDA score computed from differentially abundant features between T1_2 and T3_4. The criterion for feature selection was log(LDA score) >3.2. i Receiver operating characteristic (ROC) curve of the relative abundance of taxa as a predictor of clinical T stage status. The eight differential bacteria (genus) were tested in aggregate as average. The confidence interval threshold was set as 0.95 (p = 0.00003).

To further investigate the specific changes in the microbiome in the tumors of patients with PTC at different stages, the relative abundance of taxa was assessed. At the genus level, three genera, namely Pseudomonas (p = 0.0017), Rhodococcus (p = 0.02644), and Sphingomonas (p = 0.0073), displayed a difference in abundance among various stages (Fig. 2f and Supplementary Fig. S2a). Pseudomonas spp., the most abundant genus in all groups, presented a higher abundance in tumors of patients with T1 and T2 PTC than in those with T3 or T4 PTC (p = 0.0049, p = 0.0138, p = 0.010, and p = 0.0028, respectively; Fig. 2f). Rhodococcus abundance was also significantly higher in patients with T1 PTC than in those with T3 PTC (p = 0.0032; Fig. 2f), and Sphingomonas was more abundant in T1 and T2 than in T3 (p = 0.0001 and p = 0.0005, respectively; Fig. 2f, Supplementary Data 4). To further investigate these findings, we conducted high-dimensional comparisons using linear discriminant analysis of effect size (LEfSe). We detected considerable differences in the predominance of bacterial communities between T1_2 and T3_4 (Fig. 2g, h Supplementary Data 5). The T1_2 tumors exhibited a predominance of Pseudomonas (Pseudomonadales at the order level), Rhodococcus (Corynebacteriales at the order level), and Sphingomonas (Sphingomonadales at the order level) at the genus level. In contrast, T3_4 tumors were dominated by Streptococcus, Granulicatella, Haemophilus g_unclassified_o_Rhizobiales, and g_norank_f_norank_o_Coriobacteriales at the genus level (Fig. 2g, h). We then used the eight genera for the area under the curve (AUC)-receiver operating characteristic (ROC) analysis. The combination of the eight taxa (Pseudomonas, Rhodococcus, and Sphingomonas Streptococcus, Granulicatella, Haemophilus g_unclassified_o_Rhizobiales, and g_norank_f_norank_o_Coriobacteriales) resulted in an AUC of 0.83 in the T1_2 and T3_4 groups (Fig. 2i), and 0.91 in the T1 and T4 groups (Supplementary Fig. S2b), which showed higher discriminative capacity than in the group of bacteria filtered by random forest analysis (Supplementary Fig. S3). In addition, we conducted an AUC-ROC analysis with a separate set of samples using the combination of the eight taxa, which resulted in an AUC of 0.79 for TT1_2 and T3_4 groups (Supplementary Fig. S5b), although the composition of bacterial communities was evaluated at different taxonomic levels (Supplementary Fig. S5a). This confirmed the potential of the eight-genera microbiome signature to discriminate and influence PTC invasion status.

Intratumor microbial dysbiosis is related to thyroid function

Identified functional predictions differentially presented metaCyc pathways between patients with T1_2 and T3_4 PTC. The top pathways that differed between these two groups were listed according to the effect size (difference between proportions >0.2 and the ratio of proportions >10; p < 0.05; Supplementary Fig. S2c). The T1_2 group had higher proportions of most of the pathways, including fatty acid salvage, sulfate reduction (assimilatory), fatty acid and beta-oxidation I, octane oxidation, and l-tyrosine degradation. While the T1_2 group had higher proportions of super pathways of purine nucleotide de novo biosynthesis and palmitate biosynthesis II (bacteria and plants).

Furthermore, the association among the thyroid-related hormones free T4 (FT4), T4, FT3, T3, and thyroid-stimulating hormone (TSH) and microbial abundance were assessed. The levels of different thyroid hormones were related to different microbial genera (Fig. 3a and Supplementary Fig. S4a, Supplementary Data 6). A positive relationship was observed between FT4 and Neisseria and norank_f__norank_o__Chloroplast, FT3, and Treponema, whereas a negative relationship was observed between FT4 and Klebsiella, T4 and Klebsiella and EscherichiaShigella, T3 and Granulicatella, and TSH and norank_f__norank_o__Clostridia_UCG-014 and Prevotella. In addition, positive correlations between FT4 and norank_f__norank_o__Clostridia_UCG-014 and negative correlations between hormones T3 and FT3 and Granulicatella were found in the other set of samples (Supplementary Fig. S5c). However, all these microbes showed a relatively lower intratumor richness in PTC.

Fig. 3: Heatmap of Spearman’s correlation analysis between the specific PTC intratumor microbiome and the clinical factors.
figure 3

a Heatmap of Spearman’s correlation analysis between the specific bacterial and thyroid-related hormones. *p < 0.05 and **p < 0.01. Heatmap color keys indicate Spearman correlation coefficient. b Heatmap of Spearman’s correlation analysis between the specific bacterial and thyroid diseases (AITD)-related antibodies. *p < 0.05 and **p < 0.01. Heatmap color keys indicate Spearman correlation coefficient.

The association between the intratumor microbiome and autoimmune thyroid disease (AITD)-related antibodies was also investigated, as AITD affects the entire metabolism in the human body, and the immune status may also contribute to PTC progression. Autoantibodies against thyroid peroxidase (TPO) and thyroglobulin (TG), which characterize Hashimoto’s thyroiditis (HT), and thyroid-stimulating receptors (TSHR), a marker for Graves’ disease (GD), were used in this study. The relationship between microbial abundance and the levels of the three autoantibodies was assessed (Fig. 3b, Supplementary Data 7 and Supplementary Fig. S4b). The levels of anti-TSHR, anti-TG, and anti-TPO antibodies were adjusted; values <0.3 were set as 0.2, values <10 were set as 9, and values <5 were set as 4. A negative correlation was found between anti-TSHR levels and Klebsiella and BurkholderiaCaballeroniaParaburkholderia. Meanwhile, the anti-TG levels positively correlated with Sphingomonas, Rhodococcus, Ralstonia, and Brevundimonas but negatively correlated with Anaerococcus and Akkermansia. Nine genera, namely UCG-002, Streptococcus, Parvimonas, Akkermansia, Bacteroides, Haemophilus, Selenomonas, Prevotella, and Bifidobacterium, exhibited a negative relationship with the anti-TPO levels.

Correlation analysis between the abundance of these microbes and AITD-related antibodies was also conducted in the other set of samples (Supplementary Fig. S5c). Bifidobacterium and Haemophilus abundance showed a positive relationship with the anti-TPO levels, while Leptotrichia abundance showed a positive relationship with the anti-TSHR levels. This further proved the correlation between the abundance of microbes and thyroid-related hormones and AITD-related antibodies, while the interaction may differ according to individual conditions. Functional predictions identified different pathways for these specific microbes between patients with T1_2 and T3_4 PTC (Supplementary Fig. S4c). These metabolite differences may disrupt thyroid hormones and AITD-related antibody levels related to thyroid function.

Associations between clinical variables and the intratumor microbiome

An association between sex and intratumor microbiome diversity was found in patients with PTC (Fig. 4). A higher microbiome diversity (alpha-diversity (α-diversity)) was observed in females (Fig. 4a), who are expected to have a higher PTC incidence. However, no obvious difference in beta diversity was found (p = 0.3000). Of note, Rhodococcus, Ralstonia, Chryseobacterium, and Burkholderia-Caballeronia-Paraburkholderia were more abundant in female patients with PTC than in male patients (uncorrected p = 0.0413, p = 0.0092, p = 0.0275, and p = 0.0008, respectively, Wilcoxon rank-sum test; Fig. 4c). Remarkably, the sex-associated genus was the same T1/T2 PTC-enriched genus Rhodococcus. Moreover, a higher PTC microbiome diversity might have a worse association with PTC. We also noted correlations between intratumor bacteria or their predicted functions in patients with different tumor subtypes, ages, or lymphatic metastasis status. However, there was no significant difference with respect to the α- and β-diversities (Supplementary Fig. S6).

Fig. 4: Changes in the intratumor microbiomes of patients with PTC in different sexes.
figure 4

a α-Diversity at the OTU level (estimated using Sobs, Shannon, and Simpson indices) in patients of different sexes with PTC (*p < 0.05 and **p < 0.01). Error bars indicate SD. b PCoA score plots based on the Bray–Curtis distance at the OTU level in patients of different sexes with PTC. Wilcoxon rank-sum test and ANOSIM were used for the analysis of intragroup differences in α-diversity, Simpson index, and PCoA analysis, respectively. c Microbiome alterations at the genus level in patients of different sexes with PTC. Differences between groups were analyzed using the Wilcoxon rank-sum test, and the 95% confidence interval (CI) was calculated using the bootstrap method. G1: male; G2: female.

Leave a Comment