Continuous pre- and post-transplant exposure to a disease-associated gut microbiome promotes hyper-acute graft-versus-host disease in wild-type mice

ABSTRACT

Objective

The gut microbiome plays a key role in the development of acute graft-versus-host disease (GVHD) following allogeneic hematopoietic stem cell transplantation. Here we investigate the individual contribution of the pre- and post-transplant gut microbiome to acute GVHD using a well-studied mouse model.

Design

Wild-type mice were cohoused with IL-17RA^{–/ –} mice, susceptible to hyperacute GVHD, either pre- or post-transplant alone or continuously (i.e., pre- and post-transplant). Fecal samples were collected from both WT and IL-17RA^{–/ –} mice pre- and post-cohousing and post-transplant and the microbiome analyzed using metagenomic sequencing.

Results

Priming wild-type mice via cohousing pre-transplant only is insufficient to accelerate GVHD, however, accelerated disease is observed in WT mice cohoused post-transplant only. When mice are cohoused continuously, the effect of priming and exacerbation is additive, resulting in a greater acceleration of disease in WT mice beyond that seen with cohousing post-transplant only. Metagenomic analysis of the microbiome revealed pre-transplant cohousing is associated with the transfer of specific species within two as-yet-uncultured genera of the bacterial family Muribaculaceae; CAG-485 and CAG-873. Post-transplant, we observed GVHD-associated blooms of Enterobacteriaceae members Escherichia coli and Enterobacter hormaechei subsp. steigerwaltii, and hyperacute GVHD gut microbiome distinct from that associated with delayed-onset disease (>10 days post-transplant).

Conclusion

These results clarify the importance of the peri-transplant microbiome in the susceptibility to acute GVHD post-transplant and demonstrate the species-specific nature of this association.

KEYWORDS:

• GVHD
• hematopoietic stem cell transplant
• S24-7
• Muribaculaceae
• cohousing
• microbiome
• Enterobacteriaceae

View correction statement:

Correction

Introduction

Graft-versus-host disease (GVHD) is a serious complication of hematopoietic stem cell transplantation occurring in 30–50% of cases.¹ The condition is characterized by T cell-mediated tissue damage to target organs, which include the skin, liver and gastrointestinal (GI) tract. In fact, in T-cell replete-unrelated stem cell transplants (SCT), the GI tract is involved in virtually all fatal cases of acute GVHD, with an overall survival at 2 years after the onset of stage 3–4 gut GVHD of 25%.² The gut microbiome was initially implicated as a key contributor to the development of intestinal GVHD following experiments demonstrating reduced disease severity in germ-free and antibiotic-treated mice.^3,4 Conditioning regimens administered prior to transplant to ablate the immune system are necessary for successful engraftment, however, tissue injury resulting from these treatments permits translocation of microbes and microbial products (e.g., lipopolysaccharide) from the gastrointestinal lumen into circulation, triggering the release of inflammatory cytokines TNF, IL-1, and IL-6.^5–8 This systemic circulation of danger/pathogen-associated molecular patterns and inflammatory cytokines results in the activation of host antigen-presenting cells, subsequent priming, differentiation and expansion of donor T-cells, and an overwhelming immune response resulting in destruction of target tissue.⁹ Antibiotic treatment aimed at circumventing this process forms part of some conditioning regimens, however, concerns regarding antimicrobial resistance means the practice is not universal.¹⁰ Disruption of the commensal microbiome, a system critical to human and animal health, may also impose further stress on the transplant recipient potentially impacting recovery.¹¹ Interestingly, we recently demonstrated a critical role for intestinal microbiota in influencing MHC class II presentation by intestinal epithelial cells leading to the initiation of GVHD.¹² As such, there exists a need to increase our understanding of the precise role of the gut microbiome in regulating the development of intestinal GVHD with a view to designing better targeted prophylactic or treatment strategies.

Investigation into the role of the gut microbiome in GVHD is complicated by observed microbial community fluctuations occurring both prior to and concomitant with the development of GVHD. At the time of disease, bacterial diversity within the gut decreases in both mice and humans in concert with large compositional shifts.¹³ Members of the bacterial families Enterobacteriaceae, Lactobacillaceae and Enterococcaceae have all been shown to increase in abundance with disease development, often dominating the community as a whole.^13–16 Decreased diversity in patients at the time of transplant, but prior to GVHD initiation, is associated with increased disease mortality,¹⁷ in particular, decreased abundance of Blautia species has been linked to GVHD in humans.¹⁸ Acute GVHD development is associated with increased abundance of Firmicutes in patients pre-transplant,^19,20 as is an overall decreased level of diversity,¹⁹ suggesting that the composition of the gut microbiota prior to transplant may also play a role. Understanding whether there is a defining period during which the gut microbiome is critical in the development of GVHD is vital to fully elucidating the nature of the relationship.

