Skip to main content
Download PDF
  • Research article
  • Open access
  • Published:

Temporal dynamics in microbial soil communities at anthrax carcass sites

BMC Microbiology volume 17, Article number: 206 (2017) Cite this article

Abstract

Background

Anthrax is a globally distributed disease affecting primarily herbivorous mammals. It is caused by the soil-dwelling and spore-forming bacterium Bacillus anthracis. The dormant B. anthracis spores become vegetative after ingestion by grazing mammals. After killing the host, B. anthracis cells return to the soil where they sporulate, completing the lifecycle of the bacterium. Here we present the first study describing temporal microbial soil community changes in Etosha National Park, Namibia, after decomposition of two plains zebra (Equus quagga) anthrax carcasses. To circumvent state-associated-challenges (i.e. vegetative cells/spores) we monitored B. anthracis throughout the period using cultivation, qPCR and shotgun metagenomic sequencing.

Results

The combined results suggest that abundance estimation of spore-forming bacteria in their natural habitat by DNA-based approaches alone is insufficient due to poor recovery of DNA from spores. However, our combined approached allowed us to follow B. anthracis population dynamics (vegetative cells and spores) in the soil, along with closely related organisms from the B. cereus group, despite their high sequence similarity. Vegetative B. anthracis abundance peaked early in the time-series and then dropped when cells either sporulated or died. The time-series revealed that after carcass deposition, the typical semi-arid soil community (e.g. Frankiales and Rhizobiales species) becomes temporarily dominated by the orders Bacillales and Pseudomonadales, known to contain plant growth-promoting species.

Conclusion

Our work indicates that complementing DNA based approaches with cultivation may give a more complete picture of the ecology of spore forming pathogens. Furthermore, the results suggests that the increased vegetation biomass production found at carcass sites is due to both added nutrients and the proliferation of microbial taxa that can be beneficial for plant growth. Thus, future B. anthracis transmission events at carcass sites may be indirectly facilitated by the recruitment of plant-beneficial bacteria.

Background

The microbial composition of arid soils across the globe is distinct from other soil environments [1, 2]. The forces that shape arid soil microbial community composition include low water availability, temperature and UV radiation [3,4,5,6,7,8]. Especially, water restriction has a large influence since it affects several environmental factors, such as salinity, pH, and the availability of (in-) organic matter, which further modulate soil microbial diversity and activity [4, 5, 9,10,11]. The combination of the above factors may explain why soils in arid environments, such as deserts and arid savannahs, are taxonomically distinct from other soil types [1].

Most research on arid soil microbial diversity is directed to understand the community dynamics under changing environmental conditions, e.g. precipitation changes. However, soil communities can also be affected by the deposition of animal carcasses and their decomposition. This may have a profound impact on the soil microbiome as the carcass influences both biotic (adding new microbes) and abiotic (adding nutrients and moisture) factors. Moreover, after the introduction of a carcass, a succession within the microbial community will occur, changing abundances in accordance with nutrients being released during carcass decomposition [12]. For instance, bacteria belonging to Proteobacteria and Acidobacteria are the most common in soils during the initial stages of carcass decomposition, while Firmicutes are more prominent during active decomposition [13, 14]. The described succession is however, found under experimentally controlled conditions, where carcasses were secured at one location for the entire experiment. In contrast, in natural ecosystems carcasses are often consumed and/or dragged away from the site of death by scavengers [15]. Thus microbes and nutrients may only transiently enter the soil at the site of death [16], which might induce different microbial soil dynamics at natural carcass sites.

The present study investigates the effects of animal carcasses on soil microbial communities after an animal has died of an anthrax infection. The disease anthrax is caused by Bacillus anthracis, a gram-positive, rod shaped, sporulating bacterium [17]. This species belongs together with Bacillus cereus, Bacillus thuringiensis and several other Bacillus spp., to the B. cereus group [18], which is commonly found in soil as vegetative cells or as spores [19]. The B. cereus group bacteria are indistinguishable from each other using 16S rRNA gene sequences and show high genetic identity (> 99.6%) for certain housekeeping genes [18, 20, 21]. However, B. anthracis can genetically be distinguished from the other ‘species’ based on single nucleotide polymorphisms (SNPs) (e.g. in the plcR gene [22, 23]) and in most cases by the presence of two virulence plasmids specific to B. anthracis pXO1 and pXO2 [24,25,26]. The lifecycle of B. anthracis is different from other B. cereus group bacteria in that it predominantly targets, infects and kills mammalian herbivores, instead of insects as found for B. thuringiensis. Grazing by herbivores seems to be the main route for transmission in natural settings such as found in the semi-arid savannah in Etosha National Park (ENP), Namibia [27]. After causing the death of its host, B. anthracis returns to the soil where it sporulates. Hence, B. anthracis spores will often be found in high densities in the top layer of soil where haemorrhagic fluids have leaked from an anthrax carcass [27,28,29].

Long term measurements at carcass sites in the ENP have identified several processes occurring in and above the soil. It was found that B. anthracis cell counts in rhizosphere soils increased in the second year after carcass deposition, but not in surface soils [27]. Those results contrast experimental work where no multiplication of cells was observed in the second year [30]. Furthermore, it was observed that grass biomass and quality increased at the localised area of anthrax carcass sites [27]. This increase in localised plant growth is mostly likely due to nutrient release and/or rhizosphere plant-microbe interactions. The increase in above ground plant biomass results in attraction of grazers, which increases the potential exposure to the pathogen through ingestion of grasses and potentially roots / soil at carcass sites [16, 27]. It is thought that the bacterium sporulates shortly after entry into the soil and that it will stay dormant until favourable conditions arrive [31]. This suggests that there is only a short period, after the entry of B. anthracis into the soil, where the pathogen could interact with other microbes or plants to enhance its transmission [16, 30,31,32,33]. It is unclear if such interactions influence further transmission of the pathogen. In particular there is a lack of understanding of the dynamics of the (arid) soil microbial community after the influx of animal fluids (e.g. blood, gut contents) with high densities of B. anthracis vegetative cells. In order to address how interactions between microbes / plants and B. anthracis benefits pathogen transmission, it is important to understand the dynamics of the arid soil microbial community after carcass deposition and influx of B. anthracis vegetative cells into the soil.

There are several technical challenges in following both the B. anthracis and microbial community dynamics in soil. Currently, DNA-based methods such as 16S rRNA amplicon sequencing and shotgun metagenomics are established methods for the study of complete microbial communities [34]. As mentioned above B. anthracis is indistinguishable from other B. cereus group species on 16S rRNA gene sequences. Therefore, 16S rRNA amplicons are not suitable for studies aiming to distinguish members from this group. B. anthracis and many other microbes can be detected with shotgun metagenomics [35]. Again, B. anthracis is genetically highly similar to B. cereus group bacteria, which makes bioinformatic identification of these species challenging [18]. Therefore methods able to detect and distinguish various B. cereus group species are needed to complement and confirm the shotgun metagenomic results. Such methods include cultivation using selective media (e.g. PLET) and specific qPCR assays [23, 36]. The combination of these methods can then be used to monitor environmental B. anthracis populations and to generate hypotheses about possible species interactions, which subsequently can be tested with specifically designed experiments.

Here we present an analysis of soil microbial community composition following leaching of fluids into the soil from two spatially and temporally proximate zebra anthrax carcasses in ENP. The ENP is a semi-arid savannah environment, which for much of the year has sparse water availability, meagre vegetation and limited nutrient resources [37, 38]. Hence, carcass nutrients will likely have great influence on the localised soil microbial community [39]. Large herds of ungulates ensure that every year there is an abundance of carcasses in ENP, often through predation or diseases such as anthrax. It is estimated that up to 400 plains zebras (Equus quagga) per year die from anthrax infections in the ENP [40]. Due to the unpredictability of the presence of disease cases, and our interest in controlling for temporal and spatial variability in soil microbiota among study sites, our study was limited to two carcass sites situated proximately in space and time. Our aim was to describe the microbial community succession taking place in the soil after the influx of haemorrhagic fluids containing B. anthracis vegetative cells. Shotgun metagenomic sequencing of samples obtained through the first month of decomposition was employed to investigate temporal dynamics of the community structure and function. We investigate the temporal change of the taxonomic composition and the metabolic potential by identifying major metabolic pathways. Finally, we use a combination of techniques (qPCR, cultivation and metagenomic sequencing) to circumvent state (vegetative cell vs. spore) associated challenges to accurately track B. anthracis abundances.

Results

Carcass information and rainfall recording

On 03.03.2014 we identified two plains zebra (Equus quagga) carcasses less than 1 km apart in ENP. Hereafter referred to as Carcass 1 (Ca1) and Carcass 2 (Ca2). Ca1 was intact when sampled on day 0, while Ca2 was minimally scavenged (some intestine was dragged out the anus by vultures). The time of collection (≈ 14:00), the state of both carcasses and only the presence of a few avian scavengers and not mammalian, indicates that both animals were likely dead for less than 12 h and certainly fewer than 24 h. Our study contrasts with other carcass decomposition experiments, since scavenging was not restricted [12,13,14, 39, 41]. As such, the carcass nutrients and fluids will only leak onto the soil for a short time before scavengers consume soft tissue and move the remains off the site.

