Skip to main content
Download PDF

Adaptive modification of antiviral defense systems in microbial community under Cr-induced stress

Microbiome volume 13, Article number: 34 (2025) Cite this article

Abstract

Background

The prokaryotic antiviral defense systems are crucial for mediating prokaryote-virus interactions that influence microbiome functioning and evolutionary dynamics. Despite the prevalence and significance of prokaryotic antiviral defense systems, their responses to abiotic stress and ecological consequences remain poorly understood in soil ecosystems. We established microcosm systems with varying concentrations of hexavalent chromium (Cr(VI)) to investigate the adaptive modifications of prokaryotic antiviral defense systems under abiotic stress.

Results

Utilizing hybrid metagenomic assembly with long-read and short-read sequencing, we discovered that antiviral defense systems were more diverse and prevalent in heavily polluted soils, which was corroborated by meta-analyses of public datasets from various heavy metal-contaminated sites. As the Cr(VI) concentration increased, prokaryotes with defense systems favoring prokaryote-virus mutualism gradually supplanted those with defense systems incurring high adaptive costs. Additionally, as Cr(VI) concentrations increased, enriched antiviral defense systems exhibited synchronization with microbial heavy metal resistance genes. Furthermore, the proportion of antiviral defense systems carried by mobile genetic elements (MGEs), including plasmids and viruses, increased by approximately 43% and 39%, respectively, with rising Cr concentrations. This trend is conducive to strengthening the dissemination and sharing of defense resources within microbial communities.

Conclusions

Overall, our study reveals the adaptive modification of prokaryotic antiviral defense systems in soil ecosystems under abiotic stress, as well as their positive contributions to establishing prokaryote-virus mutualism and the evolution of microbial heavy metal resistance. These findings advance our understanding of microbial adaptation in stressful environments and may inspire novel approaches for microbiome manipulation and bioremediation.

Video Abstract

Introduction

Prokaryote-virus interactions can significantly influence the adaptive evolutionary and ecological processes of microbial communities [1]. Among them, antiviral defense systems carried by prokaryotes play a crucial role in countering virus infections, making them essential mediators of such interactions [2]. These systems target various stages of virus infection, from preventing virus attachment to inhibiting virus replication and release [3]. Previous studies have revealed the prevalence and diversity of prokaryotic antiviral defense systems, the underlying molecular mechanisms, and the biotic evolutionary drivers [4, 5]. Over a hundred antiviral defense systems have been identified in the genomes of cultivated and uncultivated bacteria and archaea [6,7,8]. However, the response of these prokaryotic antiviral defense systems to environmental disturbances remains largely unexplored in soil ecosystems, representing a crucial knowledge gap for illuminating the adaptation strategies of microbial communities.

Constrained by biotic and abiotic factors, microbial communities often exhibit trade-offs in physiological traits [9]. Studies on model bacteria have reported significant fitness costs associated with developing resistance to viruses and contaminants (e.g., slowed cell growth rates, decreased virulence, and impaired bacterial motility [10, 11]). Moreover, under abiotic stress conditions, reductions in biomass inevitably decrease encounters between viruses and prokaryotes, leading to a concurrent decrease in virus predation pressure on prokaryotes [12]. Thus, the resource trade-offs associated with increased abiotic resistance and reduced virus encounters may result in the loss or down-regulation of antiviral defense systems in polluted environments. On the other hand, changes in host cell physiology induced by pollution stress are highly likely to trigger prophage induction, potentially leading to cell death. Reinforcing defense may be necessary to prevent this scenario. At the community level, abiotic stress also drives changes in microbial composition [12], which may further influence the diversity and distribution of antiviral defense systems. Therefore, abiotic stress may reshape the prokaryotic antiviral defense systems, ultimately favoring microbial environmental adaptation to challenging conditions. Due to differences in defense mechanisms, the defense spectra and fitness costs of different antiviral defense systems vary considerably [2], resulting in divergent responses to abiotic stress and ecological consequences.

Several studies have revealed that elevated abiotic stress (e.g., pesticides [13] and arsenic (As) [14]) can lead to enhanced mutualism between the prokaryotes and their infecting viruses in soil. This mutualistic transition involves prokaryotes providing stable intracellular environments for viruses, while viruses significantly influence the adaptive evolution of their hosts by facilitating the transfer and expression of relevant resistance genes [15]. Chromium (Cr) is a highly toxic and persistent heavy metal pollutant widely distributed in industrial and mining areas. It poses significant environmental risks due to its high solubility and mobility potential, which can lead to severe soil and groundwater contamination [16]. Cr(VI) exposure disrupts microbial communities by inducing oxidative stress, reducing biomass, and altering ecological functions [16]. Therefore, understanding microbiomes adaptation to Cr(VI) is crucial for advancing bioremediation strategies and ecological risk assessment. Our previous research has demonstrated the adaptive mutualistic response between viruses and prokaryotes under Cr(VI) contamination [12]. However, our understanding of how prokaryotic defense systems are modified to maintain this mutualism with viruses under adverse habitats remains limited.

To address this gap, we established microcosm systems with varying Cr(VI) concentrations (0, 10, and 100 mg/kg, respectively), and inoculated microbial communities from both clean (Cr(VI) concentration: 0.094 mg/kg) and polluted (Cr(VI) concentration: 665 mg/kg) soil on site. The response of microbial defense systems to abiotic stress and its ecological consequences was investigated using hybrid metagenomic assembly, combining long-read and short-read sequencing, along with meta-analyses of public datasets from heavy metal-contaminated soils (Fig. 1). We investigated the modification of antiviral defense systems in terms of diversity, distribution, composition, and transcriptional expression in soil ecosystems under increasing and decreasing Cr(VI) contamination. The prokaryotic antiviral defense systems were explored regarding their relationship with viruses and their synchronicity with the resistance to heavy metals. Furthermore, we revealed the distribution of antiviral defense systems in plasmids and viruses and their contributions to prokaryotic immunity in stressful environments. Overall, this study advances our understanding of the complex dynamics and ecological significance of prokaryotic antiviral defense systems within soil ecosystems under abiotic stress.

Fig. 1
figure 1

Microcosm experimental setup. Soil samples were collected from Cr-contaminated hotspots (Cr(VI), 665 mg/kg) and from clean areas around the site (Cr(VI), 0.09 mg/kg). Cr(VI) was added to the clean soil at concentrations of 0, 10, and 100 mg/kg (added as K2Cr2O4) to create microcosm soil systems with different pollution levels. These soils were then sterilized three times using high-pressure sterilization. Using tangential flow filtration (TFF, membrane pore size 30 kDa), the microbiome (including prokaryotes and viruses) was enriched from both the clean and heavily Cr-polluted soils, abbreviated as CM and PM, respectively. These suspensions were then inoculated into sterilized soils of varying Cr(VI) concentrations to construct microcosm systems, named C-0, C-10, C-100, P-0, P-10, and P-100, which were incubated at room temperature for 35 days. Following incubation, DNA was extracted, and both long-read (MinION Nanopore) and short-read (Illumina Novaseq 6000) sequencing were performed

Full size image

Materials and methods

Soil collection and physicochemical analysis

Soils used in this study were collected from a Cr-polluted site (abbreviated as P) and relatively clean areas (abbreviated as C) in Xining, Qinghai Province, China. Soil sampling was conducted by using a five-point sampling method at a depth of 5-15 cm. The soil properties (e.g., pH, organic matter, total phosphorus, total sulfur, total potassium, available phosphorus, total nitrogen, and nitrate nitrogen), the concentrations of total Cr, Cr(VI), and other heavy metals were described in Table S1. Soil pH was measured using a pH meter in a 1:2.5 soil-to-water suspension. Organic matter and total nitrogen were analyzed by using the dichromate oxidation method and the Kjeldahl method, respectively. Total phosphorus, total potassium, and available phosphorus were determined by molybdenum blue spectrophotometry and flame photometry. Cr(VI) concentration was measured by using the diphenylcarbazide method, while total Cr and other heavy metals were analyzed by inductively coupled plasma mass spectrometry (ICP-MS) [12].

Microcosm setup

To prepare the soil suspension, 5L of buffer solution (68 mmol/L NaCl, 10 mmol/L MgSO₄, pH adjusted to 7.5 using Tris-Cl) was added to 5 kg of either polluted or clean soil. The mixture was evenly distributed into 500 mL bottles, manually stirred, and shaken for 3 min. Subsequently, the mixture was incubated at 4 °C for 30 min. Then, the suspension was filtered through a membrane with a pore size of 5 μm, followed by filtration using a tangential flow filtration (TFF) system (membrane pore size, 30 kDa). During the enrichment process, elution buffer was added multiple times to remove any possible residual Cr until the detection concentration was below 0.01 mg/kg. The final volume of the retentates was adjusted to 200 mL. Finally, the retentates, representing microorganisms derived from polluted soil (PM) and clean soil (CM), were stored at 4 °C and inoculated into different microcosms the following day.