We recently observed an association between fecal microbiota and the development of hyper-acute GVHD in a mouse model of GVHD.²¹ Following the description of a distinct microbial community within IL-17RA^{–/ –} mice, a mouse genotype in which hyper-acute GVHD develops consistently post-transplant (survival post-transplant ≤10 days), it was discovered that accelerated disease could be induced in WT mice via cohousing with IL-17RA^{–/ –} mice.²¹ Analysis of the fecal microbiome before and after cohousing revealed transfer of community members primarily from IL-17RA^{–/ –} to WT mice. However, as mice were cohoused throughout the experiment, it was unclear whether the pre-transplant cohousing was sufficient to induce accelerated disease or whether cohousing post-transplant also contributed to this observation. In order to ascertain whether there was an essential cohousing period, and by extension, an essential transfer of microorganisms, we undertook additional experiments using mice cohoused either pre- or post-transplant alone. These data demonstrate that sustained exposure to the disease-associated microbiota from IL-17 receptor-deficient mice is necessary for maximum induction of hyperacute GVHD in WT mice. Cohousing post-transplant only results in a rate of post-transplant disease progression intermediate between non-cohoused WT controls and continuously cohoused WT mice. Cohousing pre-transplant only is insufficient to significantly accelerate disease. The detrimental effect of pre-transplant priming is therefore only realized with post-transplant exposure to IL-17RA^{–/ –} mice.

Results

Cohousing with a dysbiotic microbiome plays both a priming and exacerbating role in the development of acute GVHD in WT mice

In order to determine whether there was a critical period during which exposure to microbiota from IL-17 receptor-deficient mice was able to induce accelerated GVHD in WT mice, we compared the outcome of transplants involving mice cohoused continuously (I), pre-transplant only (II) or post-transplant only (III) (Figure 1a). Note that we previously published survival curves and a 16S rRNA-based analysis of Experiment I,²¹ and in the present study used the DNA from the published experiment for metagenomic analysis. We hypothesized that pre-transplant cohousing primes the WT gut community for pathogen displacement and disease, and that post-transplant cohousing exacerbates disease progression by exposure to high pathogen load. Cohousing mice pre-transplant only (II) resulted in a reduced overall survival time, however, it was not significantly different to separately housed WT mice (median survival 33 versus 45 days, respectively; P = .3605) (Figure 1b). Cohousing mice post-transplant only (III) was more detrimental, with all WT mice succumbing to disease and a significantly accelerated rate of disease overall in comparison to separately housed WT mice (median survival 28 versus 45 days; P = .0015) (Figure 1b). However, the disease exacerbation in WT mice associated with post-transplant cohousing (III) was significantly delayed in comparison to that experienced by WT mice that were cohoused continuously (median survival 28 versus 9 days; P = .0057) (Figure 1b). While continuously cohoused WT mice (I) began to succumb to disease 4-days post-transplant, WT mice cohoused post-transplant but not pre-transplant survived for a minimum of 23 days. We suggest that this lag in disease initiation can be attributed to the absence of pre-transplant microbiome-dependent priming of the WT mice. Therefore, while the priming effect of cohousing pre-transplant is insufficient to accelerate disease in isolation, when combined with post-transplant exacerbation, primed mice are significantly more susceptible to accelerated GVHD than WT mice cohoused post-transplant only.

Figure 1. Continuous cohousing accelerates GVHD in WT mice post-transplant.

To obtain species-level resolution of transferred microbiota in addition to their functional potential we undertook metagenomic sequencing of fecal samples from one experimental replicate of each cohousing scenario, including the continuously cohoused mice analyzed previously by 16S rRNA amplicon sequencing.²¹ Samples from the same mice were sequenced pre-transplant and at the time of sacrifice post-transplant, as well as pre-cohousing, where applicable. Separately housed (control) mice from each sequenced replicate were also analyzed pre- and post-transplant. To assess the bacterial community composition, we first undertook genome recovery, yielding a dereplicated set of 95 metagenome-assembled genomes (MAGs) representing 15 bacterial families (Table S1). Read mapping to these MAGs and a set of publicly available genomes confirmed the development of disease was associated with a substantial shift in the bacterial community (Figures S1 & S2 & Table S2). The composition of the disease-associated microbiome of WT and IL-17RA^{–/ –} mice was distinguishable regardless of cohousing treatment (Figure S1B). Where cohousing occurred pre- or post-transplant only, ordination analysis suggested a more extensive difference in the disease-associated gut community between WT and IL-17RA^{–/ –} mice as compared to continuously cohoused mice (Figure S1B). Continuous cohousing was also associated with a significant decrease in community diversity in WT mice post-transplant, whereas there was no significant decrease when cohousing occurred pre- or post-transplant only (Figure S3). These data suggest that WT mice cohoused pre- or post-transplant experience a more moderate community shift associated with disease onset than continuously cohoused WT mice, in agreement with the delayed disease phenotype observed in these mice (Figure 1).