After three days both carcasses were completely consumed by scavengers and the bones were found approximately 5 m away from each sampling site. The sampling sites were visited at days: 0, 3, 7, 14, 21 and 30, to collect material for cultivation and DNA extraction. At day zero an uncontaminated control sample (Ctrl0) was taken at both sites to function as a reference sample. The soil in the study area has a alkaline pH around 8.7 ± 0.4 (Additional file 1: Table S1), low moisture and dominance of bacteria (Additional file 2: Figure S1), which is characteristic of arid soils [42].

There was some rainfall at Okaukuejo in the days prior to day 0 followed by little rain until day 18 when heavy rainfalls occurred (day 18–21, Additional file 3: Figure S2).

Detection of B. anthracis

DNA extraction efficiency of bacterial gram-positive cells and/or spores can be poor [43]. To control for extraction efficiency we spiked soil samples with 3.4 X 106 B. anthracis vaccine strain Sterne 34F2 spores (Onderstepoort Biological Products).

The qPCR showed that soil samples spiked with B. anthracis spores (see Additional file 4: for details) had an estimated recovery of 347,183 (Fig. 1a) and 654,419 (Fig. 1d) genomes representing 10.2% and 19.2% DNA extraction efficiency, respectively (equation in Additional file 4). QPCR on the carcass samples revealed that Ca1 had high abundance of B. anthracis at day 0 (Fig. 1a), while a similar peak occurred at day 3 for Ca2 (Fig. 1d). At both sites we find low B. anthracis abundances in the Ctrl0 samples. The presence of B. anthracis in these samples could either be due to spill-over of B. anthracis cells between contaminated and uncontaminated soils (see methods), or B. anthracis spores were present in those samples before carcass deposition. With our data it is not possible to determine which explanation is correct.

Fig. 1
figure 1

Estimation of B. anthracis abundance in soil samples using qPCR, metagenomic reads mapping and cultivation. For all panels, sampling time-points are on the x-axis. a and d qPCR results with estimated number of B. anthracis genomes per gram soil along the y-axis at different time-points (x-axis) for Carcass 1 (Ca1) (a) and Carcass 2 (Ca2) (d). The spiked sample is control soil from day 0 with 3.4 × 106 Stern34F2 B. anthracis spores added. P1 and P2 represent the two replicates at each time-point. b, c, e and f shows the results from mapping the metagenomic reads against the B. anthracis isolate genomes (K1 /K2) isolated from Ca1 and Ca2 as well as the other reference strains; B. cereus E33L, B. thuringiensis HD-771 and B. subtilis 168. Reads matching the references are on the y-axis. b and e and (c and f) shows the mapping of the metagenomes against the K1 and K2 strain, respectively as well as the other reference strains. Ca1 in panes (B and C) and Ca2 in panes (E and F). (G) B. anthracis spore counts of soil samples from Ca1 and Ca2 in spores per gram dry soil. Spore counts were estimated by culturing of heat-shocked soil samples and adjusted to dry weight of soil

Full size image

In order to distinguish B. anthracis metagenome reads from those of B. cereus and B. thuringiensis, mapping of reads was performed using the aln algorithm (Burrows-Wheeler Aligner (BWA-aln)) [44] with very strict mapping parameters (Additional file 4: Methods). In addition, we added closely related Bacillus spp. reference strains as bait for sequences not unique to B. anthracis [45]. (Note that without the usage of closely related reference strains, mapping of metagenomic reads against B. anthracis becomes highly unreliable). Such mapping against the genomes of two B. anthracis strains, K1 and K2 isolated from the carcasses studied here [46], resulted in a similar pattern as observed in B. anthracis specific qPCR experiments (Fig. 1a, d) with abundances peaking at day 0 for Ca1 and at day 3 for Ca2. For both carcasses there are no significant differences between the mapping results when using the K1 (Fig. 1b, e) or K2 (Fig. 1c, f) genomes. There are few reads (2–94) mapping to B. cereus, B. thuringiensis or B. subtilis in the Ctrl0 samples (Fig. 1, Additional file 5: Table S2). The frequencies for these three species peak at day 3 at both carcass sites, but abundances are lower than for B. anthracis (Fig. 1b–f). After day 3 the numbers of reads mapping to all Bacillus spp. decrease gradually, except for a slight increase in Ca1 at day 30. A close inspection of the mapping process revealed that many reads with low mapping qualities, e.g. poor alignments due to mismatches, were removed in our final filtering step (Additional file 4: Methods). Removal of reads with mismatches did not change the abundance pattern of our time-series for B. anthracis, while there was a significant change for the other species used (Additional file 6: Figure S3). This difference is likely due to the presence of DNA sequences in the metagenomes related, but not identical to the B. cereus group genomes available and used as references.

Culturing from soil samples that had been heated to kill vegetative cells, showed no B. anthracis spores present at days 0 and 3 (Fig. 1g) in contrast to the qPCR results. Ca1 had >100,000 spores on day 7 and 14, followed by a sharp drop in spores at day 21, and an increase to >50,000 spores at day 30. Ca2 had >2000 spores on day 7 and >60,000 spores on day 14; at the rest of the time-points Ca2 only had <130 spores. Some of these aberrations from the general trend in counts may be products of sampling biases. In addition, the reduction of B. anthracis levels on day 21 (Fig. 1g) could also be due to a response to the rainfall prior to day 21 (Additional file 3: Figure S2).

Temporal changes of GC content and average genome size

Changes in the average % guanine-cytosine (GC) content of metagenomics data indicate changes in microbial community composition [47]. The average GC-content of the time-series metagenomes ranged between 54.6 and 69.0% (Table 1). Interestingly, in both time-series there was a marked drop in average %GC content at day 3 and 7 before it increased again at day 14. This drop was more prevalent in Ca1 than Ca2. This change is due to addition of an extra %GC content peak (≈43%) in the %GC profile, which reduces the average GC-content, and the skewness of the GC-content (Table 1).

Table 1 Overview metagenome shotgun samples
Full size table

The drop in GC content corresponded with a drop in average genome size (AGS) and an increase in genome equivalents (Table 1) [48]. The AGS drop and the shift in average % GC suggests that the carcass soil community composition is changing due to actively growing microbes responding to the carcass nutrients. In addition, the drop in AGS may indicate that the metabolic capacity of the soils changes over the time-series, since AGS can be used as a measure for the metabolic complexity of an ecosystem [49].

Taxonomic classification of shotgun sequences

Many tools exist for the taxonomic classification of metagenomic shotgun sequences. Our aim was to use a tool that captures most of the diversity present in the samples. We therefore compared metaxa2 [50], which classifies rRNA sequences in the metagenomes using a modified Silva SSU /LSU database, with metaBIT (wrapper around Metaphlan2) [51, 52], Kraken [53] and MEGAN [54]. Metaphlan2 and Kraken uses a database derived from microbial whole genome sequences, with the difference that Metaphlan2 only uses signature sequences to identify taxa, instead of whole genomes. MEGAN uses the Non-Redundant protein database from NCBI (NR) for classification.

By comparing the relative abundance of reads classified by each tool we found that the tools differed in the amount of reads classified per sample (Additional file 7: Table S3; Additional file 8: Figure S4) [55]. Interestingly, Kraken showed different temporal dynamics with respect to the change in relative abundance of classified reads compared to the other tools. MetaBIT shows an especially large change between days 0 and 3 compared to the other tools (Additional file 8: Figure S4). We assume that in both cases these results are due to the databases used by metaBIT and Kraken. This is also reflected in the number of prokaryotic taxa identified, with MEGAN detecting a maximum of 180 prokaryotic orders, while the three other tools detected fewer taxa (metaBIT:18, Kraken: 121 and metaxa2: 159). This shows that detection of environmental prokaryotes is significantly influenced by database choice.

The results of the above comparison suggest that metaxa2 is the best available choice among the four tools tested here to describe the bacterial and eukaryotic composition of the soil communities in our study (Additional file 7: Table S3), since the database it relies on is both well curated and contains the widest taxonomic range available. Moreover, MEGAN is a good complement to metaxa2, since it uses most reads for taxonomic and functional classification and can provide diversity measures that support community comparisons (e.g. PCoA, clustering etc).

The metaxa2 data was used to analyse the community diversity using rarefaction and rank-abundance plots. The number of classified reads per sample varied from 19,661 to 110,996 corresponding to about 50% of the identified rRNA sequences (Additional file 2: Figure S1a). The classified reads were dominated by bacteria, with a relative abundance of 70% ± 0.13% and 69% ± 0.1% for Ca1 and Ca2, respectively (Additional file 2: Figure S1b), which stayed constant regardless of time-point. The rRNA reads clustered into 557 OTUs using 97% similarity cut-off. The two Ctrl0 samples have the lowest number of OTUs, while the highest number of OTUs was found in Ca1 day 7 (Ca1_7) (Fig. 2a). The rarefaction curves (Fig. 2a) indicate that we have captured the dominant taxa present, since the curves start levelling off at around 5400 sequences. Nonetheless, the curves do not become horizontal, which indicates that additional low abundant taxa can still be detected with additional sequencing effort. Ca2 had fewer total rRNA reads than Ca1 at all time-points, except the control samples (Additional file 9: Table S4), which results also in fewer OTUs for Ca2 time-points compared to Ca1. An OTU rank abundance plot (Fig. 2b) indicates that Ca1 has higher OTU abundances for the dominant OTUs at all the time-points compared to Ca2. For both carcasses, the soil microbial community is dominated by a few OTUs on days 3 and 7. At the other time-points the evenness of the community is higher as seen in the day 30 samples where almost all OTUs have similar abundances (Fig. 2b). The OTU abundance distribution for the different time-points suggests different soil communities for each of the samples. And indeed, when we analysed community diversity differences using PCoA, based on MEGAN shotgun read classification, we found a clear separation between the communities (Fig. 3). Here metagenomic sequences from both carcasses follow a similar trend, where they are most divergent from the Ctrl0 and day 0 samples at days 3 and 7, before clustering closer to the Ctrl0 and day 0 sample again on days 14 and 21. Ca1 however, has a deviation from the trend on day 30: this sample is close to the Ca1 day 14 sample. The day 30 sample for Ca2 is similar to both Ctrl0 samples.