Three sets of 3 kg of clean soil were taken and potassium dichromate (K2Cr2O4) was added to each to prepare soils with total Cr(VI) concentrations of 0 mg/kg, 10 mg/kg, and 100 mg/kg, representing low, moderate, and heavy contamination levels, respectively. Subsequently, the contaminated soils were subjected to high-pressure sterilization (121 °C, 20 min) and left for 1 day before the next sterilization cycle, repeated three times in total. Subsequently, soil sterility was assessed using plate culture techniques (LB medium). Once validated, the sterile soils were utilized for subsequent cultivation experiments.

The retentate PM was divided into 10 portions, each comprising 20 mL. Nine of these portions were separately added to 200 g of sterilized soil with low, moderate, or heavy contamination (three replicates). This setup aimed to investigate the response of microbial communities subjected to contamination on stress relief. The treatments were designated as P-0, P-10, and P-100. Similarly, the retentate CM obtained earlier was divided into 10 portions. Nine of these portions were separately added to 200 g of sterilized soil with different levels of contamination (named C-0, C-10, and C-100, three replicates) to investigate the response of microbial communities without contamination stress as disruption increased. The inoculated culture dishes were covered with breathable sealing film and placed in a sterile incubator. All treatments were incubated under sterile conditions at room temperature (25 °C) for 35 days, with regular addition of sterile water to maintain soil moisture content stable (Fig. 1).

DNA extraction, metagenomic sequencing, and assembly