Members of the dominant bacterial family Muribaculaceae are transferred during cohousing and have the genetic potential to contribute to the GVHD phenotype

To investigate the putative priming effect of pre-transplant cohousing we compared both the functional potential and composition of the microbiome of cohoused WT mice to both their pre-cohousing composition and to separately housed WT mice at the same time point. At the functional level, there was no obvious distinction between pre-transplant cohorts, i.e. WT pre- and post-cohousing, nor between WT mice post-cohousing and separately housed control mice (Figure S4A-C). Only one annotation was identified as consistently differential between WT mice pre- and post-cohousing in both Experiment I and II, while also being significantly different from separately housed WT mice; tetracycline resistance monooxygenase (K18221) (Tables S3–5). The limited change in the functional potential of the microbiome during cohousing suggests that the priming effect is subtle and not detectable at the level of whole community comparison. Alternatively, the functional effect lies outside of the annotated portion of the metagenomes, or is occurring at the transcriptional/translational level.

At the compositional level, five species from the family Muribaculaceae (formerly uncultured lineage S24-7²²) increased during cohousing in WT mice from both Experiment I and II (Figure 2 & Table S6). These species were also significantly more abundant in cohoused WT mice than in separately housed WT mice (Table S7). The only non-Muribaculaceae species to be enriched during cohousing in both experiments was Prevotella sp002933775; however, the average post-cohousing abundance of this species was substantially lower than that of the Muribaculaceae species (0.3% vs 1.7%; Table S6). Only two species were depleted consistently across both Experiment I and II: one from the family Muribaculaceae and one from the family Lachnospiraceae (Table S6). Notably, the pre-cohousing abundance of the depleted Lachnospiraceae species varied considerably between mice, ranging from 0.003% to 3.1%, and also did not differ in abundance between WT mice post-cohousing and separately housed WT mice (Table S7). In contrast, the depleted Muribaculaceae species was found to be significantly depleted in cohoused WT mice in comparison to separately housed WT mice (Table S7). We also undertook a broader analysis of all pre-transplant mice from all three experiments and identified a distinct microbiome composition between those typically susceptible to early-onset GVHD post-transplant (IL-17RA^{–/ –} and pre-transplant cohoused WT mice) and those typically exhibiting delayed disease (separately housed WT mice) (Figure 2b & Table S8). Comparison of MAGs associated with each disease type did not reveal any significant functional distinctions across the Pfam, KEGG and CAZy databases (Tables S9–11), further supporting the priming effect as occurring below the level of these functional categories or being potentially driven by variation in expression. Multivariate analysis of all pre-transplant WT mice also revealed a distinction between non-cohoused and cohoused WT mice (Figure 2d). Muribaculaceae species were consistently identified as critical to both this division and that distinguishing early- and delayed-onset disease (Figure 2c and e & Table S8). These findings are in agreement with the results of our previous 16S rRNA gene sequencing analysis of pre- and post-transplant cohoused WT mice that identified members of the family Muribaculaceae as the dominant taxa with altered relative abundance via cohousing.²¹

Figure 2. Species from the family Muribaculaceae are transferred during cohousing and are associated with disease susceptibility. (a) Heatmap displaying genomes identified as consistently significantly enriched or depleted in WT mice during pre-transplant cohousing across Experiments I and II (n = 4–5 per group per experiment). Complete analysis output and relative abundance values are contained in Table S6. (b) PCA displaying all pre-transplant mice across all three experiments, I, II and III, based on read mapping to a genome database composed of recovered MAGs plus genomes from NCBI (see methods). Dotted line indicates division between mice susceptible to early-onset (IL-17RA^−/− and cohoused WT) and delayed-onset (separately housed WT) disease. Week 0 mice are pre-transplant mice from Experiment III and have undergone no cohousing pre-transplant. Week 4 mice have either been cohoused for 4 weeks or maintained separately for the same time period (control mice) and originate from Experiment I and II. (c) Heatmap displaying genomes with significantly different abundance between early-onset and delayed-onset GVHD susceptible mice, using division indicated in (b). Complete analysis output and relative abundance values are contained in Table S8. (d) Multivariate sparse partial least squares discriminant analysis (sPLS-DA) based on read mapping to genome database composed of recovered MAGs plus genomes from NCBI of all pre-transplant WT mice from Experiments I, II and III. (e) Genome-based signature contributing to separation along component 1 of (D). Color indicates group in which median genome abundance is highest. Bar length corresponds to loading weight. CH: cohoused, SH: separately housed.