Fig. 2
figure 2

OTU richness of 16S /18S rRNA sequences from soil metagenomic samples. 16S /18S rRNA sequences were extracted with metaxa2 and diversity was analysed with MetaAmp. a Rarefaction curves. X-axis indicates number of rRNA sequences, and y-axis shows OTU numbers b Rank abundance of OTUs, where OTU rank is on the x-axis and the abundance per OTU is on the y-axis. The solid lines are Ca1 samples and the dotted lines are Ca2 samples. Line colour indicates time-points as indicated in the legend

Full size image
Fig. 3
figure 3

Soil microbial community relationships by time-point and carcass site. Principal coordinate analysis of a distance matrix created by normalised total counts and using Bray-Curtis dissimilarity. Ca1 is represented by the circles and Ca2 by the triangles, the time-points are visualised in different colours and the arrows are pointing in the direction of increasing days, red arrows for Ca1 and light blue for Ca2. The dotted arrows show the relationship between the day 30 to the Ctrl0 sample

Full size image

Temporal dynamics of the soil microbial communities

The 50 most abundant orders from the metaxa2 analysis were extracted and visualised for Ca1 and Ca2 to track fluctuations of the soil community over a 30-day time-period (Fig. 4). For both carcasses there is a clear relative abundance change in part of the microbial community after the influx of body fluids, while another part of the community does not seem to respond to this influx. Relative abundances may reflect an increase or decrease of a taxon due to changes in the abundance of another taxon, which can be misleading. We therefore show in Fig. 4b the log5 average fold changes between the time points of the raw counts for both carcass communities. This illustrates both the average direction of change over time and the average size of the read abundance change.

Fig. 4
figure 4

Temporal dynamics of microbial order abundances at carcass 1 and 2. a Heatmap visualisation of the relative abundance for the 50 most abundant orders at each time-point for carcass 1 and 2. b Barplots that show the log5 average fold change of the raw reads counts across the time-series for the 50 taxa shown in figure (a). The taxa in the heatmaps (a) are sorted using the order abundances in the control sample of Carcass1 at day 0 (Control). In order to visualise the relative abundances (a) values were log5 normalised and then scaled so that the sum of each column equals 1. Eukaryotic orders are marked with *, the remaining orders are bacterial. The bottom seven entries have not been classified to lower phylogeny than Kingdom, Phylum or Class

Full size image

Several orders show a sharp abundance increase after day 0 (log5 average fold change >1) at both carcass sites. These include Bacillales and Pseudomonodales species. In contrast, bacteria belonging to orders such as Frankiales, Rhizobiales and Solirubrobacteriales do not show large changes in their read abundances after the blood/nutrient influx at either carcass site. This is despite the large increase in total reads classified (Additional file 2: Figure S1). Finally, there are also orders that show limited increases. The abundance increase can be early (T = 0 days) and then decline (e.g. Clostridiales, Bacteriodales), or increase only later (T ≥ 7 days) in the time-series (e.g. Burkholderiales, Corynebacteriales).

The orders Bacillales and Pseudomonadales have the highest abundances at days 3 and 7 (Additional file 10: Figure S5), with genera such as Acinetobacter, Lysinibacillus and Kurthia being dominant (Additional file 9: Table S4). These orders are present at all time-points for both carcasses, but are especially abundant on days 3–14. Many of the genomes in the NCBI database from these genera are small (Size <3 Mbp) and their presence can explain the drop in average genome size (AGS) that we observed at day 3 (Table 1).

Between days 14 and 21 the relative rRNA read abundance for Bacillales drops at Ca1 and Ca2 by 33 and 23%, respectively. For Pseudomonadales, we identify a drop of 17 and 13% for Ca1 and Ca2. In total 9 orders show a relative rRNA read abundance drop (>10%) by day 21 in Ca1 and 12 orders in Ca2. At day 30 the abundances of both Bacillales and Pseudomondales increased again in Ca1, but not in Ca2. Interestingly, for the “stable” orders there was an overall increase in relative rRNA read abundances at day 21 and a subsequent drop at day 30.

The abundance changes observed at days 21 and 30 are likely influenced by heavy rainfall (16–103 mm) recorded at all the weather stations in the ENP in the three days prior to sampling day 21 (Additional file 3: Figure S2). Thus, we assume that the abundance increase at day 30 for many orders could be a response to this heavy rainfall. Moreover, the precipitation could explain the species abundance profile changes since it can result in short-term changes in the microbial community [56].

In addition to bacterial orders that are normal for soils, we also observed representatives of several bacterial orders likely introduced from the zebra. For instance, Fusobacteriales are present in the first days for both Ca1 and Ca2, where they are most abundant on day 0 for both samples. They are completely absent on day 14 and 30 in Ca1 and from day 7 onwards in Ca2. The increase in abundance for this order on day 0 is most likely due to the introduction of Fusobacterium spp. such as Fusobacterium equinum, which is a known inhabitant of the oral cavity and lower respiratory tract of horses [57]. In Ca1 sequences classified to the order Pasteurellales are highly abundant on day 0, but completely disappear after day 7. Their abundance at day 0 is mainly caused by Actinobacillus spp. and Pasteurella caballi (Additional file 9: Table S4). The later species is a commensal of the upper respiratory tract of horses [58].

Microbial community metabolism

The soil metagenomes at the early time-points show a decline in AGS compared to the control sample and the later time-points (Fig. 5). AGS is ecologically informative, as in general, soil bacteria have larger genomes than specialised, opportunistic, parasitic and symbiotic bacteria [59]. Thus AGS change suggests taxonomic and metabolic changes taking place within the microbial community. We therefore correlated KEGG pathway abundances of 112 pathways with AGS to study the temporal change of metabolic potential in the metagenomic time-series. For Ca1 we identified 29 pathways that were significantly correlated with AGS (p < 0.05), with 23 being specific to Ca1 (Fig. 5, Additional file 11: Table S5, Additional file 12: Table S6). Six pathways showed a negative correlation with AGS for both Ca1 and Ca2. For the Ca2 metagenomes we identified 7 pathways correlated with AGS, with only one pathway (aminobenzoate degradation) positively correlated with a large average genome size in Ca2 but not Ca1. In Ca1 six pathways show a decrease in abundance that is correlated with the decrease of AGS (Fig. 5, top panel). These pathways show an abundance drop at days 3 and 7 and then a subsequent increase. In Ca2 these pathways show a similar pattern (Additional file 11: Table S5, Additional file 12: Table S6). For the metabolic pathways: novobiocin metabolism (KEGG map00401), streptomycin metabolism (KEGG map00521) and naphthalene degradation (KEGG map00626), the genes identified are not unique for the pathways. For example, no reads were mapped to the genes involved in the final step of novobiocin production and only a few reads to the final step of streptomycin production. The later examples highlight the difficulty in functional annotation of shotgun metagenomic data and indicate the need for cautious interpretation of such results [60, 61].

Fig. 5
figure 5

Heatmap of significantly different metabolic pathways per time-point. KEGG metabolic pathways significantly correlating with AGS for microbial communities of Ca1 and Ca2. Pathways were determined by MEGAN classification. KEGG pathways abundances were normalised with DESEQ2 [96] and correlated to the AGS using the Spearman correlation method with the False Discovery Rate (FDR) test to calculate probabilities (FDR cut-off at 0.05). KEGG-pathways correlating positively or negatively with a p-value <0.05 are shown for Ca1 (blue) and Ca2 (red) for all the time-points, the AGS is shown in the top panel. # indicates pathways that are significant in both Ca1 and Ca2, § indicates pathways that are significant in Ca2, and the rest are only significant in Ca1. Pathway abundances were centred and scaled per row and positive and negative correlations were clustered based on the sign of the Spearman correlations

Full size image

Discussion

Shotgun metagenomic sequencing, which allows detection of seemingly all species in a community based on the presence of their DNA, is generally considered a superior method to culturing [1, 62, 63]. This is because it addresses “the great plate count anomaly”, which claims that only about 1% of the microorganisms seen during microscopy can be cultured [64, 65]. However, DNA-based analyses of microbial communities also have technical pitfalls, particularly when studying spore-forming microorganisms such as B. anthracis, which can result in large discrepancies between actual and estimated abundances [43]. For B. anthracis this is further confounded by the high levels of DNA sequence identity shared with other members of the B. cereus group.