After microcosm cultivation, soil DNA extraction was performed using the DNeasy PowerSoil (Qiagen, Hilden, Germany). Additionally, DNA was extracted directly from the collected polluted and clean soil samples. Due to the low microbial biomass in the polluted and oligotrophic soil samples, 10 g of soil per sample was used. Multiple extractions were performed for each sample following the kit instructions, and the resulting DNA was pooled and concentrated using freeze-drying. DNA quantification was carried out using the Qubit 3.0 fluorometer with the dsDNA system (Invitrogen, Waltham, MA). For Illumina sequencing, libraries were prepared using the TruSeq DNA PCR-Free Kit (Illumina). For ONT sequencing, libraries were prepared using the Ligation Sequencing Kit (SQK-LSK109, Oxford Nanopore Technologies), following the manufacturer’s protocols. Subsequently, short-read sequencing (Illumina Novaseq 6000, PE150, 10G) and long-read sequencing (MinION nanopore, 3G) were conducted for all samples [17]. Raw sequencing reads produced by MinKNOW were subject to base-calling utilizing Albacore (v2.1.10) to generate fasta files. The resultant reads, deemed of sufficient quality, underwent adapter trimming via PoreChop (0.2.3, https://github.com/rrwick/Porechop), wherein the parameter “–discard_middle” was employed to eliminate reads containing internal adapters. MinION sequencing data statistical assessment was conducted and visualized using NanoPack. Fastp (v0.23.2) was employed for quality control of reads (-D -5 -W 5 -M 20 -q 15 -u 40 -l 50 –thread 12). Subsequently, the hybrid assembly of short and long reads for individual samples was performed using OPERA-MS (v0.9.0) (–no-polishing –no-ref-clustering –no-gap-filling) [18]. For short-reads, metaSPAdes v3.13.0 [19] (-t 200 -m 600) and megahit v1.2.9 [20] (–k-min 21 –k-max 141 –k-step 10 –min-contig-len 500) were utilized for single-sample assembly. The assembly quality was assessed with QUAST v5.2.0. Due to the insufficient quality of DNA extracted from sample C-0-1, sequencing was not performed.

Functional analysis of prokaryotic communities

Contigs were deduplicated by clustering using MMseqs2 with the parameters –min-seq-id 0.95 and -c 0.90. Prodigal v2.6.1 was employed to predict open reading frames (ORFs) of the deduplicated contigs. Then, gene clustering was performed using Mmseqs with default parameters (-e 0.001 –min-seq-id 0.95 –c 0.90) at the nucleotide level. The longest sequence in each cluster was selected as the representative sequence, resulting in a non-redundant gene catalog (Unigenes). Clean reads were aligned to each ORF using BWA to quantify gene abundance (-k 19, transcripts per million, TPM). Functional annotation of the CDSs was performed using EggNOG-mapper with default parameters against the KEGG and eggNOG databases [21]. Heavy metal resistance genes were annotated using DIAMOND (blastp) with the BacMet EXP database (e-value < 1e − 5) [22].

Construction and classification of MAGs

The assembly contigs of a single sample were clustered to recover bins using MaxBin2 v2.2.7, CONCOCT v1.0.0, and MetaBAT2 v2.0 with default parameters as integrated in the MetaWRAP pipeline [23]. After that, the “bin_refinement” and “bin_reassembly” modules in MetaWRAP v1.2.1 were used to deduplicate and improve the quality of bins, and high-quality bins (considered metagenome-assembled genomes, MAGs) with completeness > 70% and contamination < 10% were retained [23]. The taxonomic annotation was implemented based on the Genome Taxonomy database Release 207 (GTDB) using GTDB-Tk v1.5.0 with the “classify_wf” module, and the phylogenetic trees of all bacterial MAGs and MAGs with CRISPR-Cas system were constructed using GTDB-Tk v1.5.0 with the “infer” module [24]. A total of 746 medium-quality MAGs were obtained from the hybrid assembly (completeness > 70%, contamination < 10%), with an average genome length of 3570 kb. These MAGs represented 1 archaeal and 13 bacterial phyla. CoverM was used to estimate the relative abundance of MAGs across samples (-m TPM, –trim-min 0.10, –trim-max 0.90, –min-read-percent-identity 0.95, –min-read-aligned-percent 0.80).

Identification and analyses of viral sequences

The viral sequences were identified using Virsorter2 based on viral genomic content and structural features (–min-length 3000 -j 4 all) [25]. Simultaneously, DeepVirFinder was also employed for viral identification, with a selection criterion of score > 0.9 and p-value < 0.01 [26]. Subsequently, TBtools extracted relevant viral sequences based on the combined results from Virsorter2 and DeepVirFinder. CheckV was then applied to recognize and eliminate host contamination and assess the confidence and completeness of identified viral sequences [27]. Viral sequences longer than 3 kb were retained for further analysis. Then, the viral lifestyle was determined by PhaTYP [28]. Prodigal was employed for predicting ORFs of viral sequences, followed by manual annotation of their functions according to the BacMet EXP database using DIAMOND (blastp) [29]. The TPM of viral contigs was calculated using CoverM (-m TPM, –trim-min 0.10, –trim-max 0.90, –min-read-percent-identity 0.95, –min-read-aligned-percent 0.80).

Identification of plasmid, provirus, and antiviral defense systems

GeNomad was used to identify and extract provirus and plasmid sequences with the default parameters (using the end-to-end command), due to its superior performance in recognizing mobile genetic elements (MGEs), particularly plasmids and prophages [30]. DefenseFinder was employed to identify the antiviral defense systems and defense genes with the default parameters [6]. Input files comprised assemblies of prokaryotic, plasmid, and viral contigs in fasta format. To ensure accurate defense system identification, only systems in which all required genes were located on the same contig were considered present. Additionally, CRISPRCasFinder was utilized for the identification of CRISPR arrays and Cas proteins [31]. CRISPR spacer sequences were extracted and counted using CRISPRCasTyper [32]. Notably, to avoid the loss of defense information due to redundant operations, the analysis for the antiviral defense system was conducted on all 749 MAGs that were not redundant.

Host prediction of viruses

To ensure accurate host prediction, viral contigs incorrectly binned into prokaryotic MAGs were manually removed after binning. This step was necessary because viral sequences misassigned to MAGs could lead to false associations between the virus and its host. Specifically, viral contigs were identified and removed based on the previously determined IDs of free viruses. To establish links between viral contigs and prokaryotic MAGs, multiple methods were used, including CRISPR spacer matching and tRNA matching. CRISPR spacer matching is based on the principle that CRISPR-Cas systems capture fragments of foreign DNA, such as viral sequences, and store them as spacers in the host genome, enabling targeted recognition and interference. Similarly, tRNA matching relies on the observation that viruses often acquire tRNA genes from their hosts during co-evolution, resulting in high sequence similarity between viral and host tRNA genes. Specifically, (i) CRISPR spacer matching: CRISPR spacers of MAGs were identified via CRISPRCasTyper (–prodigal meta) to detect CRISPR spacers on repeat sequences, followed by a blastdb of prokaryotic CRISPR spacers was created and BLASTn with an e-value of 1e − 10 was utilized for virus-host relationship prediction (100% identity over 100% coverage). (ii) tRNA matching: Identification of tRNA genes in viral contigs and prokaryotic MAGs was performed using tRNAscan-SE (v2.0.9) [32] with parameters “-A” and “-B”. Subsequently, a blastdb of prokaryotic tRNA genes was created and BLASTn with an e-value of 1e − 10 was utilized for matching viral and prokaryotic tRNA genes (100% identity over 100% coverage). Following this, the defense systems carried by the associated MAGs and viral contigs were compared. Furthermore, for viral contigs carrying CRISPR-Cas systems, a CRISPR spacer was extracted and matched with MAGs sequences and other viral contigs from the same sample to establish targeting connections.

Transcription analysis of microbial antiviral defense systems and heavy metal resistance gene under varying Cr(VI) concentrations

Following the construction of the microcosm systems as previously described, we conducted a second round of identical microcosm cultivation. After microcosm cultivation, soil samples were collected for RNA extraction to investigate the expression of defense systems and heavy metal resistance genes. In brief, total RNA was extracted from 30 mg soil using RNeasy Mini Kit (Qiagen, Germany), which quality and quantity were determined by Agilent 4200 system (Agilent Technologies, Waldbron, Germany) and NanoDrop2000 (Thermo Fisher, MA, USA), respectively. The rRNA transcripts were eliminated by Ribo-Zero rRNA Removal Kits. The mRNA was converted into cDNA, followed by purification and amplification to obtain the sequencing library using the TruSeq™ RNA Sample Prep Kit. Paired-end sequencing was conducted on a NovaSeq6000 system to generate 150 bp paired-end reads. Due to samples used for RNA and DNA extraction not being from the same batch, a non-referential transcriptome analysis method was adopted, as detailed below. Fastp (v0.23.2) was used for quality control of sequencing reads as metagenomics. Non-protein-coding rRNA sequences from eukaryotes and prokaryotes were filtered by Sortmerna. Afterward, the remaining clean reads were assembled into a de novo transcriptome using Trinity (v2.15.1) with the parameters: –pre_correction –mink 20 –maxk 60 –step 10. Prodigal (v2.6.3) was used to predict ORFs for assembled contigs longer than 200 bp and the unigenes were clustered using methods consistent with metagenomics. BWA was utilized to align the clean reads back to the unigenes for quantification of gene expression levels (-k 19, TPM). The identification of defense systems and genes, and the functional annotation of genes are consistent with the previous metagenomic analysis.

Meta-analyses of public datasets from various heavy metal-contaminated sites

Metagenomic sequence data were collected from the NCBI database to corroborate the findings of the microcosm system. These data encompassed several bioprojects, namely PRJNA694468, focusing on soils from Lan Mu Chang (LMC), a site in Guizhou, China with co-contamination of multiple heavy metals, and relatively clean southern soils (MS) [33]; PRJNA883072 and PRJNA886109, sourced from soil samples within the world’s largest Sb mining area in XKS, Hunan, China. The sampled areas included pollution epicenters at the surface (HH), distant low-contamination points (HL), as well as vertical sampling across high-contamination and low-contamination zones (VH and VL) [34]. Additionally, samples from Cr-polluted sites at varying concentrations in Luzhou, Sichuan (LZ), and Zhangye, Gansu (ZY), were incorporated (GSA database PRJCA003236 and PRJCA004276) [12]. The heavy metal content of these soil samples was recorded in Table S3. These sequences underwent quality control using fastp (v0.23.2), and were subsequently assembled using metaSPAdes (v3.13.0). Identification, and lifestyle predictions of viral sequences, as well as the recognition of antiviral defense systems, were conducted following previously described procedures.

Statistical analyses

All statistical analyses were conducted by R v4.0.0 and Python 3.9. The Wilcoxon rank sum test was used for pairwise comparison between different treatments, including all data presented in this article using box plots, such as the relative abundance of contigs and MAGs carrying antiviral defense systems, as well as the proportion of temperate viruses. Correlation analysis between antiviral defense systems and microbial lysogenic levels, as well as heavy metal resistance genes while employing Spearman correlation coefficients. Graphical analyses were generated using the online platform https://www.chiplot.online [35].

Results

The overall prokaryotic antiviral defense potential strengthened with increased heavy metal-induced stress

The hybrid assembly of short-read and long-read sequencing significantly improved the assembly of longer sequences, whereas Illumina short reads generated a larger number of assembly contigs compared to the hybrid assembly. Specifically, the numbers of contigs longer than 5 kb, 10 kb, 25 kb, and 50 kb from the hybrid assembly increased by 16%, 32%, 40%, and 46%, respectively, compared to metaSPAdes assembly, and by 8%, 27%, 38%, and 45%, respectively, compared to MEGAHIT assembly (Table S2). Additionally, the mean N50 values for contigs obtained from metaSPAdes, MEGAHIT, and hybrid assembly were 2410, 2288, and 5230, respectively. These results corroborate that hybrid assembly significantly enhances the quality of assembled long contigs, therefore all subsequent analyses were based on the hybrid assembly.

We obtained 745 medium-quality MAGs with completeness exceeding 70% and contamination below 10%, with only one MAG being archaea. We identified 3119 antiviral defense systems encoded by 6634 genes in 628 MAGs (accounting for 84.2%). The frequency of restriction-modification (RM) systems based on methylation-modification recognition was the highest (1004 systems in 53.5% MAGs), followed by the independent protein SoFIC with Fic domains (256 systems in 28.6% MAGs), adaptive immunity CRISPR-Cas (228 systems in 24.7% MAGs) and abortion infection AbiE (177 systems in 21.0% MAGs) (Table S4). There were 237 MAGs carrying more than 10 types of antiviral defense systems, suggesting their widespread coexistence. The distribution of most antiviral defense systems on the genome exhibited characteristics of defense islands, with the RM system most commonly appearing adjacent to other systems (Fig. 2A). Among the antiviral defense systems with the top 20 average relative abundance carried by MAGs, AbiE exhibited the broadest distribution range spanning 11 bacterial phyla, followed by RM and SoFIC carried by 10 bacterial phyla (Fig. 2B).

Fig. 2
figure 2

The distribution and abundance response of defense systems to Cr-induced stress. A Defense islands on the genomes. B The distribution patterns of defense systems within prokaryotes. C The number of defense system genes per unit length of MAGs (kb). D The number of defense systems per unit length of MAGs (kb). E Relative abundance of contigs with antiviral defense systems. (C, clean soil from the site (Cr(VI), 0.09 mg/kg); P, Cr-polluted soil from the site (Cr(VI), 665 mg/kg); C-0, C-10, and C-100 respectively refer to microcosm systems with the addition of Cr(VI) concentrations of 0, 10, and 100 mg/kg, and inoculation with microbiome from clean soil; P-0, P-10, and P-100 respectively refer to microcosm systems with the addition of Cr(VI) concentrations of 0, 10, and 100 mg/kg, and inoculation with microbiome from polluted soil. Boxplot components, center lines, medians; box limits, 25th and 75th percentiles; whiskers, 1.5 × interquartile range from the 25th and 75th percentiles; statistical significance symbols: * p ≤ 0.05, ns p > 0.05 (Wilcoxon rank sum test)

Full size image

The relative abundance of MAGs with antiviral defense systems showed no significant variations across different microcosm systems (Fig. S1). However, upon normalization for the number of antiviral defense systems and associated genes, MAGs from highly polluted soils harbored a higher number of antiviral defense systems and genes per unit length of the genome compared to those from slightly and moderately polluted soils (p < 0.05, Fig. 2C, D). Additionally, for the overall microbial community, the relative abundance of contigs carrying antiviral systems significantly increased from 0.33% in C-0 to 0.43% in C-100, and that from microbes in polluted sites significantly decreased after stress reduction (p < 0.05, Fig. 2E). This indicates that the defense potential of prokaryotes strengthens simultaneously with increasing pollution stress level.

Furthermore, microbial data from four sites (Hunan, Gansu, Sichuan, and Guizhou, China) contaminated by Sb, Cr, and multiple heavy metals, were collected to validate the changes in prokaryotic antiviral defense systems facing heavy metal concentrations. A total of 4386 antiviral defense systems were identified, encompassing 103 types and 175 subtypes. Consistent with observations in our microcosm experiments, the occurrence frequency of RM systems (1453) far exceeded that of any other antiviral defense system, followed by CRISPR-Cas (352), AbiE (261), SoFIC (245), and pAgo (218), with all other systems below 200 (Table S5). Additionally, prokaryotes facing higher heavy metal concentrations harbored more antiviral defense systems in these soils (Fig. S2), corroborating that the prokaryotic antiviral defense systems are enhanced under pollution stress.

The enhancement of prokaryotic antiviral defense systems associated with improved prokaryote-temperate virus mutualism

The lysogens level of microbial communities is an important indicator of prokaryote-virus mutualism. The proportion of temperate viruses in contaminated soil was higher than that in clean soil from fields (37% vs. 32%, Fig. 3A). Furthermore, within microcosm systems, the proportion of temperate viruses in C-0, C-10, and C-100 was 34%, 35%, and 42%, respectively. As Cr(VI) concentration decreased, the proportion of temperate viruses significantly decreased (p < 0.05) from 38% (P-100) to 34% (P-10) and 32% (P-0). Moreover, the GeNomad-based results revealed a significant decrease in the number of provirus carried by unit length contigs with decreasing Cr(VI) pollution (p < 0.05, Fig. 3B). Therefore, our findings corroborated that prokaryote-virus interactions exhibit a sensitive response to pollution stress, with mutualism strengthening as pollution intensified and mutualism diminishing as stress decreased [12].

Fig. 3
figure 3

The association of prokaryotic antiviral defense systems with prokaryote-virus symbiosis and heavy metal resistance. A The proportion of temperate viruses in the viral community (Wilcoxon rank sum test). B Number of proviruses carried by contig of unit length (Mb) (Wilcoxon rank sum test). C Correlation analysis of prokaryotic defense systems with microbial lysogeny (Spearman correlation). D The heatmap of heavy metal resistance and growth-related functional genes of prokaryotes. E Correlation analysis of prokaryotic defense systems with prokaryotic heavy metal resistance genes, as well as genes associated with cell growth (Spearman correlation). Statistical significance symbols: *p ≤ 0.05, ns p > 0.05

Full size image

The analysis of the relative abundances of contigs carrying the top 20 antiviral defense systems reflected the defense strategies of prokaryotes. There were significant positive correlations between the proportion of temperate viruses and the relative abundance of contigs harboring defense systems AbiE, RloC, Reston, dGTPase, and Shedu (p < 0.05, Fig. 3C). Moreover, the frequency of proviruses was positively correlated with the relative abundance of contigs harboring defense systems dGTPase, Lamassu-Fam, defense-associated reverse transcriptase (DRT) and dCTPdeaminase. Conversely, the relative abundances of contigs carrying systems CRISPR-Cas and MazEF exhibited significant negative correlations with the frequency of proviruses (p < 0.05, Fig. 3C).

For microbiomes from the other four heavy metal contaminated sites, virus-prokaryotic mutualism was also reinforced with increased pollution stress as reflected by the elevated proportion of temperate virus (Fig. S3), which was consistent with the discoveries in microcosm systems (Fig. 3B). Moreover, microbiomes were further analyzed to gain insights into prokaryotic antiviral defense systems from Cr-contaminated soils in Luzhou (LZ), Sichuan and Zhangye (ZY), Gansu. In LZ, a significant positive correlation was observed between the proportion of temperate virus and the relative abundance of RosmeTA, SspBCDE, and RloC systems (Fig. S4, p < 0.01). Similarly, in ZY, a significant relative positive correlation was found between the proportion of temperate viruses and the abundance of RosmeTA, RloC, and Retron systems (Fig. S5, p < 0.01). These antiviral defense systems primarily function by inhibiting viral transcription and translation, as well as inducing cell abortion during the late stages of viral infection [3]. Overall, the relative abundance of prokaryotic antiviral defense systems that promote virus-prokaryote mutualism increased with intensified heavy metal stress.

The synchronization of changes in prokaryotic virus resistance and heavy metal resistance facing elevated abiotic stress

Exploring the trade-off between biotic and abiotic resistance in prokaryotic communities can advance the understanding of their environmental adaptation. Heavy metal resistance genes primarily encompass oxidation-reduction, efflux pumps, biofilm formation, and nucleic acid repair [16]. The relative abundance of these genes in C-100 was significantly higher than that in C-0, and there was a slight decline in the relative abundance of heavy metal resistance genes in P-0 compared to P-100 (Fig. 3D). Based on transcriptomic results, microbiomes transplanted from clean soil notably upregulated expression of Cr resistance gene in the C-100 treatment compared to C-0, and those transplanted from polluted soil significantly downregulated such gene expression in the P-0 treatment compared to P-100 (p < 0.01, Fig. S6). These results corroborate that microbial communities exhibit fitness costs for developing pollution stress resistance under Cr-induced stress.

Interestingly, the relative abundance of the top 20 antiviral defense systems was mostly positively correlated with the relative abundance of heavy metal resistance genes (p < 0.05, Fig. 3E), while negatively correlated with that of genes related to cell growth (e.g., RloC, Reston and dGTPase). It suggests prokaryotic overall viral resistance and heavy metal resistance exhibited synchronized changes in environments with Cr-induced stress. However, it should be specifically mentioned that the relative abundance of contigs carrying systems such as CRISPR-Cas and MazEF were significantly negatively correlated with genes associated with heavy metal resistance (p < 0.05, Fig. 3E).

Prokaryotes carrying adaptive immune CRISPR-Cas system exhibited reduced survival ability in highly polluted soils

The CRISPR-Cas system, as the unique autoimmune system, has great potential for application in multiple fields, and can sensitively reflect the dynamics of prokaryote-virus interactions [36]. We identified CRISPR-Cas systems in 185 MAGs (accounting for 24.8% of all the MAGs). Among these, 51 MAGs contained two or more CRISPR-Cas subtype systems (Table S6). The most frequently occurring subtypes were I-C, I-G, I-E, and I-F (Fig. 4A and Fig. S7). Moreover, MAGs harboring the CRISPR-Cas system were classified into 10 phyla, mainly Pseudomonadota, Bacteroidota, and Verrucomicrobiota (Fig. 4A, B).

Fig. 4
figure 4

Distribution patterns and environmental responses of CRISPR-Cas systems. A Phylogenetic tree of MAGs harboring CRISPR-Cas systems, along with the features of CRISPR-Cas systems. B The relative abundance of contigs with CRISPR-Cas systems (Wilcoxon rank sum test, statistical significance symbols: *p ≤ 0.05, ns p > 0.05). C The distribution patterns of MAGs harboring CRISPR-Cas systems. D Composition of CRISPR-Cas subtypes carried by MAGs

Full size image

For the sites, MAGs with relative abundance exceeding 30% carried the CRISPR-Cas system in clean soil, while MAGs hosting the CRISPR-Cas system remained undetectable in Cr-contaminated soils (Fig. S8). Moreover, a decrease in the relative abundance of MAGs carrying the CRISPR-Cas system was observed during both pollution increase and alleviation treatments as the Cr(VI) concentration increased in microcosm systems (Fig. S9). Similar trends were observed at the overall community level, where the relative abundance of contigs carrying the CRISPR-Cas system was higher in low-pollution treatments (Fig. 4C). These results suggest a substantial fitness cost associated with the maintenance of spacer sequences and the complex action of Cas proteins [37]. Furthermore, the CRISPR-Cas system was distributed across a broader range of phyla in C-0 and C-10 microbiomes, while it was mainly distributed in Pseudomonadota within the C-100, P-100, P-10, and P-0L microbial communities (Fig. 4B). Regarding subtypes, MAGs carrying I-E, I-F, and II-C displayed relatively higher abundance in stress increase treatments, whereas those carrying I-C and I-D exhibited higher abundance in stress alleviation treatments (Fig. 4D). Therefore, the subtype composition and distribution of prokaryotic CRISPR-Cas systems alternated as pollution changed.

During CRISPR system activities, spacer sequences are transcribed and processed into CRISPR RNA (crRNA) for viral immunity recognition [38]. Therefore, the number of spacers can reflect the defense spectra of CRISPR systems against viruses [37]. There were 162 MAGs with CRISPR-Cas systems carrying 5085 spacers (Table S7 and Fig. 4A). The number of spacer arrays was not correlated to the genome size of MAGs, but may be influenced by the taxonomic classification of the prokaryote and the number of associated CRISPR systems (Fig. 4A). The count of spacers per MAG ranged from 2 to 246, with 7 MAGs harboring over 100 spacers (Fig. 4A and Fig. S10). These MAGs belonged to Saccharospirillum in Pseudomonadota, NIC37A-2 in Myxococcota, Imperializer in Bacteroidota, and UBA6623 in Acidobacteriota, suggesting their heightened potential for virus immunity. The matching linkages between spacers and viral genomes can reflect prokaryote-virus interactions [39]. Among the 140 linkages of prokaryote-virus established based on spacer match, 41 were associated with temperate viruses, whereas 99 were associated with virulent viruses (Fig. S11).

Plasmid-carried antiviral defense systems increased relative abundance and diversity in polluted soils

Plasmids, as pivotal vessels of horizontal gene transfer, can carry and facilitate the transfer of antiviral defense systems among microbial communities [40]. A total of 12,442 antiviral defense systems were detected in these microbial communities, with 20.4% distributed on plasmids (Table S8). Specifically, approximately 4% of plasmids (2,088) carried 2,539 antiviral defense systems, encompassing 88 types and 138 subtypes, with predominant systems including RM, RloC, Gabija, AbiE, and CBASS (Fig. 5A). The proportion of antiviral defense systems carried by plasmid among microbial community increased by about 43% from C-0 to C-100, and decreased from 26% in P-100 to 22% in P-0 (p < 0.05, Fig. 5B).

Fig. 5
figure 5

Defense systems carried by plasmids. A Overview of defense systems carried by plasmids in overall microbial communities. B The proportion of defense systems carried by plasmid (%). C Relative abundance of plasmids within microbial communities. D Genomic structure of some plasmid genomes simultaneously carrying defense systems and heavy metal resistance genes. Statistical significance symbols: ***p ≤ 0.001, **p ≤ 0.01, *p ≤ 0.05, ns p > 0.05 (Wilcoxon rank sum test)

Full size image

Furthermore, there was a significant increase in the relative abundance of plasmids from low contamination to heavy contamination (p < 0.01, Fig. 5C). Additionally, nearly 20% of plasmids carrying antiviral defense systems also harbored genes associated with heavy metal resistance, exemplified by plasmid pMOL28 (Fig. 5D). Consequently, the diverse and numerous antiviral defense systems carried by plasmids constituted a vital component of the microbial defense systems, and the augmentation of prokaryotic antiviral defense systems in heavily contaminated soils was closely related to the enrichment of plasmids.

Temperate viruses contributed more to defending against virulent viruses in stressful environments

Viruses, as the main targets of prokaryotic antiviral defense systems, can also carry antiviral defense systems [41, 42]. There were 1051 prokaryotic antiviral defense systems, consisting of 63 types. The proportion of antiviral systems carried by viruses among the overall defense systems of each sample ranged from 5.3 to 9.2%, with a higher proportion in highly polluted soil (Fig. 6A). These systems are distributed among 3.5% of viral contigs (822), with 162 contigs harboring more than two types (Table S8). Additionally, 108 viral contigs possessed both antiviral defense systems and heavy metal resistance genes. Notably, these antiviral defense systems were predominantly distributed in temperate viral contigs (56–71%) (Fig. 6B). The abundance of viruses carrying antiviral defense systems remained relatively stable at around 3.5% under stress alleviation treatments (P-100, P-10, and P-0). However, the abundance of viruses carrying antiviral defense systems increased from approximately 1% in C-10 to 3.5% in C-100 (p < 0.05, Fig. 6C). Moreover, the changes in abundance of viruses carrying antiviral defense systems exhibited strong consistency with that of viral auxiliary metabolic genes (AMGs) related to Cr resistance (Fig. 6C and Fig. S12). The RM system was the most common antiviral defense system carried by viruses (Fig. 6D). In treatments P-100, P-10, P-0, and C-100, a greater diversity of antiviral defense systems was observed among the viral populations compared to C-10 and C-0 (Fig. 6D). Therefore, the antiviral defense systems carried by viruses also constitute an important component of the defense, contributing to the evolution and propagation of prokaryotic antiviral defense systems.

Fig. 6
figure 6

Defense systems carried by viruses. A The proportion of antiviral defense systems carried by viruses among all microbial antiviral defense systems. B The lifestyle of viruses carrying antiviral systems. C The relative abundance of viral contigs carrying defense systems. D The type composition of defense systems carried by viruses. Statistical significance symbols: *p ≤ 0.05, ns p > 0.05 (Wilcoxon rank sum test)

Full size image

Host prediction established 93 links between viruses carrying antiviral defense systems and prokaryotic MAGs. Surprisingly, only 24 viruses harbored antiviral defense systems identical to those of their hosts, while 69 viruses carried antiviral defense systems absent in their associated MAGs (Table S9). Furthermore, spacers of 27 virus-carried CRISPR-Cas systems were matched against other virus and host genomes, and the results showed that 9 viruses from treatments P-100, P-10, and P-0 directly target other viruses within the same treatment, establishing a total of 18 links (Table S10). However, no spacer sequences targeting prokaryotic MAGs were found. This suggests that under stress conditions, viruses could not only facilitate the horizontal transfer of antiviral defense systems but also directly participate in the defense of prokaryotes against other viruses.

Discussion

Although prokaryotic antiviral defense systems are ecologically significant, our understanding of prokaryotic defense systems in complex environments remains limited. Therefore, it is crucial to elucidate whether and how environmental abiotic stresses drive the evolution of prokaryotic antiviral defense systems and understand their ecological roles in shaping microbial adaptability. Our study addresses this gap by investigating the distribution patterns, dynamic response, and ecological impacts of the prokaryotic defense systems in soil ecosystems under heavy metal-induced stress. We used a combination of hybrid metagenomic assembly, microcosm systems, and site investigations to provide comprehensive insights into microbial defense strategies within soil ecosystems.

Our metagenomic and metatranscriptomic analyses revealed increased density, diversity, and activity of antiviral defense systems harbored by prokaryotes in heavily polluted soils (Fig. 2E and Fig. S13), which suggests that Cr(VI)-induced stress could actually foster the enrichment of antiviral defense systems. This aligns with the findings of Wu et al. [43], who demonstrated the plasticity of bacterial defensomes in response to environmental pressures. Note that changes in the defense systems are closely related to microbial composition shifts driven by pollution stress. As shown in Fig. 2B, the distribution of some defense systems exhibits species preferences, indicating that community composition changes influence the overall defense phenotype. Correspondingly, the regulation of phage-host interactions by specific defense systems also contributes to the dominance of particular taxa in polluted environments. Moreover, there was a significant positive correlation between the relative abundance of prokaryotes carrying these antiviral defense systems (such as dGTPase, RloC, Retron, and dCTPdeaminase) and the prevalence of temperate viruses and provirus within microbial communities (Fig. 3C). The defensive mechanisms of dGTPase, RloC, and Retron involve the inhibition of viral replication and translation but do not prevent the injection of viral DNA into the host cell, indicating that these systems may promote the cessation of temperate viruses in their lysogenic state (proviruses) [44, 45]. Similarly, DRT strongly suppressed the expression of viral late genes, such as capsid proteins, whereas early and middle genes were not substantially affected [46]. Thus, increasing the presence of these systems might promote and maintain the interaction between prokaryotes and temperate viruses under stressful environments, advancing our understanding of the ecological mechanisms behind the formation of prokaryote-virus mutualism under stress conditions.

The adaptive immune CRISPR-Cas systems, are widely distributed in bacterial and archaeal genomes [47], and Wang et al. predicted that the presence of CRISPR-Cas systems would decrease as the demand for beneficial genes increases within cells [48]. Our results confirmed this prediction at the community level, whereby prokaryotes carrying CRISPR-Cas systems decreased in highly polluted soils (Fig. 4C and Fig. S9), which favored cells acquiring heavy metal resistance genes carried by plasmids or proviruses. However, in low-pollution soils where viruses become the primary mortality pressure, prokaryotes carrying these antiviral defense systems were reselected (Fig. 4C). Two additional reasons for this result may be that CRISPR immunity systems (except for type III targeting RNA) were likely to be incompatible with the integration of proviruses, as the CRISPR response would create an auto-immune problem [49]. Additionally, the complex mechanism of Cas protein’s action and the acquisition and possession of spacers sequence [37], require greater fitness cost and make it difficult to achieve immediate virus defense compared to other antiviral defense systems [50]. Therefore, in highly polluted soils, there was a decrease in prokaryotes carrying CRISPR-Cas systems which favored establishing mutualism between the overall prokaryotic community and proviruses, thereby reinforcing adaptive phenotypes to Cr(VI)-induced stress.

The relationship between biotic resistance and prokaryotic abiotic resistance is a crucial aspect of understanding their adaptive mechanisms in adverse environments. For example, the evolution of resistance to environmental pressures, such as antibiotics and other abiotic stresses, can lead to pleiotropic interactions that may constrain virus resistance [51]. Our previous studies have shown decreased soil microbial biomass due to enhanced Cr(VI) contamination [12]. Interestingly, contrary to a simple trade-off where cell growth slowed as heavy metal resistance increased, the majority of antiviral defense systems exhibit a concurrent increase with the heavy metal resistance of microbial communities under stressful conditions (Fig. 3E). Furthermore, transcriptome analysis also indicated that with increased Cr(VI) concentration, the expression of antiviral defense systems was upregulated concurrently with Cr resistance genes (Fig. S6 and Fig. S13, Table S11). However, defense systems with high fitness, such as CRISPR-Cas, exhibited an opposite trend to heavy metal resistance (Fig. 3E). Therefore, the reinforcement of antiviral defense systems conducive to prokaryote-virus mutualism and abiotic stress resistance may be attributed to cellular sacrifices in proliferation and the loss of some systems with high adaptive costs. Additionally, a substantial proportion of plasmids and temperate viruses simultaneously carry antiviral defense systems and heavy metal resistance genes (Fig. 5D), suggesting synchronicity in transferring these genes mediated by MGEs. This also indicated that the metabolic burden of antiviral defense systems might been compensated by other beneficial genes. Among them, as a representative MGE with both antiviral defense systems and heavy metal resistance genes, plasmid pMOL28 has significant potential in future research on heavy metal remediation ecological engineering and microbial genetic evolution. Moreover, previous studies have unveiled the non-defensive roles of prokaryotic antiviral defense systems, such as microbial quorum sensing and biofilm formation [38]. The augmentation of these antiviral defense systems may be attributed to their direct involvement not only in modulating the interaction between prokaryotes and viruses but also in the adaptive regulation of prokaryotic communities to adverse conditions.

Although the crucial roles of viruses and plasmids in the horizontal transfer of antiviral defense systems are well recognized [52], the diversity of antiviral defense systems within plasmid and viral genomes in environmental microbial communities and their responses to abiotic stress remain poorly elucidated. Our results revealed the diversity and distribution patterns of antiviral defense systems carried by plasmids and viruses in Cr-contaminated soil. We quantified that plasmids and viruses harbor comparable antiviral defense systems (4% vs 3.5%), providing valuable resources for exploring the phenomenon of plasmid- and virus-mediated viral defense. As environmental stress escalates, prokaryotes could enhance the uptake of mobile elements harboring beneficial genes, as evidenced by the rise in the abundance of provirus and plasmid (Figs. 3B and 5C). As expected, compared to microbiomes in clean or slightly polluted soils, those in highly contaminated soils showed a higher prevalence of antiviral defense systems carried by MGEs, which facilitated horizontal gene transfer and resource sharing of defense systems, which are important factors in enhancing the defense of the microbial community.

Viruses are the primary targets of antiviral defense systems, and increasing evidence suggests that the provirus can directly assist the host in blocking virulent virus growth, promoting bacterial survival, and enabling efficient lysogeny [53,54,55]. For instance, Dedrick et al. demonstrated that (p)ppGpp synthetase coded by prophage Phrann can defend against phage Tweety infection [41]. Moreover, Staphylococcus aureus prophages can block the reproduction of the induced Tha-positive prophages by coding the tail-activated defense gene Tha [54]. Our research provides support for establishing compelling conclusions at the community level, including the finding that viruses carrying antiviral defense systems are mainly temperate viruses (Fig. 6B), and the antiviral defense systems they carry exhibit significant distinctions from the corresponding host’s antiviral defense systems. Additionally, the spacers of CRISPR-Cas systems carried by viruses target other virulent viruses from the same environment rather than their host. It was elucidated that the frequency of viruses carrying antiviral defense systems in soil habitats increased with escalating stress levels (Fig. 6), emphasizing the significant impact of virus-carried antiviral defense systems on the adaptive evolution of prokaryotes in adverse environments.

To summarize, we depicted the adaptive modification of microbial defense systems in soil ecosystems under heavy metal disturbances, revealing how these communities enhance their adaptability by adjusting prokaryotic immune strategies. This includes illustrating the impact of antiviral defense systems on prokaryote-virus mutualism, recognizing the synchronous changes between microbial antiviral defense and heavy metal resistance, and emphasizing the critical role of MGEs in the evolution of prokaryotic antiviral defense systems in stressful environments. Our research advances the understanding of prokaryote-virus interactions and microbial adaptation strategies under stress conditions. However, further investigations are needed to explore the mechanisms of collaboration among prokaryotic antiviral defense systems, and the potential non-defensive functions of certain antiviral defense systems, as they may directly contribute to the evolution of microbial resistance to abiotic stress.

Data availability

Datasets are publicly available at the Sequence Read Archive – NCBI (PRJNA1104006).

References

  1. Chevallereau A, Pons BJ, van Houte S, Westra ER. Interactions between bacterial and phage communities in natural environments. Nat Rev Microbiol. 2022;20(1):49–62.

    Article  CAS  PubMed  Google Scholar 

  2. Georjon H, Bernheim A. The highly diverse antiphage defence systems of bacteria. Nat Rev Microbiol. 2023;21(10):686–700.

    Article  CAS  PubMed  Google Scholar 

  3. Doron S, Melamed S, Ofir G, Leavitt A, Lopatina A, Keren M, Amitai G, Sorek R. Systematic discovery of antiphage defense systems in the microbial pangenome. Science. 2018;359(6379):eaar4120.

    Article  PubMed  PubMed Central  Google Scholar 

  4. Hampton HG, Watson BNJ, Fineran PC. The arms race between bacteria and their phage foes. Nature. 2020;577(7790):327–36.

    Article  CAS  PubMed  Google Scholar 

  5. Millman A, Melamed S, Leavitt A, Doron S, Bernheim A, Hor J, Garb J, Bechon N, Brandis A, Lopatina A, et al. An expanded arsenal of immune systems that protect bacteria from phages. Cell Host Microbe. 2022;30(9):1556–69.

    Article  CAS  PubMed  Google Scholar 

  6. Tesson F, Hervé A, Mordret E, Touchon M, d’Humières C, Cury J, Bernheim A. Systematic and quantitative view of the antiviral arsenal of prokaryotes. Nat Commun. 2022;13(1):2561.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  7. Zhen X, Xu X, Ye L, Xie S, Huang Z, Yang S, Wang Y, Li J, Long F, Ouyang S. Structural basis of antiphage immunity generated by a prokaryotic argonaute-associated SPARSA system. Nat Commun. 2024;15(1):450.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  8. Beavogui A, Lacroix A, Wiart N, Poulain J, Delmont TO, Paoli L, Wincker P, Oliveira PH. The defensome of complex bacterial communities. Nat Commun. 2024;15(1):2146.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  9. Huang W, Traulsen A, Werner B, Hiltunen T, Becks L. Dynamical trade-offs arise from antagonistic coevolution and decrease intraspecific diversity. Nat Commun. 2017;8:2059.

    Article  PubMed  PubMed Central  Google Scholar 

  10. Burmeister AR, Turner PE. Trading-off and trading-up in the world of bacteria-phage evolution. Curr Biol. 2020;30(19):R1120–4.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  11. Jessup CM, Bohannan BJM. The shape of an ecological trade-off varies with environment. Ecol Lett. 2008;11(9):947–59.

    Article  PubMed  Google Scholar 

  12. Huang D, Yu PF, Ye M, Schwarz C, Jiang X, Alvarez PJJ. Enhanced mutualistic symbiosis between soil phages and bacteria with elevated chromium-induced environmental stress. Microbiome. 2021;9(1):1.

    Article  Google Scholar 

  13. Zheng X, Jahn MT, Sun M, Friman V-P, Balcazar JL, Wang J, Shi Y, Gong X, Hu F, Zhu YG. Organochlorine contamination enriches virus-encoded metabolism and pesticide degradation associated auxiliary genes in soil microbiomes. ISME J. 2022;16(5):1397–408.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  14. Tang X, Zhong LR, Tang L, Fan CZ, Zhang BW, Wang MR, Dong HR, Zhou CY, Rensing C, Zhou SG, et al. Lysogenic bacteriophages encoding arsenic resistance determinants promote bacterial community adaptation to arsenic toxicity. ISME J. 2023;17(7):1104–15.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  15. Feiner R, Argov T, Rabinovich L, Sigal N, Borovok I, Herskovits AA. A new perspective on lysogeny: prophages as active regulatory switches of bacteria. Nat Rev Microbiol. 2015;13(10):641–50.

    Article  CAS  PubMed  Google Scholar 

  16. Xia X, Wu SJ, Zhou ZJ, Wang GJ. Microbial Cd(II) and Cr(VI) resistance mechanisms and application in bioremediation. J Hazard Mater. 2021;401:123685.

    Article  CAS  PubMed  Google Scholar 

  17. Che Y, Xia Y, Liu L, Li AD, Yang Y, Zhang T. Mobile antibiotic resistome in wastewater treatment plants revealed by Nanopore metagenomic sequencing. Microbiome. 2019;7(1):44.

    Article  PubMed  PubMed Central  Google Scholar 

  18. Bertrand D, Shaw J, Kalathiyappan M, Ng AHQ, Kumar MS, Li CH, Dvornicic M, Soldo JP, Koh JY, Tong CX, et al. Hybrid metagenomic assembly enables high-resolution analysis of resistance determinants and mobile elements in human microbiomes. Nat Biotechnol. 2019;37(8):937–44.

    Article  CAS  PubMed  Google Scholar 

  19. Nurk S, Meleshko D, Korobeynikov A, Pevzner PA. metaSPAdes: a new versatile metagenomic assembler. Genome Res. 2017;27(5):824–34.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  20. Li D, Liu CM, Luo R, Sadakane K, Lam TW. MEGAHIT: an ultra-fast single-node solution for large and complex metagenomics assembly via succinct de Bruijn graph. Bioinformatics. 2015;31(10):1674–6.

    Article  CAS  PubMed  Google Scholar 

  21. Cantalapiedra CP, Hernández-Plaza A, Letunic I, Bork P, Huerta-Cepas J. eggNOG-mapper v2: functional annotation, orthology assignments, and domain prediction at the metagenomic scale. Mol Biol Evol. 2021;38(12):5825–9.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  22. Pal C, Bengtsson-Palme J, Rensing C, Kristiansson E, Larsson DGJ. BacMet: antibacterial biocide and metal resistance genes database. Nucleic Acids Res. 2013;42(D1):D737–43.

    Article  PubMed  PubMed Central  Google Scholar 

  23. Uritskiy GV, DiRuggiero J, Taylor J. MetaWRAP-a flexible pipeline for genome-resolved metagenomic data analysis. Microbiome. 2018;6:158.

    Article  PubMed  PubMed Central  Google Scholar 

  24. Chaumeil PA, Mussig AJ, Hugenholtz P, Parks DH. GTDB-Tk: a toolkit to classify genomes with the Genome Taxonomy Database. Bioinformatics. 2020;36(6):1925–7.

    Article  CAS  Google Scholar 

  25. Guo JR, Bolduc B, Zayed AA, Varsani A, Dominguez-Huerta G, Delmont TO, Pratama AA, Gazitúa MC, Vik D, Sullivan MB, et al. VirSorter2: a multi-classifier, expert-guided approach to detect diverse DNA and RNA viruses. Microbiome. 2021;9(1):37.

    Article  PubMed  PubMed Central  Google Scholar 

  26. Ren J, Song K, Deng C, Ahlgren NA, Fuhrman JA, Li Y, Xie XH, Poplin R, Sun FZ. Identifying viruses from metagenomic data using deep learning. Quant Biol. 2020;8(1):64–77.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  27. Nayfach S, Camargo AP, Schulz F, Eloe-Fadrosh E, Roux S, Kyrpides NC. CheckV assesses the quality and completeness of metagenome-assembled viral genomes. Nat Biotechnol. 2021;39(5):578–85.

    Article  CAS  PubMed  Google Scholar 

  28. Shang JY, Tang XB, Sun YN. PhaTYP: predicting the lifestyle for bacteriophages using BERT. Brief Bioinform. 2023;24(1):bbac487.

    Article  PubMed  Google Scholar 

  29. Pal C, Bengtsson-Palme J, Rensing C, Kristiansson E, Larsson DGJ. BacMet: antibacterial biocide and metal resistance genes database. Nucleic Acids Res. 2014;42(D1):D737–43.

    Article  CAS  PubMed  Google Scholar 

  30. Camargo AP, Roux S, Schulz F, Babinski M, Xu Y, Hu B, Chain PSG, Nayfach S, Kyrpides NC. Identification of mobile genetic elements with geNomad. Nat Biotechnol. 2023;42:1303–12.

    Article  PubMed  PubMed Central  Google Scholar 

  31. Couvin D, Bernheim A, Toffano-Nioche C, Touchon M, Michalik J, Néron B, Rocha EPC, Vergnaud G, Gautheret D, Pourcel C. CRISPRCasFinder, an update of CRISRFinder, includes a portable version, enhanced performance and integrates search for Cas proteins. Nucleic Acids Res. 2018;46(W1):W246–51.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  32. Russel J, Pinilla-Redondo R, Mayo-Muñoz D, Shah SA, Sorensen SJ. CRISPRCasTyper: automated identification, annotation, and classification of CRISPR-Cas loci. CRISPR J. 2020;3(6):462–9.

    Article  CAS  PubMed  Google Scholar 

  33. Xu R, Sun X, Haggblom MM, Dong Y, Zhang M, Yang Z, Xiao E, Xiao T, Gao P, Li B, et al. Metabolic potentials of members of the class Acidobacteriia in metal-contaminated soils revealed by metagenomic analysis. Environ Microbiol. 2022;24(2):803–18.

    Article  CAS  PubMed  Google Scholar 

  34. Liu S, Zeng J, Yu H, Wang C, Yang Y, Wang J, He Z, Yan Q. Antimony efflux underpins phosphorus cycling and resistance of phosphate-solubilizing bacteria in mining soils. ISME J. 2023;17(8):1278–89.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  35. Xie J, Chen Y, Cai G, Cai R, Hu Z, Wang H. Tree Visualization By One Table (tvBOT): a web application for visualizing, modifying and annotating phylogenetic trees. Nucleic Acids Res. 2023;51(W1):W587–92.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  36. Wu R, Chai B, Cole JR, Gunturu SK, Guo X, Tian R, Gu JD, Zhou J, Tiedje JM. Targeted assemblies of cas1 suggest CRISPR-Cas’s response to soil warming. ISME J. 2020;14(7):1651–62.

    Article  PubMed  PubMed Central  Google Scholar 

  37. McGinn J, Marraffini LA. Molecular mechanisms of CRISPR-Cas spacer acquisition. Nat Rev Microbiol. 2019;17(1):7–12.

    Article  CAS  PubMed  Google Scholar 

  38. Devi V, Harjai K, Chhibber S. CRISPR-Cas systems: role in cellular processes beyond adaptive immunity. Folia Microbiol. 2022;67(6):837–50.

    Article  CAS  Google Scholar 

  39. Dion MB, Plante P-L, Zufferey E, Shah SA, Corbeil J, Moineau S. Streamlining CRISPR spacer-based bacterial host predictions to decipher the viral dark matter. Nucleic Acids Res. 2021;49(6):3127–38.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  40. Chevallereau A, Westra ER. Bacterial immunity: mobile genetic elements are hotspots for defence systems. Curr Biol. 2022;32(17):R923–6.

    Article  CAS  PubMed  Google Scholar 

  41. Dedrick RM, Jacobs-Sera D, Bustamante CAG, Garlena RA, Mavrich TN, Pope WH, Reyes JCC, Russell DA, Adair T, Alvey R, et al. Prophage-mediated defence against viral attack and viral counter-defence. Nat Microbiol. 2017;2:16251. https://doiorg.publicaciones.saludcastillayleon.es/10.1038/nmicrobiol.2016.251.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  42. Huang D, Yuan MM, Chen J, Zheng X, Wong D, Alvarez PJJ, Yu P. The association of prokaryotic antiviral systems and symbiotic phage communities in drinking water microbiomes. ISME Commun. 2023;3(1):46.

    Article  PubMed  PubMed Central  Google Scholar 

  43. Wu Y, Garushyants SK, van den Hurk A, Aparicio-Maldonado C, Kushwaha SK, King CM, Ou Y, Todeschini TC, Clokie MRJ, Millard AD, et al. Bacterial defense systems exhibit synergistic anti-phage activity. Cell Host Microbe. 2024;32(4):557-572.e556.

    Article  CAS  PubMed  Google Scholar 

  44. Tal N, Millman A, Stokar-Avihail A, Fedorenko T, Leavitt A, Melamed S, Yirmiya E, Avraham C, Brandis A, Mehlman T, et al. Bacteria deplete deoxynucleotides to defend against bacteriophage infection. Nat Microbiol. 2022;7(8):1200–9.

    Article  CAS  PubMed  Google Scholar 

  45. Davidov E, Kaufmann G. RloC: a wobble nucleotide-excising and zinc-responsive bacterial tRNase. Mol Microbiol. 2008;69(6):1560–74.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  46. Gao L, Altae-Tran H, Böhning F, Makarova KS, Segel M, Schmid-Burgk JL, Koob J, Wolf YI, Koonin EV, Zhang F. Diverse enzymatic activities mediate antiviral immunity in prokaryotes. Science. 2020;369(6507):1077–84.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  47. Mohanraju P, Makarova KS, Zetsche B, Zhang F, Koonin EV, van der Oost J. Diverse evolutionary roots and mechanistic variations of the CRISPR-Cas systems. Science. 2016;353(6299):aad5147.

    Article  PubMed  Google Scholar 

  48. Jiang W, Maniv I, Arain F, Wang Y, Levin BR, Marraffini LA. Dealing with the evolutionary downside of CRISPR immunity: bacteria and beneficial plasmids. Plos Gene. 2013;9(9):e1003844.

    Article  CAS  Google Scholar 

  49. Nobrega FL, Walinga H, Dutilh BE, Brouns SJJ. Prophages are associated with extensive CRISPR-Cas auto-immunity. Nucleic Acids Res. 2020;48(21):12074–84.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  50. Westra ER, Levin BR. It is unclear how important CRISPR-Cas systems are for protecting natural populations of bacteria against infections by mobile genetic elements. Proc Natl Acad Sci. 2020;117(45):27777–85.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  51. Burmeister AR, Fortier A, Roush C, Lessing AJ, Bender RG, Barahman R, Grant R, Chan BK, Turner PE. Pleiotropy complicates a trade-off between phage resistance and antibiotic resistance. Proc Natl Acad Sci. 2020;117(21):11207–16.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  52. Vassallo CN, Doering CR, Littlehale ML, Teodoro GIC, Laub MT. A functional selection reveals previously undetected anti-phage defence systems in the E. coli pangenome. Nat Microbiol. 2022;7(10):1568–79.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  53. Bondy-Denomy J, Qian J, Westra ER, Buckling A, Guttman DS, Davidson AR, Maxwell KL. Prophages mediate defense against phage infection through diverse mechanisms. ISME J. 2016;10(12):2854–66.

    Article  PubMed  PubMed Central  Google Scholar 

  54. Rostøl JT, Quiles-Puchalt N, Iturbe-Sanz P, Lasa Í, Penadés JR. Bacteriophages avoid autoimmunity from cognate immune systems as an intrinsic part of their life cycles. Nat Microbiol. 2024;9(5):1312–24.

    Article  PubMed  PubMed Central  Google Scholar 

  55. Fillol-Salom A, Rostøl JT, Ojiogu AD, Chen J, Douce G, Humphrey S, Penadés JR. Bacteriophages benefit from mobilizing pathogenicity islands encoding immune systems against competitors. Cell. 2022;185(17):3248-3262.e3220.

    Article  CAS  PubMed  Google Scholar 

Download references

Funding

This work was financially supported by the National Key Research and Development Program of China (2022YFC3704700), the National Natural Science Foundation of China (42277418, 42177113, and 42477131), and the Youth Innovation Promotion Association, Chinese Academy of Sciences (Y2022084).

Author information

Authors and Affiliations

  1. College of Environmental and Resource Sciences, Zhejiang University, Hangzhou, 310058, China

    Dan Huang, Dongsheng Wang & Pingfeng Yu

  2. Department of Civil and Environmental Engineering, Virginia Tech, Blacksburg, VA, 24060, USA

    Jingqiu Liao

  3. Catalan Institute for Water Research (ICRA-CERCA), Girona, 17003, Spain

    Jose Luis Balcazar

  4. Key Laboratory of Soil Environment and Pollution Remediation, Institute of Soil Science, Chinese Academy of Sciences, Nanjing, 210008, China

    Mao Ye

  5. Earth and Biological Sciences Directorate, Pacific Northwest National Lab, Richland, WA, 99352, USA

    Ruonan Wu

  6. Department of Civil and Environmental Engineering, Rice University, Houston, TX, 77005, USA

    Pedro J. J. Alvarez

Authors
  1. Dan Huang
    View author publications

    You can also search for this author inPubMed Google Scholar

  2. Jingqiu Liao
    View author publications

    You can also search for this author inPubMed Google Scholar

  3. Jose Luis Balcazar
    View author publications

    You can also search for this author inPubMed Google Scholar

  4. Mao Ye
    View author publications

    You can also search for this author inPubMed Google Scholar

  5. Ruonan Wu
    View author publications

    You can also search for this author inPubMed Google Scholar

  6. Dongsheng Wang
    View author publications

    You can also search for this author inPubMed Google Scholar

  7. Pedro J. J. Alvarez
    View author publications

    You can also search for this author inPubMed Google Scholar

  8. Pingfeng Yu
    View author publications

    You can also search for this author inPubMed Google Scholar

Contributions

All authors contributed intellectually to and agreed to this submission. HD, YPF and YM designed the study. HD collected and analyzed the data. HD and YPF wrote the initial draft of the manuscript while LQJ, YM, JLB, WRN, WDS and PA provided substantial critical feedback and editing.

Corresponding authors

Correspondence to Mao Ye or Pingfeng Yu.

Ethics declarations

Ethics approval and consent to participate

Not applicable.

Consent for publication

Not applicable.

Competing interests

The authors declare no competing interests.

Additional information

Publisher’s Note

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

Supplementary Information

40168_2025_2030_MOESM1_ESM.xlsx

Supplementary Material 1: Table S1. Physiochemical properties of field soils. Table S2. Metagenomic assembly quality metrics (QUAST evaluation). Table S3. Metadata collection information. Table S4. Distribution of defense systems in MAGs. Table S5. Information on defense systems in meta-analysis samples. Table S6. Distribution of CRISPR-Cas systems in MAGs. Table S7. Spacer information for CRISPR-Cas systems in MAGs. Table S8. Distribution of defense systems across viruses, plasmids, and prokaryotic genomes. Table S9. Host prediction for virus sequences carrying defense systems. Table S10. Virus sequences targeted by phage-carried CRISPR-Cas system spacers. Table S11. Transcriptional activity of defense systems. Table S12. Metagenome accessions and associated metadata. Table S13. Distribution of heavy metal resistance genes.

40168_2025_2030_MOESM2_ESM.docx

Supplementary Material 2: Figure S1. The relative abundance of MAGs carrying defense systems. Figure S2. Abundance of contigs carrying defense systems in four heavy metal-polluted soils. Figure S3. Proportion of temperate phages in microbial communities in four heavy metal-polluted soils. Figure S4. Correlation analysis between the top 20 defense systems and the proportion of temperate phages in Cr-contaminated soil from LZ. Figure S5. Correlation analysis between the top 20 defense systems and the proportion of temperate phages in Cr-contaminated soil from ZY. Figure S6. Transcriptional expression of Cr resistance genes. Figure S7. Composition of CRISPR-Cas subtype systems. Figure S8. The coverage of CRISPR-Cas system in MAGs. Figure S9. The relative abundance of contigs carring CRISPR-Cas system. Figure S10. CRISPR spacer sequences carried by MAGs. Figure S11. The life history of viruses associated with the CRISPR-Cas system. Figure S12. Relative abundance of viruses carrying auxiliary metabolic genes (AMGs). Figure S13. Transcription of contigs with defense systems.

Rights and permissions

Open Access This article is licensed under a Creative Commons Attribution-NonCommercial-NoDerivatives 4.0 International License, which permits any non-commercial use, sharing, distribution and reproduction in any medium or format, as long as you give appropriate credit to the original author(s) and the source, provide a link to the Creative Commons licence, and indicate if you modified the licensed material. You do not have permission under this licence to share adapted material derived from this article or parts of it. The images or other third party material in this article are included in the article’s Creative Commons licence, unless indicated otherwise in a credit line to the material. If material is not included in the article’s Creative Commons licence and your intended use is not permitted by statutory regulation or exceeds the permitted use, you will need to obtain permission directly from the copyright holder. To view a copy of this licence, visit http://creativecommons.org/licenses/by-nc-nd/4.0/.

Reprints and permissions

About this article

Check for updates. Verify currency and authenticity via CrossMark

Cite this article

Huang, D., Liao, J., Balcazar, J.L. et al. Adaptive modification of antiviral defense systems in microbial community under Cr-induced stress. Microbiome 13, 34 (2025). https://doiorg.publicaciones.saludcastillayleon.es/10.1186/s40168-025-02030-z

Download citation

  • Received:

  • Accepted:

  • Published:

  • DOI: https://doiorg.publicaciones.saludcastillayleon.es/10.1186/s40168-025-02030-z

Keywords

Microbiome

ISSN: 2049-2618

Contact us