Of the enriched Muribaculaceae species, four are members of the genus CAG-873 and one is a member of the genus CAG-485, both genera currently lacking cultured representatives (Figure S5). The Muribaculaceae species depleted during pre-transplant cohousing also belonged to genus CAG-873 (GVHD27), suggesting priming effects must be at least species-specific. Two species (CAG-485 sp. GVHD19 and CAG-873 sp. GVHD29) displayed a relative abundance pattern of particular interest; both increased consistently in WT mice during cohousing experiments and were also observed in WT mice following disease onset (Figure S6). To assess whether these species displayed similar functional profiles we compared their representative genomes to others from the family Muribaculaceae using functional annotation of predicted proteins, including all species enriched in WT mice during cohousing (Figure S7A). The enriched genomes did not appear distinct amongst the family, indicating their broad metabolic profile was typical of Muribaculaceae. Based on their carbohydrate active enzyme complement, both GVHD19 and GVHD29 are part of the α-amylase trophic guild,²³ a guild defined by a reduced set of glycoside hydrolases (Figure S7B). No Muribaculaceae species enriched in WT mice during cohousing were members of the host glycan guild, characterized by increased abundance of enzymes associated with the degradation of mucin,²³ suggesting that they are not priming the gut for disease by promoting barrier dysfunction and subsequent pathogen translocation.²⁴

We subsequently undertook a comparison of orthologous proteins amongst Muribaculaceae genomes, specifically seeking those present in species enriched in WT mice during pre-transplant cohousing but rare in other species and with a predicted function contributing to capacity for disease priming such as adhesion or protein interaction domains, peptidases or surface antigens (Table S12). Multiple C10 peptidases were identified in the enriched genomes that displayed limited sequence conservation within the family Muribaculaceae (Table S12). The C10 protease family includes streptococcal pyrogenic exotoxin B (SpeB) from Streptococcus pyogenes capable of degrading multiple components of the host immune system.²⁵ Within GVHD19, several C10 peptidases are encoded directly adjacent to a putative serine protease inhibitor (serpin) that may act as a regulator of peptidase expression.²⁶ GVHD19 also encodes two-class C25 peptidases, of which gingipain produced by Porphyromonas gingivalis is the reference enzyme, and both GVHD19 and GVHD29 encode M6 peptidases, metalloendopeptidases shown to be important in degrading immune system elements.²⁵ Other proteins of interest present in the enriched genomes include homologs of potential adhesins such as the leucine-rich repeat encoding internalin J, most closely resembling that of Porphyromonas spp. (Figure S8A), which may play a role in biofilm formation,²⁷ and the cleaved adhesin domain carrying hemagglutinin, resembling that of Prevotella spp. (Figure S8B-D) which confers capacity to adhere to host cells.²⁸ Eukaryotic-like domains such as leucine-rich repeats, tetratricopeptide repeats and fibronectin type III are present in multiple additional proteins representing further candidates with potential for host interaction.²⁹

Continuous cohousing with a dysbiotic microbiome promotes a disease-associated bloom of Enterobacteriaceae

Due to the large shift in gut microbiome composition associated with the onset of GVHD (Figure S1A), we initially examined the disease-associated microbiome within each group of WT mice in comparison to their pre-transplant composition at the bacterial family level. The family Enterobacteriaceae was significantly enriched in cohoused WT mice at the time of sacrifice regardless of when cohousing occurred (Table S13), with continuous cohousing associated with a substantially greater enrichment of this family (Figure 3 & Figure S2). WT mice undergoing continuous cohousing also displayed a significant enrichment in the family Acutalibacteraceae as well as a significant depletion of Muribaculaceae upon onset of GVHD (Table S13). No significant change in either family was observed in WT mice displaying GVHD from the alternative cohousing scenarios, although a similar trend was visible (Table S13, Figure 3a). Separately housed WT mice that developed GVHD also displayed a significant increase in Enterobacteriaceae with disease (Table S14) and there was no significant difference in the abundance of Enterobacteriaceae between cohoused and separately housed WT mice with GVHD in any experiment (Table S15). Enterobacteriaceae blooms therefore appear characteristic of disease regardless of the cohousing scenario.

Figure 3. Increased abundance of Enterobacteriaceae is associated with disease. (a) Relative abundance of bacterial families Enterobacteriaceae, Acutalibacteraceae and Muribaculaceae across each experiment, before and after transplant. (b) Multivariate sparse partial least squares discriminant analysis (sPLS-DA) based on read mapping to genome database composed of recovered MAGs plus genomes from NCBI (see methods) of post-transplant WT mice. (c) Genome-based signature contributing to component 1 of (B). Genomes with frequency >0.75 across validation replicates displayed. Color indicates group in which median genome abundance is highest. Bar length corresponds to loading weight. CH: cohoused, SH: separately housed.