The technical pitfalls of DNA-based analyses are evident from our data, where B. anthracis spores could readily be cultured from our samples, but were not as easily detected in metagenomes without using a very strict mapping approach. The onset of B. anthracis sporulation occurs within the first 72 h after carcass deposition and can continue up to eight days post-mortem [29, 66]. The qPCR data (Fig. 1) indicates that over time B. anthracis abundance peaks before almost completely disappearing from the samples. These results, together with the cultivation data, suggest that the reduced quantities of B. anthracis observed in the qPCR experiments are due to reduced DNA extraction efficiency (because of sporulation) [67]. This also suggests that abundances of B. anthracis and similar organisms are likely underestimated in the metagenomes at several time-points due to low extraction efficiency. Nonetheless, the metagenomes still captured a large part of the microbial community present in the carcass site soils, and can therefore be used to study the temporal changes taking place.

Why does a carcass promote plant growth?

The data generated reflects the temporal dynamics of the soil microbial community in a natural ecosystem after inoculation of a pathogen (B. anthracis), other host-associated microbes and nutrient influx from zebra carcasses. Turner et al. [16] showed that carcass sites in the ENP have higher quality and more abundant vegetation, which could be related to higher levels of phosphate, nitrogen and lower pH than the surrounding soils.

The microbial communities at both carcass sites show a clear response to the influx of bodily fluids, with similar community composition (at order level) of the dominant taxa for the two carcasses with only minor variation in abundances (Fig. 4). Shifts in the community are indicated by several observations: changes in relative abundance of classified reads between the samples (Additional file 3: Figure S2), changes of the AGS of the community studied (Fig. 5) and by variation in OTU diversity and abundance for each of the time-points (Fig. 2). The AGS decreases in the first week of the sample period before increasing again, suggesting a community shift towards copiothrophic bacteria in the first week (Table 1). This is expected as copiothrophic bacteria are better adapted to local increases in nutrient availability than typical soil bacteria and therefore can increase drastically in a short amount of time as observed here [68].

In contrast to the copiothrophic bacteria, the typical soil microbial community consisted of orders like Frankiales, Rhizobiales and Solirubrobacterales (Fig. 4), which remain relatively stable throughout the sample period. This suggests that these bacteria do not need, or are unable to use, the nutrients from the carcass. Moreover, Rhizobiales are nitrogen-fixing bacteria living in symbiosis with plants that thrive under nitrogen poor conditions [69]. The extra nitrogen provided by an animal carcass would make them less competitive in soils. Frankiales spp. can be plant symbionts, but have also been found to grow on rock beds, which indicate a lifestyle specialized for oligothrophic or extreme conditions [70, 71]. In that respect, it is reasonable that these oligothrophic taxa do not react dramatically to the nutrients of a decomposing animal.

Among the copiotrophic bacteria we see an increase in genera like Acinetobacter, Lysinibacillus and Kurthia in the first week, which is similar to observations in other decomposition studies [12]. These genera belong to the orders Pseudomonadales and Bacillales and are known to contain species that can be either pathogenic or beneficial for both animals and plants [72,73,74,75,76].

Interestingly, many of the orders reacting quickly to the nutrient influx remain at a relatively high level throughout the sampling period with variations in abundance of genera within orders (Fig. 4). For example, Xanthomonadales is abundant at relatively stable levels from day 3 onwards, but the genus Wohlfahrtiimonas within this order is only abundant on day 7 at Ca1. Wohlfahrtiimonas is a known parasite of different Wohlfahrtia spp. (flesh flies) and other such insects. The abundance of Wohlfahrtiimonas seen on day 7 for Ca1 may be a result of an increased number of insect hosts appearing through the decomposition process. This example indicates that short-term species composition changes might be due to local conditions in the soil and are not necessarily due to direct competition or sample variation.

Nonetheless, after day 7 other bacterial orders, such as Burkholderiales, increase in abundance. This order is typically found in higher abundances in the rhizosphere of plants than in the surrounding soils because they actively feed on plant root exudates [77,78,79]. Furthermore, Burkholderiales are known to have members that are plant growth-promoting bacteria (PGPB), like the Bacillales and Pseudomonadales genera described above [80,81,82]. Turner et al. [16] showed that decomposition of carcasses improves soil fertility and enhances plant growth. Our results suggest that this process is promoted by the soil microbial community, which shows a sharp increase in the abundance of orders known to have PGPBs. The promoting functions of PGPBs are many, such as nitrogen fixation, siderophore and phytohormone production, phosphorus solubilisation, and the suppression of plant diseases [79]. Bacillus anthracis is an obligate-killer pathogen, able to transmit only by killing its host. Interestingly, by taking advantage of the nutrient- and microbial-driven stimulation of plant growth occurring at carcass sites, transmission of B. anthracis is enhanced since this plant growth is an attractive food source for herbivores [16].

Microbial metabolism in arid environments

To get a functional overview of the microbial communities in our study we investigated how metabolic pathways correlate with changes in AGS over time (Fig. 3). Several of the positively correlating pathways can be linked to environmental stressors typical for arid soils, such as desiccation, high temperatures and (UV) radiation. For instance, starch and sucrose (S&S) metabolism is essential for the synthesis of the disaccharides sucrose and trehalose, which can play a role in prokaryotic desiccation resistance [83]. The importance of DNA repair in the surface soil exposed to high radiation levels is reflected in the positive correlation of the non-homologous end joining (NHEJ) repair system, which is also needed for survival during quiescent states such as in the spore stage [84,85,86].

The presence of other positively correlating pathways can be explained by their role in nutrient acquisition. Arid soils have low organic matter content and most of it comes from decaying plant material [56]. The positive correlation of the aminobenzoate degradation pathway probably reflects its involvement in the degradation of recalcitrant compounds such as lignin.

Many of the pathways negatively correlated with AGS are essential pathways involved in membrane transport, DNA metabolism (transcription, translation, replication, repair and recombination), amino-acid, lipid, carbohydrate and cofactor metabolism (Additional file 13: Figure S6). This reflects the dominance of opportunistic bacteria, especially at day 3–7, with small genomes that lack the arid soil ‘specialist’ pathways discussed above.

Overall, the analysis of the functional metabolism shows a clear change between the different time-points. The “normal” microbial soil community metabolism in ENP is represented by species that have genes that are involved in dealing with stress encountered in arid soils. The carcass fluids and nutrients activate many dormant or low abundant copiothrophic lineages present in the soil, which is reflected by the pathways needed for proliferation especially present in the first week of our experiment. These two states of the soil community, dormant vs. active, are typical for arid environments where water input stimulates in a pulse-like manner microbial activity and if sufficient even plant growth [56]. Rainfall gives enough moisture to activate both plant growth and the microbial decomposition compartment of the soil, providing plants with necessary nutrients [56, 87, 88]. However, in contrast to rainfall, animal carcasses not only provide moisture to the soil, but also additional nutrients. This pulse of nutrients and moisture allows the microbial community at very local scales to break down carcass nutrients, laying the foundation for enhanced plant growth directly and in the future [16, 56].

Conclusion

Our work is the first study to use shotgun metagenomics to describe the temporal changes in the soil microbial community after deposition of anthrax carcasses. We observed that the metagenomic approach is not the best way to study B. anthracis in a natural setting as it sporulates quickly making detection via DNA extraction difficult. In addition, it is too genetically similar to its close relatives such that most currently available tools to analyse (meta-) genomic sequences are not able to differentiate between such similar organisms. Nevertheless, the metagenomic data allows us to assess how this pathogen influences the microbial soil community.

The metagenomic time-series analysis indicated that carcass introduction stimulated the increase in abundance of part of the microbial soil community. This could play an indirect role in the future transmission of B. anthracis by stimulating the growth of taxa known to have PGPBs. The taxonomic shift in the microbial community was accompanied by shifts in the metabolic potential. The typical semi-arid soil community was involved in stress evasion, while the carcass introduction increased pathways involved in proliferation. Nonetheless, the time-series showed that within a month the semi-arid bacterial soil community was similar to the community at T = 0 with respect to taxonomic and metabolic composition.

Our work provides a background to study the ecology of spore-forming pathogens in a natural setting using culture-independent methods. It highlights the difficulty of using DNA based approaches to study spore-forming organisms such as B. anthracis in arid soils where they are dormant and poorly detectable.

Methods

Carcass and sample information

Sample collection in ENP was done after receiving permission from the Ministry of Environment and Tourism of Namibia (permit number: 1857/2013). Two fresh (likely <24 h since death) plains zebra (Equus quagga) anthrax carcasses (Ca1 & Ca2) were found on 03.03.2014 between Sprokieswoud and Charl Marais Dam, 835 m apart at coordinates S19.031/E015.548 (Ca1) and S19.037/E015.553 (Ca2), respectively. Information about ENP, carcass sites and sample collection can be found in the Additional file 4. Soil samples were taken from the area of blood-spill at six time-points (T), starting at T = 0 days (03.03.2014), T = 3, 7, 14, 21, 30 days. A control sample (Ctrl0) was also taken from the grids at T = 0 days, avoiding areas covered in blood (Additional file 8: Figure S4).

DNA isolation and purification

A total of 14 soil samples were collected (7 per carcass). DNA was isolated in two replicates from each sample (P1 and P2) resulting in a total of 14 DNA samples per carcass site. DNA was isolated using the FastDNA® spin kit for soil (MP Biomedicals, Santa Ana, California, USA), following the manufacturer’s protocol with adjustments specified in Additional file 4. The DNA samples were filter sterilised using an Ultrafree® Durapore PVDF 0.1 μM spinfilter (Millipore, Darmstadt, Germany).