At the species level, Escherichia coli comprised the majority of the significant increase in abundance of Enterobacteriaceae in cohoused WT mice post-transplant in all three cohousing scenarios (Table S16). Other members of the family were present at low abundance, with the exception of Enterobacter hormaechei subsp. steigerwaltii, which increased significantly in all experiments and markedly more so with continuous cohousing (maximum 14% relative abundance). The relative abundance of Bacteroides vulgatus, Lactobacillus murinus and a member of the genus CAG-180 from the family Acutalibacteraceae also increased significantly with disease progression in cohoused WT mice in all three cohousing scenarios (Table S16). Separately housed WT mice also experienced significant increases in each of these species with the exception of CAG-180, which only increased significantly in one experimental replicate (Table S17). Despite these similar microbiome transitions amongst WT mice from different cohousing scenarios, multivariate comparison of the disease-associated microbiome of all WT mice clearly indicated separation of continuously and partially (pre- or post-transplant) cohoused mice based on higher relative abundance of Enterobacteriaceae species in the former group (Figure 3b and c). Continuous cohousing may therefore promote acute disease through increased expansion of Enterobacteriaceae species beyond that experienced by mice under the alternative cohousing scenarios.

Disease-associated Enterobacteriaceae have no conspicuous GVHD-specific virulence factors

To determine whether the association observed between E. coli and E. hormaechei subsp. steigerwaltii and GVHD was driven by a specific set of virulence factors we undertook comparative analysis between all MAGs representing these species (non-dereplicated) and 21 high quality reference genomes not known to be associated with GVHD. A total of 28 E. coli and two E. hormaechei subsp. steigerwaltii MAGs (>90%-estimated completeness) were recovered across all experiments (Figure S9A-B). We additionally recovered four E. coli MAGs from publicly available datasets from patients undergoing hematopoietic stem cell transplantation; one adult that developed acute GVHD and three pediatric patients that did not.^30,31 A further E. hormaechei subsp. steigerwaltii MAG was also obtained from a pediatric patient that developed acute GVHD.³¹ Ordination analysis suggested that the virulence profile was driven more by phylogeny than by disease association, as MAGs from the current dataset clustered with their phylogenetic counterparts (Figure S9C). Inspection of individual virulence factors in the GVHD-associated E. coli MAGs also support differences being driven by phylogeny rather than disease association, e.g. iron acquisition capability and type VI secretion systems differed between the phylotypes (Figures S9A-B & S10, Table S18). Virulence factors that were shared across the majority of GVHD MAGs did not appear GVHD-specific; e.g., FdeC, capable of mediating adhesion to mammalian cells and found at higher prevalence in pathogenic subtypes,³² was also found in non-GVHD reference genomes and was absent from the clinical GVHD E. coli MAG (Figure S10). Unbinned contigs homologous to a public set of E. coli genomes were also analyzed to identify potential plasmid elements not initially clustered with the MAGs; however, no further virulence factors were identified (Table S18). The lack of a clear GVHD-associated virulence profile supports E. coli and E. hormaechei subsp. steigerwaltii as having an opportunistic role in disease progression.

Survival post-transplant is associated with microbiome composition

Across all three cohousing scenarios we encountered a range of survival times post-transplant amongst WT mice; some succumbed to GVHD rapidly, similar to IL-17RA deficient mice, while others did not develop disease for over a month (Figure 1). Using this distribution of length of survival post-transplant, we identified a negative correlation between the abundance of Enterobacteriaceae and survival of WT mice post-transplant, and positive correlation of the abundance of Muribaculaceae with survival (Figure 4a, Table S19). These trends were also observable at the genome level; however, there was clear variation in abundance between mice with similar survival times (Table S20). We subsequently defined three disease types for further comparison based on observable clustering of survival times post-transplant: hyper-acute (survival ≤10 days), intermediate (survival 20–35 days) and delayed (survival 45+ days) (Figure 4a). Ordination analysis suggested a distinction between the microbiome composition with hyper-acute disease (both WT and IL-17RA^–/–) and those with intermediate or delayed disease onset (Figure 4b and c). Hyper-acute GVHD was associated with a significantly higher abundance of E. coli and B. vulgatus (Figure 4d, Table S21). Intermediate or delayed disease was associated with higher abundance of members of the families Muribaculaceae and Lachnospiraceae, consistent with a greater retention of commensal organisms.

Figure 4. Length of survival post-transplant is associated with altered microbiome composition. (a) Correlation between relative abundance of Muribaculaceae and Enterobacteriaceae in WT mice at the time of sacrifice due to GVHD with length of survival post-transplant. Benjamini-Hochberg adjusted P values shown. (b) PCA displaying all mice sacrificed due to the development of GVHD across the three experiments (n = 4–5 per group per experiment) based on read mapping to genome database composed of recovered MAGs plus genomes from NCBI (see methods). (c) As per (B) with colored lines connecting mice within each survival timeframe. (d) Heatmap displaying genomes identified as significantly different between hyper-acute and delayed-onset GVHD at the time of sacrifice. Genomes with relative abundance >1% in at least one mouse displayed. Full list contained in Table S21. CH: cohoused, SH: separately housed.

Discussion

Here we investigate the effect of the pre- and post-transplant gut microbiome on the development of hyper-acute GVHD. Whilst we had previously observed accelerated disease development in WT mice following cohousing with IL-17RA^{–/ –} mice,²¹ it was unknown whether the phenotype was induced by exposure to the microbiome from IL-17RA^{–/ –} mice before or after transplant. Pre-transplant priming of the gut microbiome via cohousing was insufficient to significantly accelerate disease in isolation, while cohousing post-transplant resulted in a significantly higher morbidity of WT mice versus non-cohoused controls. However, disease onset following cohousing post-transplant and not pre-transplant was delayed in comparison to WT mice cohoused continuously, indicating pre-transplant priming is also necessary for generating the hyper-acute phenotype. These novel data support both a priming and exacerbating role for the gut microbiome in disease development. We observed transfer of members of the bacterial family Muribaculaceae during cohousing pre-transplant from IL-17RA^{–/ –} mice to WT mice and a bloom of members of the family Enterobacteriaceae post-transplant that associated with disease. Continuous cohousing was associated with an exaggeration of the expansion of Enterobacteriaceae, including increased abundance of both E. coli and E. hormaechei subsp. steigerwaltii.

As a dominant family in fecal microbiomes of laboratory mice,²³ it was unsurprising to see transition occurring amongst members of Muribaculaceae during cohousing. This family is regularly noted to vary in abundance in perturbation studies,^33,34 however, it is currently unclear whether the group plays a role in disease development. Functional studies have been impeded by a lack of cultured Muribaculaceae representatives with the first isolates only recently being reported.^22,35 Our use of metagenomic sequencing in the current study enabled us to identify, with species-level resolution, which members of the family were changing in abundance during cohousing. The improved resolution of metagenomic sequencing in comparison with amplicon sequencing means that species assignment will be directly comparable with future studies, including with new isolates as they become available. We used comparative genomics to compare the species of interest with available genomes from the family to compile a list of candidate genes with the potential to play a role in priming the gut to drive acute GVHD. A number of peptidases in this list suggest one possible priming mechanism is active degradation of components of the immune system rendering mice more susceptible to blooms of opportunistic pathogens post-transplant, a mechanism with potential to also operate in the human gut. While members of Muribaculaceae are found in the human gut, their prevalence is low (~7%).³⁵ However, the functions identified here as GVHD-associated are not exclusive to Muribaculaceae and therefore may be playing a similar role in other gut species. We identified peptidases from the families C10, C25 and M6, each of which has been demonstrated to have immunomodulatory activity and are identified in a variety of bacterial species.²⁵ We also observed a potential regulatory arrangement between peptidase and inhibitor similar to that described in Bacteroides fragilis where the gene pair is co-transcribed and responsive to environmental stimuli, particularly oxygen.²⁶ Increased oxygen availability associated with epithelial damage induced by transplant conditioning may therefore provide stimulus for peptidase expression while also promoting expansion of facultative anaerobes such as E. coli.³⁶

The bloom of E. coli we observed associated with the development of acute GVHD is consistent with previous observations both in mice and humans.^14,30 Expansion of E. hormaechei subsp. steigerwaltii has not been described previously, although unclassified Enterobacter spp. have been associated with bloodstream infection in patients post-transplant in connection with acute GVHD.^37,38 Using genome level analysis of the associated populations we confirmed potential pathogenicity of the identified E. coli supporting an active role in disease progression. While we observed a correlation between the abundance of Enterobacteriaceae and length of survival post-transplant, one limitation of these data is that they are only taken at the end-point of disease. There is a possibility that the species responsible for exacerbating acute GVHD bloomed early post-transplant and then diminished. In contrast, GVHD itself may directly modulate the microbiome and potentially hide preceding species of interest. For example, Paneth cells decrease in number and are damaged during GVHD, potentially contributing to altered microbiota via decreased production of the antimicrobial peptide α-defensin.^13,14 However, sampling more frequently post-transplant is difficult due to the severity of GVHD and the consequence of this on the overall health of the mice.

In the broader disease context, these data support both the pre- and post-transplant gut microbiome as critical factors in the susceptibility to acute GVHD. Multiple factors affect the gut microbiome during stem cell transplant complicating the design and application of approaches for management of the bacterial community. Prior to transplant, antibiotics may be administered prophylactically; however, despite clear beneficial effects in randomized studies,^39,40 some negative effects have been also been suggested in retrospective cohort analysis,⁴¹ perhaps consistent with the fact that the microbiome varies significantly between individuals.⁴² Conditioning regimens also affect pre-transplant microbiota, however, to what degree is difficult to establish due to prior and/or concomitant antibiotic administration.^43–45 Antibiotic treatment is often initiated post-transplant in response to infectious complications and may negatively impact overall prognosis depending on the antibiotics used.^46,47 As an alternative to microbiome depletion, prophylactic post-transplant manipulation of the microbiota via fecal microbiota transplantation (FMT) is currently under investigation as a means of restoring commensal organisms, thereby improving bacterial diversity. Whilst this research is only emerging and still largely experimental, early clinical data showing feasibility, safety and some efficacy of FMT to treat refractory acute GVHD in the GI tract is encouraging.^48–50 Restoration and/or maintenance of a diverse microbiome may increase resistance to blooms of opportunistic pathogens such as E. coli.⁵¹ Our results demonstrate association of specific species with disease outcome, where other members of the same genus appear unassociated, indicating a necessary resolution for manipulation that is unlikely to be achievable using antibiotics. Therefore, promoting bacterial diversity using methods such as FMT may be a more successful strategy for reducing acute GVHD development following stem cell transplantation.