DNA was shipped dry, following an ethanol precipitation, from ENP to Norway (see Additional file 4: for details). Upon arrival DNA was re-suspended in DES buffer from FastDNA® spin kit for soil. Samples were concentrated and purified using Agencourt® AMPure® XP beads (Beckman Coulter, Beverly, Massachusetts, USA), and treated with a PowerClean® Pro DNA clean-Up Kit, (MO BIO Laboratories, Carlsbad, California, USA) to remove any remaining inhibitors.

Metagenome sequencing and quality control

Environmental DNA samples were sequenced at the Norwegian Sequencing Centre (NSC) using Illumina MiSeq® (Illumina Inc., San Diego, California, USA). A 250 bp paired-end sequencing library (400 bp insert length) was generated using a Regular TruSeq® adapter ligation kit (Illumina Inc., San Diego, California, USA). Fastq files were processed using cut-adapt (v1.8) [89] and prinseq-lite (v0.20.4) [90] (See Additional file 4: for settings). Metagenomic read GC-content was calculated using infoseq (EMBOSS v. 6.5.7 [91]). The infoseq output was used to calculate average GC-content and GC-content skewness (timeData package version: 3012.100) using R-Studio (v3.3.0).

The average genome size (AGS) of each metagenome was estimated using MicrobeCensus on clean paired-end files (Additional file 4: Methods). The likelihood of sampling a universal single copy is inversely correlated with the AGS of a metagenomes [48, 92]. Thus universal single copy abundance will show significant differences when AGS between metagenomes differ.

Bacillus anthracis quantification from soils

Two methods were used to quantify B. anthracis in the soil samples; serial dilution culturing and quantitative PCR (qPCR) (see Additional file 4: for details).

The genome coverages of the B. anthracis isolate K1 and K2 genomes (Accession numbers: LBBZ00000000 and LBCA00000000 [46]) isolated from the Ca1 and Ca2 sites respectively, were investigated by mapping metagenomic sequences using the Burrows-Wheeler Aligner (BWA v0.7.8) [44]. The genomes of B. cereus E33L, B. thuringiensis HD-771 and Bacillus subtilis 168 were also used as references (see Additional file 4: for details).

Taxonomic and functional classification of metagenomic reads

Taxonomic composition of the metagenomes was determined using metaxa2 (v2.0.1) which extracts SSU rRNA sequences [50], MEGAN (v5.10.15) [54], Kraken (0.10.5-beta) [53] and metaBIT [51]. For the details on the comparison of the taxonomic classifiers see Additional file 4. The metaxa2 taxonomic profiles were visualised with a heatmap using ggplot2 in R-studio (see Additional file 4: for settings). Metaxa2 extracted SSU reads were run through the MetaAmp 1.1 pipeline to obtain operational taxonomical unit (OTU) abundances.

Shotgun sequences from each dataset were compared to the NCBI nr database (accessed on 20.01.2016) using Diamond (v0.7.11) [93] and the output was used for taxonomic and functional classification in MEGAN (v5.10.15) [54]. Temporal species composition dynamics of metagenomic samples was examined by principle coordinate analysis (PCoA) in R-Studio (v3.3.0). Changes in the metabolic potential of the community over time were analysed by correlation of KEGG KO-terms and pathways [94, 95] abundances with AGS size (For details see Additional file 4).

Data availability

The datasets supporting the conclusions of this article are available in the SRA repository, SRP076706. All sequence data was submitted to Genbank as part of bioproject: PRJNA281298.