Methods

Mice and allogeneic stem cell transplantation

Wild-type C57BL/6 mice (referred to as WT herein) were purchased from the Animal Resource Center (Perth, Western Australia). IL-17RA^{–/ –} mice on a C57BL/6 background (Amgen, Washington, USA) were bred in-house. Animal procedures were undertaken using protocols approved by the QIMR Berghofer Animal Ethics committee. Mice were transplanted and monitored as described previously.^21,52 Briefly, recipient mice received total body irradiation of 1000 cGy. Recombinant human Granulocyte Colony Stimulating Factor (G-CSF; Amgen Inc., Thousand Oaks, CA, USA) was administered to donor mice subcutaneously at 10 µg/dose/animal for 6 days. Mice were transplanted with either 25 × 10⁶ T-cell replete or 20 × 10⁶ T-cell depleted (TCD) G-CSF mobilized splenocytes. All transplanted mice were housed in sterilized microisolator cages and received acidified autoclaved water. GVHD was assessed using an established scoring system⁵³ and mice with clinical scores ≥6 were sacrificed in accordance with institutional guidelines. The data from the continuous cohousing experiments (Experiment I in Figure 1a) has been published previously.²¹

Metagenomic data processing, assembly and MAG recovery

Mouse reads were removed by mapping against the Mus musculus genome (GRCm38.p5) using BWA v0.7.12⁵⁴ requiring a minimum alignment length of 30 bases and maximum of 15 clipped bases for reads to be considered of mouse origin. Non-mouse reads from each sample were assembled independently using metaSPAdes v3.12.0.⁵⁵ Coverage based binning was performed by read mapping of a subset of six samples to each resulting assembly using BamM v1.7.3 (https://github.com/ecogenomics/BamM) with bins recovered using MetaBAT v2.12.1⁵⁶ with a minimum contig length of 2000 bases. Contamination and completeness of resulting bins was assessed using CheckM v1.0.12.⁵⁷ Bins with completeness >75% and contamination <7% were dereplicated using dRep v2.2.3⁵⁸ at 95% minimum secondary average nucleotide identity (-sa 0.95). The taxonomic affiliation of recovered MAGs was determined using GTDB-Tk v0.3.0 with the Genome Taxonomy Database (GTDB) Release 04-RS89.⁵⁹

Metagenomic community profiling

Reads for each sample were mapped to a dereplicated set (dRep, 95%) of 16,958 bacterial and archaeal genomes from GTDB Release 03-RS86 (>80% complete, <7% contamination) using BamM with a minimum seed length of 25. Genomes with >1x coverage of >1% of the genome and overall coverage of 0.01X, as determined using Mosdepth v0.2.3,⁶⁰ were retained and combined with the dereplicated set of recovered MAGs for assessment of community composition. The final genome set of 237 comprised 95 MAGs recovered from this study plus 142 genomes from NCBI (NCBI accessions are contained in Table S23). Read counts for the final genome set were determined for each sample via mapping using BamM with minimum seed length of 25 bases and subsequent filtering for minimum mapping percentage identity of 95%. Per genome read counts were scaled to account for genome size whilst maintaining the raw unmapped read percentage for each sample as a reflection of unrepresented diversity. Relative abundance was calculated using scaled read counts as a fraction of total non-host reads per sample. Alpha-diversity was calculated using QIIME v1.8.0⁶¹ with counts normalized using the size factor method implemented within the R package DESeq2 v1.20.0.⁶² MAGs were annotated using Prokka v1.12.⁶³

Functional annotation

Functional annotation of raw reads and predicted proteins from MAGs was undertaken via alignment with HMMER v3.1b2⁶⁴ to the hidden Markov model databases dbCAN CAZy v6,⁶⁵ Pfam r32⁶⁶ and TIGRFAM v15⁶⁷ (MAGs only) with maximum e-value cutoff of 1e-10. KEGG orthology was determined via BLAST v2.8.1⁶⁸ alignment to UniProt UniRef100 database downloaded July 2017⁶⁹ with maximum e-value of 1e-10 and subsequent extraction of associated KO terms.

Muribaculaceae MAG analysis

A Muribaculaceae genome tree was constructed using recovered MAGs plus publicly available genomes based on alignment of 120 single-copy marker genes.⁵⁹ A bootstrapped maximum-likelihood tree was inferred using IQ-TREE v1.6.9⁷⁰ (100 replicates, non-parametric) using model LG+C10 + F + G with posterior mean site frequency approximation⁷¹ based on alignment positions containing a residue within ≥50% of sequences (38,690 positions). Trophic guilds were assigned using the method described in ref.²³ Orthologous proteins were identified using Proteinortho v5.16b⁷² with an e-value cutoff 1e-6. Gene trees were constructed using homologues identified within the recovered Muribaculaceae MAGs and the GTDB Release 03-RS86 by GeneTreeTK v0.0.14 (https://github.com/dparks1134/GeneTreeTk) with default settings. Bootstrapped maximum-likelihood trees (100 replicates, non-parametric) were constructed using IQ-TREE v1.6.9⁷⁰ with ModelFinder⁷³ used for model selection. Alignment sites with a minimum similarity of 30% were used for phylogenetic inference. ARB⁷⁴ was used for alignment filtering and tree curation.

Escherichia coli MAG analysis

Virulence factor prediction in binned and unbinned E. coli contigs was undertaken via BLAST search of recovered MAGs and public genomes against the VFDB core database⁷⁵ (accessed February 2019) with maximum e-value cutoff 1e-10, minimum alignment identity of 30% and minimum alignment fraction of 70%. Virulence factor classes were obtained from reference genomes downloaded from VFDB VFanalyzer. Unbinned contigs from all mice at sacrifice were clustered using CD-HIT⁷⁶ using the cd-hit-est command with 99% identity (-c 0.99) and 90% coverage of the shorter sequence (-aS 0.9) thresholds. Representative contigs homologous to E. coli were identified via BLAST search of the NCBI nt database (downloaded January 2019) with maximum e-value cutoff 1e-10 followed by filtering for sequence description containing ‘Escherichia coli’, minimum alignment identity of 97% and minimum alignment fraction of 50% of the query contig. Additional E. coli MAGs were recovered using the same methodology described above from SRA experiments SRR3340629, SRR5050584, SRR5050585, and SRR5050587. Bootstrapped maximum-likelihood tree was inferred using IQ-TREE as described above.

Statistical analysis and graphical presentation

Survival curves were plotted using Kaplan-Meier estimates and compared by log-rank analysis. The Wilcoxon rank-sum test was used for the statistical analysis of alpha-diversity values with Benjamini-Hochberg adjustment for multiple comparisons. P < .05 was considered statistically significant. Dot plots are presented as mean ± standard error of the mean. Principal component analysis was conducted using the R package vegan v2.5–1⁷⁷ on data normalized using log cumulative-sum-scaling (log-CSS) implemented within metagenomeSeq v1.22.0.⁷⁸ Differential abundance of bacterial taxa and functional annotations (raw reads) between sample groups was assessed using the Wald test within DESeq2 v1.20.0⁶² based on total annotated read counts or read counts per genome scaled to account for genome size with the Benjamini-Hochberg adjustment for multiple comparisons. P < .001 was considered statistically significant. sPLS-DA analysis was conducted using the R package mixOmics v6.3.2⁷⁹ using centered log-ratio transformed relative abundance values (pseudocount 1e-07) with 50xM-fold cross-validation (fivefold pre-transplant WT, threefold post-transplant WT). Comparison of functional annotations within MAGs associated with early-onset and delayed-onset GVHD was undertaken using EnrichM v0.5.0 (https://github.com/geronimp/enrichM) with Fisher’s exact test employed to compare the number of MAGs encoding/not-encoding each functional category (i.e., each individual CAZy, KO, or Pfam classification) associated with each disease type (with Benjamini-Hochberg adjustment for multiple comparisons). Spearman’s rho calculated using ‘corr.test’ function within R package psych v1.8.12⁸⁰ with the Benjamini-Hochberg adjustment for multiple comparisons using centered log-ratio transformed relative abundance values. Bacterial families filtered for those with minimum relative abundance ≥0.05%. P value of <0.05 was considered statistically significant. Heatmaps were produced using the R package pheatmap v1.0.10,⁸¹ box plots produced using the R package ggplot2 v2.2.1⁸² with extension GGally v1.4.0⁸³ and dot plots of family abundance generated within GraphPad Prism v8.0.1 (GraphPad Software).

Data availability

Raw sequencing data and recovered MAGs are available via NCBI BioProject PRJNA544874. Sample accessions are provided in Table S22.

Disclosure of Potential Conflict of Interest

No potential conflicts of interest were disclosed.

Supplementary Material

Supplemental data for this article can be accessed on the publisher’s website.

Additional information

Funding

This work was supported by a Cancer Council Queensland Project Grant [1129969];National Health and Medical Research CouncilProgram Grant [APP1071822].