References

  1. Fierer N, Leff JW, Adams BJ, Nielsen UN, Bates ST, Lauber CL, Owens S, Gilbert JA, Wall DH, Caporaso JG. Cross-biome metagenomic analyses of soil microbial communities and their functional attributes. Proc Natl Acad Sci U S A. 2012;109(52):21390–5.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  2. O'Brien SL, Gibbons SM, Owens SM, Hampton-Marcell J, Johnston ER, Jastrow JD, Gilbert JA, Meyer F, Antonopoulos DA. Spatial scale drives patterns in soil bacterial diversity. Environ Microbiol. 2016;18(6):2039–51.

    Article  PubMed  PubMed Central  Google Scholar 

  3. Rousk J, Bååth E, Brookes PC, Lauber CL, Lozupone C, Caporaso JG, Knight R, Fierer N. Soil bacterial and fungal communities across a pH gradient in an arable soil. ISME J. 2010;4(10):1340–51.

    Article  PubMed  Google Scholar 

  4. Bell CW, Tissue DT, Loik ME, Wallenstein MD. Acosta - Martinez V, Erickson RA, Zak JC: soil microbial and nutrient responses to 7 years of seasonally altered precipitation in a Chihuahuan Desert grassland. Glob Chang Biol. 2014;20(5):1657–73.

    Article  PubMed  Google Scholar 

  5. Johnson SL, Kuske CR, Carney TD, Housman DC, Gallegos-Graves LV, Belnap J. Increased temperature and altered summer precipitation have differential effects on biological soil crusts in a dryland ecosystem. Glob Chang Biol. 2012;18(8):2583–93.

    Article  Google Scholar 

  6. Garcia-Pichel F, Loza V, Marusenko Y, Mateo P, Potrafka RM. Temperature drives the continental-scale distribution of key microbes in topsoil communities. Science. 2013;340(6140):1574–7.

    Article  CAS  PubMed  Google Scholar 

  7. Niu F, He J, Zhang G, Liu X, Liu W, Dong M, Wu F, Liu Y, Ma X, An L, et al. Effects of enhanced UV-B radiation on the diversity and activity of soil microorganism of alpine meadow ecosystem in Qinghai–Tibet plateau. Ecotoxicology. 2014;23(10):1833–41.

    Article  CAS  PubMed  Google Scholar 

  8. Bell CW, Acosta-Martinez V, McIntyre NE, Cox S, Tissue DT, Zak JC. Linking microbial community structure and function to seasonal differences in soil moisture and temperature in a Chihuahuan desert grassland. Microb Ecol. 2009;58(4):827–42.

    Article  CAS  PubMed  Google Scholar 

  9. Bi J, Zhang NL, Liang Y, Yang HJ, Ma KP. Interactive effects of water and nitrogen addition on soil microbial communities in a semiarid steppe. J Plant Ecol. 2012;5(3):320–9.

    Article  Google Scholar 

  10. Rath KM, Maheshwari A, Rousk J. The impact of salinity on the microbial response to drying and rewetting in soil. Soil Biol Biochem. 2017;108:17–26.

    Article  CAS  Google Scholar 

  11. Van Horn DJ, Okie JG, Buelow HN, Gooseff MN, Barrett JE, Takacs-Vesbach CD. Soil microbial responses to increased moisture and organic resources along a salinity gradient in a polar desert. Appl Environ Microbiol. 2014;80(10):3034–43.

    Article  PubMed  PubMed Central  Google Scholar 

  12. Howard GT, Duos B, Watson-Horzelski EJ. Characterization of the soil microbial community associated with the decomposition of a swine carcass. Int Biodeterior Biodegradation. 2010;64(4):300–4.

    Article  Google Scholar 

  13. Cobaugh KL, Schaeffer SM, DeBruyn JM. Functional and structural succession of soil microbial communities below decomposing human cadavers. PLoS One. 2015;10(6):e0130201.

    Article  PubMed  PubMed Central  Google Scholar 

  14. Weiss S, Carter DO, Metcalf JL, Knight R. Carcass mass has little influence on the structure of gravesoil microbial communities. Int J Legal Med. 2016;130(1):253–63.

    Article  PubMed  Google Scholar 

  15. Bump JK, Peterson RO, Vucetich JA. Wolves modulate soil nutrient heterogeneity and foliar nitrogen by configuring the distribution of ungulate carcasses. Ecology. 2009;90(11):3159–67.

    Article  PubMed  Google Scholar 

  16. Turner WC, Kausrud KL, Krishnappa YS, Cromsigt JP, Ganz HH, Mapaure I, Cloete CC, Havarua Z, Küsters M, Getz WM. Fatal attraction: vegetation responses to nutrient inputs attract herbivores to infectious anthrax carcass sites. Proc R Soc Lond B Biol Sci. 2014;281(1795):20141785.

    Article  Google Scholar 

  17. Sharp RJ, Roberts AG. Anthrax: the challenges for decontamination. J Chem Technol Biotechnol. 2006;81(10):1612–25.

    Article  CAS  Google Scholar 

  18. Liu Y, Lai Q, Göker M, Meier-Kolthoff JP, Wang M, Sun Y, Wang L, Shao Z. Genomic insights into the taxonomic status of the Bacillus Cereus group. Sci Rep. 2015;5:14082.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  19. Jensen GB, Hansen BM, Eilenberg J, Mahillon J. The hidden lifestyles of Bacillus cereus and relatives. Environ Microbiol. 2003;5(8):631–40.

    Article  CAS  PubMed  Google Scholar 

  20. Helgason E, Økstad OA, Caugant DA, Johansen HA, Fouet A, Mock M, Hegna I, Kolstø A-B. Bacillus anthracis, Bacillus cereus, and Bacillus thuringiensis—one species on the basis of genetic evidence. Appl Environ Microbiol. 2000;66(6):2627–30.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  21. Zwick ME, Joseph SJ, Didelot X, Chen PE, Bishop-Lilly KA, Stewart AC, Willner K, Nolan N, Lentz S, Thomason MK. Genomic characterization of the Bacillus cereus sensu lato species: backdrop to the evolution of Bacillus anthracis. Genome Res. 2012;22(8):1512–24.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  22. Mignot T, Mock M, Robichon D, Landier A, Lereclus D, Fouet A. The incompatibility between the PlcR-and AtxA-controlled regulons may have selected a nonsense mutation in Bacillus anthracis. Mol Microbiol. 2001;42(5):1189–98.

    Article  CAS  PubMed  Google Scholar 

  23. Easterday WR, Ert MN, Simonson TS, Wagner DM, Kenefic LJ, Allender CJ, Keim P. Use of single nucleotide polymorphisms in the plcR gene for specific identification of Bacillus anthracis. J Clin Microbiol. 2005;43:1995–7.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  24. Green BD, Battisti L, Koehler TM, Thorne CB, Ivins BE. Demonstration of a capsule plasmid in Bacillus anthracis. Infect Immun. 1985;49(2):291–7.

    CAS  PubMed  PubMed Central  Google Scholar 

  25. Mikesell P, Ivins BE, Ristroph JD, Dreier TM. Evidence for plasmid-mediated toxin production in Bacillus anthracis. Infect Immun. 1983;39(1):371–6.

    CAS  PubMed  PubMed Central  Google Scholar 

  26. Toby IT, Widmer J, Dyer DW. Divergence of protein-coding capacity and regulation in the Bacillus cereus sensu lato group. BMC Bioinformatics. 2014;15(11):1.

    CAS  Google Scholar 

  27. Turner WC, Kausrud KL, Beyer W, Easterday WR, Barandongo ZR, Blaschke E, Cloete CC, Lazak J, Van Ert MN, Ganz HH. Lethal exposure: an integrated approach to pathogen transmission via environmental reservoirs. Sci Rep. 2016;6:27311.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  28. Dragon DC, Bader DE, Mitchell J, Woollen N. Natural dissemination of Bacillus anthracis spores in northern Canada. Appl Environ Microbiol. 2005;71(3):1610–5.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  29. Bellan SE, Turnbull PC, Beyer W, Getz WM. Effects of experimental exclusion of scavengers from carcasses of anthrax-infected herbivores on Bacillus anthracis sporulation, survival, and distribution. Appl Environ Microbiol. 2013;79(12):3756–61.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  30. Ganz HH, Turner WC, Brodie EL, Kusters M, Shi Y, Sibanda H, Torok T, Getz WM. Interactions between Bacillus anthracis and plants may promote anthrax transmission. PLoS Negl Trop Dis. 2014;8(6):e2903.

    Article  PubMed  PubMed Central  Google Scholar 

  31. Schuch R, Fischetti VA. The secret life of the anthrax agent Bacillus anthracis: bacteriophage-mediated ecological adaptations. PLoS One. 2009;4(8):e6532.

    Article  PubMed  PubMed Central  Google Scholar 

  32. Saile E, Koehler TM. Bacillus anthracis multiplication, persistence, and genetic exchange in the rhizosphere of grass plants. Appl Environ Microbiol. 2006;72(5):3168–74.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  33. Dey R, Hoffman PS, Glomski IJ. Germination and amplification of anthrax spores by soil-dwelling amoebas. Appl Environ Microbiol. 2012;78(22):8075–81.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  34. Logares R, Haverkamp TH, Kumar S, Lanzén A, Nederbragt AJ, Quince C, Kauserud H. Environmental microbiology through the lens of high-throughput DNA sequencing: synopsis of current platforms and bioinformatics approaches. J Microbiol Methods. 2012;91(1):106–13.

    Article  CAS  PubMed  Google Scholar 

  35. Be NA, Thissen JB, Gardner SN, McLoughlin KS, Fofanov VY, Koshinsky H, Ellingson SR, Brettin TS, Jackson PJ, Jaing CJ. Detection of bacillus anthracis DNA in complex soil and air samples using next-generation sequencing. PLoS One. 2013;8(9):e73455.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  36. Dragon DC, Rennie RP. Evaluation of spore extraction and purification methods for selective recovery of viable Bacillus anthracis spores. Lett Appl Microbiol. 2001;33(2):100–5.

    Article  CAS  PubMed  Google Scholar 

  37. Le Roux C, Grunow J, Morris J, Bredenkamp G, Scheepers J. A classification of the vegetation of the Etosha national park. S Afr J Bot. 1988;54(1):1–10.

    Article  Google Scholar 

  38. Engert S. Spatial variability and temporal periodicity of rainfall in the Etosha National Park and surrounding areas in northern Namibia. Modoqua. 1997;20(1):115–20.

    Google Scholar 

  39. Lauber CL, Metcalf JL, Keepers K, Ackermann G, Carter DO, Knight R. Vertebrate decomposition is accelerated by soil microbes. Appl Environ Microbiol. 2014;80(16):4920–9.

    Article  PubMed  PubMed Central  Google Scholar 

  40. Bellan SE, Gimenez O, Choquet R, Getz WM. A hierarchical distance sampling approach to estimating mortality rates from opportunistic carcass surveillance data. Methods Ecol Evol. 2013;4(4):361–9.

    Article  Google Scholar 

  41. Hyde ER, Haarmann DP, Petrosino JF, Lynne AM, Bucheli SR. Initial insights into bacterial succession during human decomposition. Int J Legal Med. 2015;129(3):661–71.

    Article  PubMed  Google Scholar 

  42. Lauber CL, Hamady M, Knight R, Fierer N. Pyrosequencing-based assessment of soil pH as a predictor of soil bacterial community structure at the continental scale. Appl Environ Microbiol. 2009;75(15):5111–20.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  43. Albertsen M, Karst SM, Ziegler AS, Kirkegaard RH, Nielsen PH. Back to basics–the influence of DNA extraction and primer choice on phylogenetic analysis of activated sludge communities. PLoS One. 2015;10(7):e0132783.

    Article  PubMed  PubMed Central  Google Scholar 

  44. Li H, Durbin R. Fast and accurate long-read alignment with burrows–wheeler transform. Bioinformatics. 2010;26(5):589–95.

    Article  PubMed  PubMed Central  Google Scholar 

  45. Rasmussen S, Allentoft ME, Nielsen K, Orlando L, Sikora M, Sjögren K-G, Pedersen AG, Schubert M, Van Dam A, Kapel CMO. Early divergent strains of Yersinia pestis in Eurasia 5,000 years ago. Cell. 2015;163(3):571–82.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  46. Valseth K, Nesbø CL, Easterday WR, Turner WC, Olsen JS, Stenseth NC, Haverkamp THA. Draft genome sequences of two bacillus anthracis strains from Etosha National Park, Namibia. Genome Announc. 2016;4(4):e00861–16.

    Article  PubMed  PubMed Central  Google Scholar 

  47. Reichenberger ER, Rosen G, Hershberg U, Hershberg R. Prokaryotic nucleotide composition is shaped by both phylogeny and the environment. Genome Biol Evol. 2015;7(5):1380–9.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  48. Nayfach S, Pollard KS. Average genome size estimation improves comparative metagenomics and sheds light on the functional ecology of the human microbiome. Genome Biol. 2015;16(1):51.

    Article  PubMed  PubMed Central  Google Scholar 

  49. Raes J, Korbel JO, Lercher MJ, von Mering C, Bork P. Prediction of effective genome size in metagenomic samples. Genome Biol. 2007;8(1):R10.

    Article  PubMed  PubMed Central  Google Scholar 

  50. Bengtsson-Palme J, Hartmann M, Eriksson KM, Pal C, Thorell K, Larsson DGJ, Nilsson RH. Metaxa2: improved identification and taxonomic classification of small and large subunit rRNA in metagenomic data. Mol Ecol Resour. 2015;15(6):1403–14.

    Article  CAS  PubMed  Google Scholar 

  51. Louvel G, Der Sarkissian C, Hanghøj K, Orlando L. metaBIT, an integrative and automated metagenomic pipeline for analysing microbial profiles from high-throughput sequencing shotgun data. Mol Ecol Resour. 2016;16(6):1415–27.

    Article  CAS  PubMed  Google Scholar 

  52. Truong DT, Franzosa EA, Tickle TL, Scholz M, Weingart G, Pasolli E, Tett A, Huttenhower C, Segata N. MetaPhlAn2 for enhanced metagenomic taxonomic profiling. Nat Methods. 2015;12(10):902–3.

    Article  CAS  PubMed  Google Scholar 

  53. Wood DE, Salzberg SL. Kraken: ultrafast metagenomic sequence classification using exact alignments. Genome Biol. 2014;15(3):R46.

    Article  PubMed  PubMed Central  Google Scholar 

  54. Huson DH, Auch AF, Qi J, Schuster SC. MEGAN analysis of metagenomic data. Genome Res. 2007;17(3):377–86.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  55. Lindgreen S, Adair KL, Gardner PP. An evaluation of the accuracy and speed of metagenome analysis tools. Sci Rep. 2016;6:19233.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  56. Collins SL, Sinsabaugh RL, Crenshaw C, Green L, Porras-Alfaro A, Stursova M, Zeglin LH. Pulse dynamics and microbial processes in aridland ecosystems. J Ecol. 2008;96(3):413–20.

    Article  Google Scholar 

  57. Dorsch M, Lovet D, Bailey GD. Fusobacterium equinum sp. nov., from the oral cavity of horses. Int J Syst Evol Microbiol. 2001;51(6):1959–63.

    Article  CAS  PubMed  Google Scholar 

  58. Schlater LK, Brenner D, Steigerwalt A, Moss CW, Lambert M, Packer R. Pasteurella caballi, a new species from equine clinical specimens. J Clin Microbiol. 1989;27(10):2169–74.

    CAS  PubMed  PubMed Central  Google Scholar 

  59. Konstantinidis KT, Tiedje JM. Trends between gene content and genome size in prokaryotic species with larger genomes. Proc Natl Acad Sci U S A. 2004;101(9):3160–5.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  60. Prosser JI. Dispersing misconceptions and identifying opportunities for the use of 'omics' in soil microbial ecology. Nat Rev Microbiol. 2015;13(7):439–46.

    Article  CAS  PubMed  Google Scholar 

  61. Teeling H, Glöckner FO. Current opportunities and challenges in microbial metagenome analysis--a bioinformatic perspective. Brief Bioinform. 2012;13(6):728–42.

    Article  PubMed  PubMed Central  Google Scholar 

  62. Fierer N, Bradford MA, Jackson RB. Towards an ecological classification of soil bacteria. Ecology. 2007;88(6):1354–64.

    Article  PubMed  Google Scholar 

  63. Tan J, Zuniga C, Zengler K. Unraveling interactions in microbial communities - from co-cultures to microbiomes. J Microbiol. 2015;53(5):295–305.

    Article  PubMed  Google Scholar 

  64. Staley JT, Konopka A. Measurement of in situ activities of nonphotosynthetic microorganisms in aquatic and terrestrial habitats. Annu Rev Microbiol. 1985;39(1):321–46.

    Article  CAS  PubMed  Google Scholar 

  65. Torsvik V, Øvreås L. Microbial diversity and function in soil: from genes to ecosystems. Curr Opin Microbiol. 2002;5(3):240–5.

    Article  CAS  PubMed  Google Scholar 

  66. Lindeque PM, Turnbull PCB. Ecology and epidemiology of anthrax in Etosha National Park, Namibia. Onderstepoort J Vet Res. 1994;61(1):71–83.

    CAS  PubMed  Google Scholar 

  67. Sedlackova V, Dziedzinska R, Babak V, Kralik P. The detection and quantification of Bacillus thuringiensis spores from soil and swabs using quantitative PCR as a model system for routine diagnostics of Bacillus anthracis. J Appl Microbiol. 2017; https://0-doi-org.brum.beds.ac.uk/10.1111/jam.13445.

  68. Ketola T, Mikonranta L, Laakso J, Mappes J. Different food sources elicit fast changes to bacterial virulence. Biol Lett. 2016;12(1):20150660.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  69. Erlacher A, Cernava T, Cardinale M, Soh J, Sensen CW, Grube M, Berg G. Rhizobiales as functional and endosymbiontic members in the lichen symbiosis of Lobaria pulmonaria L. Front Microbiol. 2015;6:53.

    Article  PubMed  PubMed Central  Google Scholar 

  70. Goodfellow M, Kämpfer P, Busse H-J, Trujillo ME, Suzuki K-I, Ludwig W, Whitman WB. Bergey's manual of systematic bacteriology : Volume 5: The Actinobacteria. In: vol. 5, 2 edn. New York: Springer Science & Business Media; 2012: 509–546.

  71. Foesel B, Geppert A, Rohde M, Overmann J. Parviterribacter kavangonensis and Parviterribacter multiflagellatus a novel genus and two novel species within the order Solirubrobacterales and emended description of the classes Thermoleophilia and Rubrobacteria and its orders and families. Int J Syst Evol Microbiol. 2016;66:652–65.

    Article  CAS  Google Scholar 

  72. Weber BS, Harding CM, Feldman MF. Pathogenic Acinetobacter: from the cell surface to infinity and beyond. J Bacteriol. 2015;198(6):880–7.

    Article  PubMed  Google Scholar 

  73. Berry C. The bacterium, Lysinibacillus sphaericus, as an insect pathogen. J Invertebr Pathol. 2012;109(1):1–10.

    Article  PubMed  Google Scholar 

  74. Vos P, Garrity G, Jones D, Krieg NR, Ludwig W, Rainey FA, Schleifer K-H, Whitman W. Bergey's manual of systematic bacteriology: Volume 3: The Firmicutes. In: vol. 3, 2 edn. New York: Springer Science & Business Media; 2011: 364–370.

  75. Andrade LF, de Souza G, Nietsche S, Xavier AA, Costa MR, Cardoso AMS, Pereira MCT, Pereira D. Analysis of the abilities of endophytic bacteria associated with banana tree roots to promote plant growth. J Microbiol. 2014;52(1):27–34.

    Article  CAS  PubMed  Google Scholar 

  76. Batool F, Rehman Y, Hasnain S. Phylloplane associated plant bacteria of commercially superior wheat varieties exhibit superior plant growth promoting abilities. Front Life Sci. 2016;9(4):313–22.

    Article  CAS  Google Scholar 

  77. Uroz S, Buée M, Murat C, Frey-Klett P, Martin F. Pyrosequencing reveals a contrasted bacterial diversity between oak rhizosphere and surrounding soil. Environ Microbiol Rep. 2010;2(2):281–8.

    Article  CAS  PubMed  Google Scholar 

  78. Vandenkoornhuyse P, Mahé S, Ineson P, Staddon P, Ostle N, Cliquet JB, Francez AJ, Fitter AH, Young JP. Active root-inhabiting microbes identified by rapid incorporation of plant-derived carbon into RNA. Proc Natl Acad Sci U S A. 2007;104(43):16970–5.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  79. Bulgarelli D, Schlaeppi K, Spaepen S, Ver Loren van Themaat E, Schulze-Lefert P. Structure and functions of the bacterial microbiota of plants. Annu Rev Plant Biol. 2013;64:807–38.

    Article  CAS  PubMed  Google Scholar 

  80. Baldani JI, Rouws L, Cruz LM, Olivares FL, Schmid M, Hartmann A. The family Oxalobacteraceae. In: The prokaryotes: Alphaproteobacteria and Betaproteobacteria. Vol. 9783642301971; 2014. p. 919–74.

    Chapter  Google Scholar 

  81. Pieterse CM, Zamioudis C, Berendsen RL, Weller DM, Van Wees SC, Bakker PA. Induced systemic resistance by beneficial microbes. Annu Rev Phytopathol. 2014;52:347–75.

    Article  CAS  PubMed  Google Scholar 

  82. Choudhary DK, Johri BN. Interactions of Bacillus spp. and plants--with special reference to induced systemic resistance (ISR). Microbiol Res. 2009;164(5):493–513.

    Article  CAS  PubMed  Google Scholar 

  83. Potts M. Desiccation tolerance of prokaryotes. Microbiol Rev. 1994;58(4):755–805.

    CAS  PubMed  PubMed Central  Google Scholar 

  84. Shuman S, Glickman MS. Bacterial DNA repair by non-homologous end joining. Nat Rev Microbiol. 2007;5(11):852–61.

    Article  CAS  PubMed  Google Scholar 

  85. Moeller R, Stackebrandt E, Reitz G, Berger T, Rettberg P, Doherty AJ, Horneck G, Nicholson WL. Role of DNA repair by nonhomologous-end joining in Bacillus subtilis spore resistance to extreme dryness, mono-and polychromatic UV, and ionizing radiation. J Bacteriol. 2007;189(8):3306–11.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  86. Garcia-Gonzalez A, Vicens L, Alicea M, Massey S. The distribution of recombination repair genes is linked to information content in bacteria. Gene. 2013;528(2):295–303.

    Article  CAS  PubMed  Google Scholar 

  87. Du Plessis W. Effective rainfall defined using measurements of grass growth in the Etosha National Park, Namibia. J Arid Environ. 2001;48(3):397–417.

    Article  Google Scholar 

  88. Clark JS, Campbell JH, Grizzle H, Acosta-Martìnez V, Zak JC. Soil microbial community response to drought and precipitation variability in the Chihuahuan Desert. Microb Ecol. 2009;57(2):248–60.

    Article  PubMed  Google Scholar 

  89. Martin M. Cutadapt removes adapter sequences from high-throughput sequencing reads. EMBnet journal. 2011;17(1):10–2.

    Article  Google Scholar 

  90. Schmieder R, Edwards R. Quality control and preprocessing of metagenomic datasets. Bioinformatics. 2011;27(6):863–4.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  91. Rice P, Longden I, Bleasby A. EMBOSS: the European molecular biology open software suite. Trends Genet. 2000;16(6):276–7.

    Article  CAS  PubMed  Google Scholar 

  92. Beszteri B, Temperton B, Frickenhaus S, Giovannoni SJ. Average genome size: a potential source of bias in comparative metagenomics. ISME J. 2010;4(8):1075–7.

    Article  PubMed  Google Scholar 

  93. Buchfink B, Xie C, Huson DH. Fast and sensitive protein alignment using DIAMOND. Nat Methods. 2015;12(1):59–60.

    Article  CAS  PubMed  Google Scholar 

  94. Kanehisa M, Goto S. KEGG: kyoto encyclopedia of genes and genomes. Nucleic Acids Res. 2000;28(1):27–30.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  95. Kanehisa M, Sato Y, Kawashima M, Furumichi M, Tanabe M. KEGG as a reference resource for gene and protein annotation. Nucleic Acids Res. 2016;44(D1):D457–62.

    Article  CAS  PubMed  Google Scholar 

  96. Love MI, Huber W, Anders S. Moderated estimation of fold change and dispersion for RNA-seq data with DESeq2. Genome Biol. 2014;15(12):550.

    Article  PubMed  PubMed Central  Google Scholar 

Download references

Acknowledgements

We thank Namibia’s Ministry of Environment and Tourism for permission to conduct this research (permit number: 1857/2013). In addition, we thank Dr. Mari Espelund for help optimizing DNA extractions, Claudine C. Cloete and Zoe Barandongo for help with sample collecting and DNA extractions, the Etosha Ecological Institute and the Norwegian Defence Research Institute for providing laboratory facilities and the Department of Biosciences student funding at University of Oslo for travel money. We also thank Eric de Muinck and Pål Trosvik for critical feedback on our manuscript.

Funding

This work was supported by the Research Council of Norway (Project No. 180444/V40) awarded to CLN and NSF OISE-1103054 awarded to WCT.

Availability of data and materials

The datasets supporting the conclusions of this article are available in the SRA repository, SRP076706. All sequence data was submitted to Genbank as part of bioproject: PRJNA281298.

Additional file 4 accompanies the paper on the BMC Microbiology journal website.

Author information

Authors and Affiliations

  1. Department of Biosciences, Centre for Ecological and Evolutionary Synthesis (CEES), University of Oslo, The Kristine Bonnevie Building, UiO, campus Blindern, Blindern, Oslo, Norway

    Karoline Valseth, Camilla L. Nesbø, W. Ryan Easterday, Nils Chr. Stenseth & Thomas H. A. Haverkamp

  2. Norwegian Defence Research Establishment, Kjeller, Norway

    Karoline Valseth & Jaran S. Olsen

  3. Department of Biological Sciences, University of Alberta, Edmonton, Alberta, Canada

    Camilla L. Nesbø

  4. Department of Biological Sciences, University at Albany, State University of New York, Albany, New York, USA

    Wendy C. Turner

Authors
  1. Karoline Valseth
    View author publications

    You can also search for this author in PubMed Google Scholar

  2. Camilla L. Nesbø
    View author publications

    You can also search for this author in PubMed Google Scholar

  3. W. Ryan Easterday
    View author publications

    You can also search for this author in PubMed Google Scholar

  4. Wendy C. Turner
    View author publications

    You can also search for this author in PubMed Google Scholar

  5. Jaran S. Olsen
    View author publications

    You can also search for this author in PubMed Google Scholar

  6. Nils Chr. Stenseth
    View author publications

    You can also search for this author in PubMed Google Scholar

  7. Thomas H. A. Haverkamp
    View author publications

    You can also search for this author in PubMed Google Scholar

Contributions

All authors contributed to the design of the experiment. KV, WRE, WCT conducted the environmental sampling and qPCR experiments, bacterial cultivation and DNA extractions. KV prepared samples for metagenomic sequencing and performed bioinformatic quality control. KV, CLN and THA executed bioinformatics analysis of shotgun data. KV prepared the manuscript. All authors contributed to critical revision of the manuscript. All authors read and approved the final manuscript.

Corresponding author

Correspondence to Thomas H. A. Haverkamp.

Ethics declarations

Ethics approval and consent to participate

Not applicable.

Competing interests

The authors declare that they have no competing interests.

Publisher’s Note

Springer Nature remains neutral with regard to jurisdictional claims in published maps and institutional affiliations.

Additional files

Additional file 1: Table S1.

Soil nutrients and soil composition at carcass sites. (DOCX 44 kb)

Additional file 2: Figure S1.

Abundance and relative abundance of 16S/18S rRNA taxonomically classified into kingdoms. (a) Raw abundance and (b) relative abundance of metaxa2 classified reads. These are the 16S/18S reads extracted and classified by metaxa2. Roughly 50% of the extracted rRNA sequenced did not get assigned to any taxonomy and were assumed to be false positives. (EPS 1341 kb)

Additional file 3: Figure S2.

Rainfall at the Okaukuejo field station. The precipitation is measured in mm at Okaukuejo, which is approximately 40 km from the sampling site. Because of the very local rainfalls in the area, this can only be used as an indication of rainfall at our sample sites. The graph shows rainfall from the week before sampling started and to the end of the sampling period, the 2nd of April 2014. The blue arrow indicates the sampling period with sampling days marked in red. (EPS 7 kb)

Additional file 4:

A text document contained supplementary methods, equations, results and figure (mentioned only in the supplementary information). (DOCX 896 kb)

Additional file 5: Table S2.

Counts of mapped metagenomic reads of each sample against reference genomes. (XLSX 97 kb)

Additional file 6: Figure S3.

Estimation of B. anthracis abundance in soil samples using qPCR, metagenomic reads mapping and cultivation. For all panels, sampling time-points are on the x-axis. (A and D) qPCR results with estimated number of B. anthracis genomes per gram soil along the y-axis at different time-points (x-axis) for Carcass 1 (Ca1) (A) and Carcass 2 (Ca2) (D). The spiked sample is control soil from day 0 with 3.4 × 106 Stern34F2 B. anthracis spores added. P1 and P2 represent the two replicates at each time-point. (B, C, E and F) shows the results from mapping the metagenomic reads against the B. anthracis isolate genomes (K1/K2) isolated from Ca1 and Ca2 as well as the other reference strains; B. cereus E33L, B. thuringiensis HD-771 and B. subtilis 168. Reads with mismatches were not filtered away. Reads matching the references are on the y-axis. (B and E) and (C and F) shows the mapping of the metagenomes against the K1 and K2 strain, respectively as well as the other reference strains. Ca1 in panes (B and C) and Ca2 in panes (E and F). (G) B. anthracis spore counts of soil samples from Ca1 and Ca2 in spores per gram dry soil. Spore counts were estimated by culturing of heat-shocked soil samples and adjusted to dry weight of soil. (EPS 1226 kb)

Additional file 7: Table S3.

Comparison of taxonomic classification with four different methods: MetaBIT, Metaxa2, MEGAN, and Kraken. (XLSX 66 kb)

Additional file 8: Figure S4.

Comparison of the percentage of classified reads by four metagenomic classification tools. Note the very different scales on the y-axis. Barplots show the relative abundances of classified reads for each time point metagenome of Ca1 (blue) and Ca 2 (orange). The tools are Metaxa2 (A), metaBIT (B), MEGAN (C), Kraken (D). (EPS 1139 kb)

Additional file 9: Table S4.

Metaxa2 taxonomic classifications on pair1 16S /18S rRNA sequences. (XLSX 325 kb)

Additional file 10: Figure S5.

Soil microbial community relationships by time-point and carcass site. Principal coordinate analysis of a distance matrix created by normalised total counts and using Bray-Curtis dissimilarity. Ca1 is represented by the circles and Ca2 by the triangles, the time-points are visualised in different colours and the arrows are pointing in the direction of increasing days, red arrows for Ca1 and light blue for Ca2. The dotted arrows show the relationship between the day 30 to the Ctrl0 sample. The black arrows indicate the 10 most abundant orders and their effect on the soil microbial communities per time-point. Arrows were fitted using envfit with 1000 permutations using the vegan package in R-studio. (EPS 8 kb)

Additional file 11: Table S5.

Overview of the Spearman correlation coefficients of all identified KEGG pathways and their FDR corrected p-values. (XLSX 18 kb)

Additional file 12: Table S6.

KEGG KO-term abundances for significantly correlating KEGG pathways. Abundances are normalized using DESEQ2 and sorted on carcass 1, time point T = 0. (XLSX 368 kb)

Additional file 13: Figure S6.

Sample area. (a) Ca1, soil samples were taken from within the 30 × 30 cm grid, (b) Ca2, samples were taken from within the red square (resembling the 30 × 30 cm metal grid shown in (a)). (PNG 749 kb)

Rights and permissions

Open Access This article is distributed under the terms of the Creative Commons Attribution 4.0 International License (http://creativecommons.org/licenses/by/4.0/), which permits unrestricted use, distribution, and reproduction in any medium, provided you give appropriate credit to the original author(s) and the source, provide a link to the Creative Commons license, and indicate if changes were made. The Creative Commons Public Domain Dedication waiver (http://creativecommons.org/publicdomain/zero/1.0/) applies to the data made available in this article, unless otherwise stated.

Reprints and permissions

About this article

Check for updates. Verify currency and authenticity via CrossMark

Cite this article

Valseth, K., Nesbø, C.L., Easterday, W.R. et al. Temporal dynamics in microbial soil communities at anthrax carcass sites. BMC Microbiol 17, 206 (2017). https://0-doi-org.brum.beds.ac.uk/10.1186/s12866-017-1111-6

Download citation

  • Received:

  • Accepted:

  • Published:

  • DOI: https://0-doi-org.brum.beds.ac.uk/10.1186/s12866-017-1111-6

Keywords

BMC Microbiology

ISSN: 1471-2180

Contact us