This is the accepted version of the article: Chen, Qing-Lin; An, Xin-Li; Zheng, Bang-Xiao; [et al.]. «Loss of soil microbial diversity exacerbates spread of antibiotic resistance». Soil ecology letters, Vol. 1, issue 1-2 (June 2019), p. 3-13. DOI 10.1007/s42832-019-0011-0 This version is avaible at https://ddd.uab.cat/record/218194 under the terms of the license Loss of Microbial Diversity exacerbates the Spread of Antibiotic Resistome in Soil 1 Abstract: Biodiversity loss is considered as a major threat because of its importance 2 for ecosystem processes and services on which society depends. Natural ecosystems 3 show variable resistance to invasion by alien species, and this resistance can be related 4 to the microbial diversity in the system. For centuries, little has been known regarding 5 the role of soil microbial diversity in the biological barrier against microbial invasion. 6 The increasing prevalence of antibiotic resistant bacteria is one of the most serious 7 threats to public health in the 21st century. Here we explore whether reductions in soil 8 microbial diversity have consequences on the dissemination of antibiotic resistance. By 9 using a dilution-to-extinction approach coupled with a high-capacity quantitative PCR 10 arrays, we investigated the relationship between soil microbial diversity and the barrier 11 effect against antibiotic resistance invasion. We observed a negative correlation 12 between microbial diversity and the abundance of antibiotic resistance genes (ARGs), 13 and this correlation was maintained even when accounting simultaneously for other 14 drivers (incubation time, microbial abundance). Our results demonstrated that higher 15 microbial diversity can act as a biological barrier against antibiotic resistance spread. 16 Together with previous work, our results fill a critical gap in understanding the role of 17 soil microbial diversity in ecosystem health 18 Introduction 19 It is well documented that biodiversity is the foundation of the maintenance of 20 ecosystems, i.e. a large species pool is required to sustain a health ecosystem since they 21 play a key role in multiple ecosystem functions and services simultaneously, including 22 earth’s biogeochemical cycling, primary production, litter decomposition as well as 23 climate regulation 1-4. The rapid and continued development of molecular biology and 5-8 24 genomic techniques has unveiled immense microbial diversity in soil and ocean . 25 However, the roles of microbial diversity and associated traits in controlling ecosystem 26 functioning remain unclear 9-11. 27 In recent years, a few studies have shown that anthropogenic activities (such as 28 agricultural intensification, land use and nitrogen enrichment), and climate change may 29 reduce microbial diversity, a response that will likely impact ecosystem functions 12-15. 30 There is an extensive debate about relationships between microbial diversity and 16, 17 31 ecosystem functioning . Currently there are three different views about the 32 relationships. It has been proposed that soil microorganisms are key components that 33 determine life supporting functions, but functional redundancy in soil microorganisms 34 are prevalent, therefore an initial loss in microbial diversity is unlikely to substantially 35 affect ecosystem functions16, 18, 19. Studies have also demonstrated that the relationships 36 between microbial diversity and ecosystem functions are saturated 20, suggesting losses 37 of a few species at high richness levels could have minimal consequences on the 21 38 ecosystem functioning . For instance, microbial decomposer communities often 39 exhibit high redundancy for a single function, such as microbial respiration and biomass 20, 22 40 production . However, if the function of concern was more specific than general 41 functions, such as pesticide degradation are known to be limited to only a few 42 specialized functional groups23, the relationship between biodiversity and ecosystem 43 functioning is more linear than saturating 24. Moreover multifunctional redundancy was 44 generally lower, i.e. the degree of multiple functional dependence on diversity was 25 45 higher than single-functional redundancy . Therefore, the buffering capacity of 46 ecosystem functions against biodiversity loss may be limited that a moderate loss of 47 diversity may substantially impair key specialized functions 23. Additionally, Jung et al., 48 observed that genes related to the nitrogen cycling was significantly reduced, 49 contradictory results have also been obtained that the efficiency of diesel 26 50 biodegradation was increased in the low-diversity community . These results 51 indicated that the relationship between microbial diversity and ecological function 52 involves trade-offs among ecological processes, and should not be generalized as a 53 positive, neutral, or negative relationship 26. Each point of view sound reasonable, and 54 this could be, at least partly, the reason why the consequences of decline in biodiversity 55 for ecosystem processes and functioning have long been of considerable interest27-30. 56 Resistance to invasion by alien species represents a major life support function of 57 terrestrial ecosystems31, 32. The biological invasion, developed in the 1950s, which 58 indicated that ecosystems that contain a higher level of biodiversity are less vulnerable 59 to disturbances 33-35. Whether this theory can be applied to microbiology level remains 60 largely unresolved. To answer this question, some pioneering studies have investigated 61 the relationship between invasibility and the microbial diversity of ecosystems, by 19, 30, 36-38 62 manipulating the microbial community . For example, van Elsas et al., 63 observed that a positive correlation between the inoculant survival rate of Escherichia 64 coli O157:H7 and the soil fumigation depths; additionally, by using a dilution– 65 reinoculation approach they obtained a similar result that a negative correlation between 66 the soil microbial diversity and survival of the invader 36, 37, 39. An inverse relationship 67 was also found between the survival rate of the inoculant and the degree of microbial 68 diversity with respect to Pseudomonas aeruginosa, Ralstonia solanacearum and 69 Listeria monocytogenes39-41. These results shed lights on the complexity–invasiveness 70 relationship within microbial communities. However, with investigating limited species, 71 we still lack a sound evidences to achieve a universal phenomenon. To address this 72 question in an alternative way, we focus on the relationship between invasibility and 73 microbial diversity at gene-level. 74 Antibiotic resistance genes (ARGs) were chosen as the invaded genes. On one hand, 75 antibiotic-treatment failure is typically attributed to the “weapon-shield” role that ARGs 76 played in clinical settings 42. On the other, the pace of development of novel antibiotics 77 is now alarmingly low 43. Now the emergence and spread of antibiotic resistance has 78 become a global health threat 44-47. Nevertheless, factor in regulating the environmental 79 pathways of antibiotic resistance, especially the biotic factors, have not yet been 80 directly addressed. This hampers our ability to predict changes and risks in antibiotic 81 resistance under anthropogenic activities and global environmental change. 82 The main objective of this study was to investigate if, as proposed by the theory of 83 biological invasion, soil microbial species diversity could act as a biological barrier 84 preventing invasion by antibiotic resistance. In other words, we addressed the 85 consequences of soil diversity loss on the fate of antibiotic resistance genes in the soil 86 environment. Pig manure was used as the source of ARGs, and microbial diversity was 87 manipulated by employing a dilution–reinoculation approach (i.e. inoculating sterile 88 soil microcosms with serial dilutions of a soil microbial suspension) (Figure S1). To 89 test our hypothesis, we characterized bacterial communities and antibiotic resistome 90 using Illumina Miseq profiling of 16S ribosomal genes and high-capacity quantitative 91 PCR arrays with 296 primer sets targeting almost all major classes of ARGs 48, 49. Our 92 study provides empirical evidence that that microbial diversity negatively relates to 93 invasibility of antibiotic resistance in soil ecosystems; and further suggests that the loss 94 in microbial diversity will likely facilitate the spread of antibiotic resistance 95 Results and Discussion 96 Microbial community assembly and diversity 97 After the OTUs were classified according to the Ribosomal Database Project (RDP) 50 98 database , the soil microbial community assemble into 17 phyla (Figure 1A) and 99 dominated by Actinobacteria, Firmicutes, Proteobacteria, Bacteroidetes, 100 Planctomycetes and Chloroflexi, which are common bacterial phyla observed in soils 12, 51 101 worldwide . Firmicutes were significantly depleted along with the time, which is 102 contrasted with Actinobacteria and Bacteroidetes. Compared with the impact of time, 103 minor differences in the dominant bacterial phyla were observed when assessing 104 differences among the dilution depths. Nevertheless, Proteobacteria are significantly 105 (P < 0.01) enriched after high dilution. 106 Rarefaction curves were constructed for each individual sample showing the number of 107 observed OTUs (Figure 1B), defined at a 97% sequence similarity cut-off in QIIME 52, 53 108 . To assess the sequencing depth, we calculated Good’s coverage scores (Figure 1B) 109 and they were highly comparable for all samples ranging from 94.25 to 98.51% 110 indicating that the sequencing depth was adequate to reliably describe the bacterial 111 microbiome. As expected, rarefaction curve confirmed that the dilution treatment had 112 resulted in progressively decreasing species richness and diversity, which indicated that 113 dilution-to-extinction method applied here is one of the few available methods to 38 114 manipulate microbial biodiversity of complex natural ecosystems such as the soil . 115 Furthermore, the impacts of incubation time exhibited a higher degree of variation in 116 the shape of their rarefaction curves as compared to the dilution treatments. The 117 observed OTUs (richness) and diversity of OTU were significantly depleted with the 118 increasing incubation time (Figure S2). Similarly, principal coordinates (PCoA) 119 analysis of weighted UniFrac distances between samples revealed (Figure S3) that 120 microbial communities clustering together according to their initial dilution levels 121 along with PC2 which explained 10.41% of the total variations. Whereas microbial 122 communities of different incubation time are separated along with PC1 (explained 123 53.55% of the total variations). The bacterial abundance was also measured using Real- 124 time Quantitative PCR (qPCR) based on the 16S rRNA gene. The results showed that 125 in the first ten days, the bacterial abundance was decreased in treatment of S and D0 126 (which had a relatively higher diverse native microbiota), whereas, to a lower extent, a 127 slighted increase was observed in rest of treatments (which had a relatively lower 128 diverse native microbiota). After 10 days inoculation, the bacterial abundance was 129 stable over time in all constructed microcosms, being similar between treatments. These 130 findings are consistent with those from studies showing that some of manure-derived 131 bacteria could not thrive in soil environment, and gradually decreased after manure 132 treatment, which can be attributed to the competition with resident soil bacteria and the 133 differences in environmental conditions between soil and animal gut 54-56. Competition 134 for nutrients is an important mechanism that may limit invasions in highly diverse 57 135 communities . Based on the resource-based niche theories, the establishment of 136 invading species is dependent on the amount of (limiting) resources that are left 137 unconsumed by native species, as well as by the rate at which native and invader species 58 138 consume the existing resources . Additionally, most of the gut microbiota from 139 animals and humans are restricted to growth under anaerobic conditions 59, which are 140 highly different from the aerobic conditions in our microcosms. Although competition 141 for resources and shifting in oxygen condition likely affect the fate of manure-borne 142 microbiota in soil, they are only pieces of the puzzle. Other mechanisms i.e. predation 143 and negative species interactions, might also alert the survival of manure-borne species 37 144 . Figure 1 A: Phylum distribution of the OTUs (relative sequence abundance of bacterial phyla); B: Average Good’s coverage estimates (%) and rarefaction curves of different dilution treatments. Good’s coverage estimates represent averages of 20 replicates ± standard deviation. Antibiotic resistance genes and correlation with microbial diversity HT-qPCR was performed to investigate the abundance and diversity of ARGs. A total of 195 genes including 186 unique ARGs, class 1 integron-integrase gene and 8 transposase genes, were detected (Figure S?). The numbers of ARGs detected in each samples ranged from 47 to 116. The normalized abundance of ARG was ranging between 0.05 and 1.18 copies per 16S rRNA gene. As expected, the dilution treatment affected ARG abundance, with abundance of ARG increasing as dilution increased (Figure 2). At the begging of incubation (Day 0), because of treatments of S and D0 had a higher bacterial abundance (Figure S?), which were proposed to caused a lower normalized abundance of ARG. While the normalized abundance of ARGs between treatments became evident since day 10 and these differences were attribute to the effect of dilution, because of the bacterial abundance become stable and showed no significant different between treatments. Figure 2. Normalized abundance of ARG associated with the six dilution treatments during the 90 days incubation. Box plots display the first (25%) and third (75%) quartiles, the median and the maximum and minimum observed values within each data set. The biggest obstacle to understanding of the importance of microbial biodiversity for the functioning of ecosystems is the lack of sound experimental approaches to make 38 directed and predictable changes in the diversity of microbial communities in soil . Both the present and prior studies demonstrated that dilution-to-extinction method was an effective way to manipulate the soil microbial community to achieve a gradient in 30, 37 microbial diversity and species richness . By using ordinary least-squares (OLS) regression model, to our knowledge, we first explored the relationship between microbial diversity, estimated with the Inverse Simpson diversity indices and the spread of antibiotic resistance evaluated using the ARG abundance. Loss in microbial diversity (P = 0.0003, R2 = 0.1043) were linearly associated with increase in ARG abundance (Figure 3). We also found biodiversity components such as phylogenetical diversity (P = 0.0003, R2 = 0.0946) and species richness (P = 0.0016, R2 = 0.0810) were also highly and negatively related to the dissemination of antibiotic resistance. In contrast, nonsignificant relationship between evenness and ARG abundance was found (P = 0.8122, R2 = 0.0005). Further analyses provided evidence that Inverse Simpson diversity was positively and strongly related to phylogenetical diversity (P < 0.0001, R2 = 0.7227) and species richness (P < 0.0001, R2 = 0.4264). Additionally, phylogenetical diversity was also positively related to the species richness, whereas, evenness, was an exception, had no significant correlation with Inverse Simpson diversity or any of biodiversity components (Figure s?). Albeit our results are correlative in nature, and thus, the results reported here are cannot be taken as a definitive proof of causation, nevertheless, these findings are consistent with those from theoretical and experimental studies showing that a higher diverse soil microbial communities can act as a biological barrier against invasion 32, 39, 60, 61 . For example, a recent study showed that a significant negative correlation was detected between the survival rate of Listeria monocytogenes L9 and the Inverse Simpson metric (ρ = -0.817, P < 0.05) and suggested that erosion of microbial diversity may have damaging effects 39 regarding circulation of pathogenic microorganisms in the soil environment . Our results provide the first empirical evidence at gene-level to show that natural ecosystems could show variable resistance to invasion by aliens, and this resistance was relate to the species diversity in the system. Finally, our results support the hypothesis that microbial diversity is the foundation for the maintenance of ecosystems 1. Although, functional redundancy in soil microorganisms are considered prevalent, and has been thought to overwhelm any type of diversity–function relationship, i.e. initial loss in 19, 37 microbial diversity was unlikely to substantially affect ecosystem functions . However, if we take the biological barrier effect of microbial diversity into consideration, the functional redundancy could be overestimated in the previous studies. Figure 3 Ordinary least squares (OLS) regression model showing the relationships between ARG abundance, and microbial diversity, biodiversity components and evenness. The solid blue lines indicate statistical significance for the relationships, while the dashed lines indicate no statistical significance for the relationships. The shaded areas show the 95% confidence interval of the fit. Figure 4. Relationships between ARG abundance and incubation time. The solid lines represent the fitted OLS model and the shaded areas show the 95% confidence interval of the fit. Box plots display the first (25%) and third (75%) quartiles, the median and the maximum and minimum observed values within each data set (n=24). In addition to the impact of dilution treatment on ARG abundance, we also found that ARG abundance depleted with the increasing of incubation time (P<0.0001) (Figure 4). To further investigate the direct and indirect effects of dilution treatment and incubation time on profile of ARG, we generated structural equation models (SEMs) based on the known effects and relationships. Besides, SEM is an a priori approach offering the ability to separate multiple pathways of influence and take them as a system, and is useful to explore the complex networks of relationships found in ecosystems 62. Except for bacterial diversity and abundance (16S rRNA copies), MGEs was also included in our models, given their effects on profiles of ARG. Our model explained 86% of the variance found in the patterns of ARGs. Incubation time and dilution treatment can directly impact the patterns of ARGs or indirectly by strongly affecting the diversity and abundance of bacteria and abundance of MGEs. Bacterial diversity posed significant directly effect on ARG (ρ=0.09, P < 0.05), and indirectly impacted the patterns of ARGs by strongly affecting the MGEs abundance (ρ=-0.07, P < 0.05). While bacterial abundance on the other hand showed non-significant directly impacts on ARG (ρ=0.02, P > 0.05). Despite these results, the standardized effects from SEMs revealed that MGEs abundance had a direct positive effect on the pattern of ARG and was the major contributors to this model, which was align with field studies 63, indicating that MGEs were the most dominant factor altering the ARG profiles. Both the direct and indirect effects of dilution treatment showed a positive correlation with ARG abundance, and the indirect effect (through affecting the diversity and abundance of bacteria and abundance of MGEs) contributing larger than the direct effect. The incubation time had a negative effect on the pattern of ARG and the decrease in microbial diversity and MGEs abundance were the major driver. Compared with bacterial abundance, the diversity had a higher impacts in shaping the ARG profiles, although these impacts was much smaller than MGEs, dilution treatment and incubation time. Figure 5. Structural equation models showing the direct and indirect effects of dilution treatment, time, bacterial abundance, bacterial diversity, and MGEs on the ARG patterns. Black and red arrows indicate positive and negative relationships, respectively. Continuous and dashed arrows indicate significant and nonsignificant relationships, respectively. Numbers adjacent to arrows are path coefficients, and width of the arrows is proportional to the strength of path coefficients. R2 denotes the proportion of variance explained. Significance levels are indicated: *P < 0.05, **P < 0.01, and ***P < 0.001. Standardized effects (total ,direct, and ndirect effects) derived from the structural equation models. The hypothetical models fit our data well, as suggested by χ2 = 0.49, P = 0.48, df=1, GFI =0.99, and RMSEA = 0.00. Conclusions The provisioning of ecosystem services essential for human development heavily relies on the diversity of soil microorganisms. By using the microcosm model, our results provide empirical evidence that loss in microbial diversity will likely facilitate the proliferation and spread of antibiotic resistance. Our results fill a critical gap in our understanding the role of microbial diversity and provide additional insights that soil microbial communities can act as a biological barrier against invasion and the erosion of diversity may alert the circulation of antibiotic resistance genes in the environment. Altogether, ecosystems might have a limited buffering capacity of multiple ecosystem functions against biodiversity loss. Materials and Methods Microcosms and Experimental Setup Surface soil (0 to 20 cm) was collected from a cropland used for planting rice in Jiaxing, Zhejiang, south China (30°50′7.7″ N, 120°43′5.7″ E). Pig manure was obtained from a local commercial pig farm. Soil and pig manure properties are described in Table S1. The soil was sieved to < 2mm, part of the soil was sterilized by γ irradiation (35 kGy), and then 3600g irradiated soil were divided into 120 flasks (200-ml) as matrix of soil microcosms, whereas the rest of soil was used, in the dilution-to-extinction experiment, as an inoculum. To test the sterility 0.5 g of irradiated soil was spreading onto tryptic soy broth (TSB) and potatodextrose agar (PDA) media []. No bacterial and fungal growth on agar plates was observed after six days. An initial soil suspension was prepared by mixing 30g soil with 50 ml autoclaved demineralized water using a ultraviolet-sterilized blender at the maximum speed for 5 minutes[], equivalent 0.6g soil/ml approximately. Totally, four levels of dilution of the soil suspension were used as inocula to create a diversity gradient of soil microorganisms ranging from undiluted (D0) to diluted 10−2 (D1), 10−4 (D2) and 10−8 (D3) suspensions, and 5ml of suspension was subsequently inoculated into 200 ml flasks containing 30 g dry of sterile soil, equivalent to 10−1, 10−3, 10−5, 10−9 g of non-sterile soil/g sterile soil. Additionally, a positive control and negative control by using non-sterile soil (S) and sterile soil (SS) were also performed. For each treatment, four replicates were established. At last, 1.2 g pig manure as the source of antibiotic resistance was added to each flask. The flasks were then closed with sterile lids and the inoculated microcosms incubated at 20 °C, moisture content was maintained at 70% water-holding capacity (WHC), by addition of autoclaved demineralized water. At days 10, 20, 30, 60, and 90 four replicates of each combination of soil type were collected for analysis of total bacterial abundance, diversity and composition, and determination of antibiotic resistance. Assessing microbial diversity For each replicate microcosm from each dilution treatment, total DNA was extracted from 0.5g soil using a FastDNA®Spin Kit for soil according to the supplier's manual (MP Biomedical, Santa Ana, California, USA). The quality of the DNA was checked by using ND-1000 spectrophotometer (NanoDrop Technology, Wilmington, DE, USA). The concentration of DNA was determined using QubitTM dsDNA HS Assay kit through a fluorometer (QubitTM 3.0, USA). PCR was performed using 1μl of each forward (515F) and reverse (907R) bar-coded primers, 25μl of 2x ExTaq polymerase (TAKARA BIO INC, Japan) and 1μl of sample DNA as the template and 22 ul nuclease- free PCR-grade water in a total volume of 50 μl with the following PCR program of 95 °C for 5 min, followed by 30 cycles each of 95 s for 30 s, 58 °C 30 s, 72 °C for 30 s. To detect any contamination during PCR preparation, negative controls (template DNA was replaced with water) were included for all PCR reactions. PCR products of each subsample from the bar-coded primers were generated in four replicates and purified using the Wizard SV Gel and PCR Clean-Up System (TIANGEN Biotech, Beijing, China). The purified PCR products that were quantified and pooled at the same concentration, and then submitted to Illumina Hiseq2500 platform (Novogene, Beijing, China) for sequencing. To guarantee the quality of downstream analysis, raw pair-end reads were filtered to discard raw reads containing three or more ambiguous nucleotides, or with a low (< 20) average quality score, or with a short (< 100 nt) length, and barcode to generate clean joined reads capturing the complete V4-V5 region of the 16S rRNA gene by Novogene. The generated high quality sequences were processed and analyzed using QIIME pipeline []. The operational taxonomic unit (OTU) was identified using the UCLUST algorithm (Edgar, 2010) with a phylotype defined at the 97% sequence similarity level. Chimeric sequences, chloroplast and mitochondrial OTUs (around 1%), and singleton OTUs were discarded from the final OTU table. Taxonomic classification and quantification of OTUs were aligned against the Ribosomal Database Project database. Raw sequences were deposited in the National Center for Biotechnology Information (NCBI) Sequence Read Archive (SRA) under the access number SRP108158. Assessing antibiotic resistance genes A total of 296 primer sets were used to interrogate the soil DNA, these primers have been used, and validated in a previous study, which targeted resistance for all major classes of antibiotics (285 primer sets), transposase genes (8 primer sets), the universal class 1 integron-integrase gene (intI1), the clinical class 1 integron-integrase gene (cintI1) and the 16S rRNA gene (Supplementary Data 2). Before amplification with Wafergen SmartChip Real-time PCR system, all soil DNA was diluted to 50 ng μl−1 using sterile water. Amplification was conducted in a 100 nL reaction system and all qPCR reactions were conducted in triplicate and for each primer set, a non-template negative control was included. A more detailed description about the experimental procedure can be found in previous studies. HT-qPCR data was analyzed using SmartChip qPCR software (V 2.7.0.1). Reactions with poor melting curve analysis and reactions with amplification efficiency beyond the range (90%-110%) were discarded. Then screened with conditions that (1) a threshold cycle (CT) must be < 29 and (2) positive samples should have three replicates simultaneously. Relative copy number was calculated according to methods described in ref. 43: relative gene copy number =10(29−CT )/(10/3), where CT refers to quantitative PCR results and 29 refers to the detection limit. To differentiate the variations in the bacterial abundance, absolute 16S rRNA copy numbers were quantified by the standard curve (SC) method with Roche 480 system (Roche Inc., USA). Each 20 μl qPCR mixture consisted of 10 μl 2× LightCycle 480 SYBR Green I Master (Roche Applied Sciences), 0.5 μg μl–1 bovine serum albumin (BSA), 1 μM each primer, 1 ng μl–1 DNA as template and 6 μl nuclease-free PCR-grade water. The thermal cycle was set according to the previous descriptions (). A plasmid control containing a cloned and sequenced 16S rRNA gene fragment (1.18 × 1010 copies per microlitre) was used to generate eight-point calibration curves from ten-fold dilutions for standard calculation. All qPCRs were performed in technical triplicates with non-template negative controls. Statistical analysis By using ordinary least squares (OLS) regression models we explored the relationships between pattern of antibiotic resistance and microbial diversity estimated with the Inverse Simpson diversity, evenness and biodiversity components including phylogenetical diversity and species richness. OLS regression was conducted with R and visualized with package of “ggplot 2”. We used Structural equation model (SEM) (ref. 32) to evaluate the direct and indirect relationships between dilution treatments, incubation time, bacterial diversity, bacterial abundance, MGEs abundance, and pattern of antibiotic resistance. SEM is an a priori approach offering the capacity to visualize the casual relationships between variables by fitting data to the models representing causal hypotheses. Thus, the first step in SEM requires establishing an a priori model based on the known effects and relationships among the drivers shifting antibiotic resistance. The theoretical model assumptions were as follows (Figure S?): (i) Dilution treatment and incubation time might have direct influences on MGEs and the ARG patterns; (ii) and also they could indirectly influence the patterns of MGEs and ARGs by changing the bacterial abundance and diversity; (iii) bacterial abundance and diversity might have direct influences on ARG patterns and they could indirectly influence the patterns of ARGs by changing the MGEs abundance and detected number. Before modeling, we examined the distributions of all of our variables and tested their normality. All data was standardized and using Z-score in the downstream analysis. Bivariate correlation was performed to examine the pairwise correlations among these variables using SPSS 20 (IBM, Armonk, NY, USA), and the covariance matrix was imported into AMOS 21 (SPSS Inc., Chicago, IL, USA) for SEMs construction using the maximum-likelihood estimation (MLE) method. We parameterized our model using our dataset and tested its overall goodness of fit. There is no single universally accepted test of overall goodness of fit for SEMs, so we used multiple goodness-of-fit criteria. We used the chi-square (χ2) test (the model has a good fit when χ2 is low and P - value is high, traditionally P > 0.05), the root mean square error of approximation (RMSEA; the model has a good fit when RMSEA is near 0 and the probability is high, traditionally P > 0.05), high goodness-of-fit index (> 0.90) and low Akaike information criteria (AIC). With a reasonable model fit we were free to interpret the path coefficients of the model and their associated P-values. A path coefficient is analogous to a partial correlation coefficient and describes the strength and sign of the relationships between two variables. The solid and dot line represent a significant (P < 0.05) and nonsignificant correlation (P > 0.05), respectively. And we used black and red line to display the positive and negative correlation. Additionally, SEM is capable to partition direct and indirect effects that one variable may have on another and estimate the strengths of these multiple effects. To aid final interpretation in light of this capability of SEM, we calculated the standardized direct, indirect and total effects of dilution treatment, incubation time, bacterial abundance, bacterial diversity, and abundance of MGEs. The net effect of one factor upon another was calculated by summing all direct and indirect pathways between the two factors. If the model fits the data well, the total effect should approximate the bivariate correlation coefficient for that pair of factors 1. Wagg, C.; Bender, S. F.; Widmer, F.; van der Heijden, M. G., Soil biodiversity and soil community composition determine ecosystem multifunctionality. Proc Natl Acad Sci U S A 2014, 111, (14), 5266-70. 2. Wardle, D. A.; Bardgett, R. D.; Klironomos, J. N.; Setala, H.; van der Putten, W. H.; Wall, D. H., Ecological linkages between aboveground and belowground biota. Science 2004, 304, (5677), 1629-1633. 3. Bardgett, R. D.; van der Putten, W. H., Belowground biodiversity and ecosystem functioning. Nature 2014, 515, (7528), 505-511. 4. van der Heijden, M. G. A.; Bardgett, R. D.; van Straalen, N. M., The unseen majority: soil microbes as drivers of plant diversity and productivity in terrestrial ecosystems. Ecol Lett 2008, 11, (3), 296-310. 5. Torsvik, V.; Ovreas, L., Microbial diversity and function in soil: from genes to ecosystems. Curr Opin Microbiol 2002, 5, (3), 240-5. 6. DeLong, E. F.; Preston, C. M.; Mincer, T.; Rich, V.; Hallam, S. J.; Frigaard, N. U.; Martinez, A.; Sullivan, M. B.; Edwards, R.; Brito, B. R.; Chisholm, S. W.; Karl, D. M., Community genomics among stratified microbial assemblages in the ocean's interior. Science 2006, 311, (5760), 496-503. 7. Torsvik, V.; Daae, F. L.; Sandaa, R. A.; Ovreas, L., Novel techniques for analysing microbial diversity in natural and perturbed environments. J Biotechnol 1998, 64, (1), 53-62. 8. Sogin, M. L.; Morrison, H. G.; Huber, J. A.; Mark Welch, D.; Huse, S. M.; Neal, P. R.; Arrieta, J. M.; Herndl, G. J., Microbial diversity in the deep sea and the underexplored "rare biosphere". P Natl Acad Sci USA 2006, 103, (32), 12115- 12120. 9. Miki, T.; Yokokawa, T.; Matsui, K., Biodiversity and multifunctionality in a microbial community: a novel theoretical approach to quantify functional redundancy. Proc Biol Sci 2014, 281, (1776), 20132498. 10. Hillebrand, H.; Matthiessen, B., Biodiversity in a complex world: consolidation and progress in functional biodiversity research. Ecol Lett 2009, 12, (12), 1405- 1419. 11. Peter, H.; Beier, S.; Bertilsson, S.; Lindstrom, E. S.; Langenheder, S.; Tranvik, L. J., Function-specific response to depletion of microbial diversity. ISME J 2011, 5, (2), 351-61. 12. Maestre, F. T.; Delgado-Baquerizo, M.; Jeffries, T. C.; Eldridge, D. J.; Ochoa, V.; Gozalo, B.; Quero, J. L.; Garcia-Gomez, M.; Gallardo, A.; Ulrich, W.; Bowker, M. A.; Arredondo, T.; Barraza-Zepeda, C.; Bran, D.; Florentino, A.; Gaitan, J.; Gutierrez, J. R.; Huber-Sannwald, E.; Jankju, M.; Mau, R. L.; Miriti, M.; Naseri, K.; Ospina, A.; Stavi, I.; Wang, D. L.; Woods, N. N.; Yuan, X.; Zaady, E.; Singh, B. K., Increasing aridity reduces soil microbial diversity and abundance in global drylands. P Natl Acad Sci USA 2015, 112, (51), 15684-15689. 13. Mader, P.; Fliessbach, A.; Dubois, D.; Gunst, L.; Fried, P.; Niggli, U., Soil fertility and biodiversity in organic farming. Science 2002, 296, (5573), 1694-1697. 14. Gans, J.; Wolinsky, M.; Dunbar, J., Computational improvements reveal great bacterial diversity and high metal toxicity in soil. Science 2005, 309, (5739), 1387- 1390. 15. Bodelier, P. L. E., Toward understanding, managing, and protecting microbial ecosystems. Front Microbiol 2011, 2. 16. Bell, T.; Newman, J. A.; Silverman, B. W.; Turner, S. L.; Lilley, A. K., The contribution of species richness and composition to bacterial services. Nature 2005, 436, (7054), 1157-1160. 17. Allison, S. D.; Martiny, J. B. H., Resistance, resilience, and redundancy in microbial communities. P Natl Acad Sci USA 2008, 105, 11512-11519. 18. Wertz, S.; Degrange, V.; Prosser, J. I.; Poly, F.; Commeaux, C.; Freitag, T.; Guillaumaud, N.; Roux, X. L., Maintenance of soil functioning following erosion of microbial diversity. Environ Microbiol 2006, 8, (12), 2162-9. 19. Wertz, S.; Degrange, V.; Prosser, J. I.; Poly, F.; Commeaux, C.; Guillaumaud, N.; Le Roux, X., Decline of soil microbial diversity does not influence the resistance and resilience of key soil microbial functional groups following a model disturbance. Environmental Microbiology 2007, 9, (9), 2211-2219. 20. Reich, P. B.; Tilman, D.; Isbell, F.; Mueller, K.; Hobbie, S. E.; Flynn, D. F. B.; Eisenhauer, N., Impacts of Biodiversity Loss Escalate Through Time as Redundancy Fades. Science 2012, 336, (6081), 589-592. 21. Cardinale, B. J.; Matulich, K. L.; Hooper, D. U.; Byrnes, J. E.; Duffy, E.; Gamfeldt, L.; Balvanera, P.; O'Connor, M. I.; Gonzalez, A., The Functional Role of Producer Diversity in Ecosystems. Am J Bot 2011, 98, (3), 572-592. 22. Reinthaler, T.; Winter, C.; Herndl, G. J., Relationship between bacterioplankton richness, respiration, and production in the southern North Sea. Appl Environ Microb 2005, 71, (5), 2260-2266. 23. Singh, B. K.; Quince, C.; Macdonald, C. A.; Khachane, A.; Thomas, N.; Al-Soud, W. A.; Sorensen, S. J.; He, Z.; White, D.; Sinclair, A.; Crooks, B.; Zhou, J.; Campbell, C. D., Loss of microbial diversity in soils is coincident with reductions in some specialized functions. Environ Microbiol 2014, 16, (8), 2408-20. 24. Obernosterer, I.; Lami, R.; Larcher, M.; Batailler, N.; Catala, P.; Lebaron, P., Linkage Between Bacterial Carbon Processing and the Structure of the Active Bacterial Community at a Coastal Site in the NW Mediterranean Sea. Microb Ecol 2010, 59, (3), 428-435. 25. Gamfeldt, L.; Hillebrand, H.; Jonsson, P. R., Multiple functions increase the importance of biodiversity for overall ecosystem functioning. Ecology 2008, 89, (5), 1223-1231. 26. Jung, J.; Philippot, L.; Park, W., Metagenomic and functional analyses of the consequences of reduction of bacterial diversity on soil functions and bioremediation in diesel-contaminated microcosms. Scientific reports 2016, 6, 23012. 27. Chapin, F. S.; Zavaleta, E. S.; Eviner, V. T.; Naylor, R. L.; Vitousek, P. M.; Reynolds, H. L.; Hooper, D. U.; Lavorel, S.; Sala, O. E.; Hobbie, S. E.; Mack, M. C.; Diaz, S., Consequences of changing biodiversity. Nature 2000, 405, (6783), 234-242. 28. Loreau, M., Does functional redundancy exist? Oikos 2004, 104, (3), 606-611. 29. Wardle, D. A.; Bardgett, R. D.; Callaway, R. M.; Van der Putten, W. H., Terrestrial Ecosystem Responses to Species Gains and Losses. Science 2011, 332, (6035), 1273-1277. 30. Philippot, L.; Spor, A.; Henault, C.; Bru, D.; Bizouard, F.; Jones, C. M.; Sarr, A.; Maron, P. A., Loss in microbial diversity affects nitrogen cycling in soil. ISME J 2013, 7, (8), 1609-19. 31. Vitousek, P. M.; DAntonio, C. M.; Loope, L. L.; Westbrooks, R., Biological invasions as global environmental change. Am Sci 1996, 84, (5), 468-478. 32. Kennedy, T. A.; Naeem, S.; Howe, K. M.; Knops, J. M. H.; Tilman, D.; Reich, P., Biodiversity as a barrier to ecological invasion. Nature 2002, 417, (6889), 636-638. 33. Levine, J. M., Species diversity and biological invasions: Relating local process to community pattern. Science 2000, 288, (5467), 852-854. 34. Ricciardi, A.; MacIsaac, H. J., The ecology of invasions by animals and plants. Nature 2008, 452, (7183), 34-34. 35. Elton, C., The ecology of invasions by animals and plants. London: University of Chicago Press. 1958, 196pp. 36. van Elsas, J. D.; Hill, P.; Chronakova, A.; Grekova, M.; Topalova, Y.; Elhottova, D.; Kristufek, V., Survival of genetically marked Escherichia coli O157:H7 in soil as affected by soil microbial community shifts. ISME J 2007, 1, (3), 204-14. 37. van Elsas, J. D.; Chiurazzi, M.; Mallon, C. A.; Elhottova, D.; Kristufek, V.; Salles, J. F., Microbial diversity determines the invasion of soil by a bacterial pathogen. Proc Natl Acad Sci U S A 2012, 109, (4), 1159-64. 38. Yan, Y.; Kuramae, E. E.; Klinkhamer, P. G.; van Veen, J. A., Revisiting the dilution procedure used to manipulate microbial biodiversity in terrestrial systems. Appl Environ Microbiol 2015, 81, (13), 4246-52. 39. Vivant, A. L.; Garmyn, D.; Maron, P. A.; Nowak, V.; Piveteau, P., Microbial Diversity and Structure Are Drivers of the Biological Barrier Effect against Listeria monocytogenes in Soil. PLoS ONE 2013, 8, (10), e76991. 40. Irikiin, Y.; Nishiyama, M.; Otsuka, S.; Senoo, K., Rhizobacterial community-level, sole carbon source utilization pattern affects the delay in the bacterial wilt of tomato grown in rhizobacterial community model system. Appl Soil Ecol 2006, 34, (1), 27-32. 41. Matos, A.; Kerkhof, L.; Garland, J. L., Effects of microbial community diversity on the survival of Pseudomonas aeruginosa in the wheat rhizosphere. Microb Ecol 2005, 49, (2), 257-264. 42. Martinez, J. L., Antibiotics and antibiotic resistance genes in natural environments. Science 2008, 321, (5887), 365-367. 43. Walsh, C., Where will new antibiotics come from? Nature Reviews Microbiology 2003, 1, (1), 65-70. 44. Pruden, A.; Larsson, D. G.; Amezquita, A.; Collignon, P.; Brandt, K. K.; Graham, D. W.; Lazorchak, J. M.; Suzuki, S.; Silley, P.; Snape, J. R.; Topp, E.; Zhang, T.; Zhu, Y. G., Management options for reducing the release of antibiotics and antibiotic resistance genes to the environment. Environmental health perspectives 2013, 121, (8), 878-85. 45. Bush, K.; Courvalin, P.; Dantas, G.; Davies, J.; Eisenstein, B.; Huovinen, P.; Jacoby, G. A.; Kishony, R.; Kreiswirth, B. N.; Kutter, E.; Lerner, S. A.; Levy, S.; Lewis, K.; Lomovskaya, O.; Miller, J. H.; Mobashery, S.; Piddock, L. J.; Projan, S.; Thomas, C. M.; Tomasz, A.; Tulkens, P. M.; Walsh, T. R.; Watson, J. D.; Witkowski, J.; Witte, W.; Wright, G.; Yeh, P.; Zgurskaya, H. I., Tackling antibiotic resistance. Nature reviews. Microbiology 2011, 9, (12), 894-6. 46. Zhu, Y.-G.; Johnson, T. A.; Su, J.-Q.; Qiao, M.; Guo, G.-X.; Stedtfeld, R. D.; Hashsham, S. A.; Tiedje, J. M., Diverse and abundant antibiotic resistance genes in Chinese swine farms. P Natl Acad Sci USA 2013, 110, (9), 3435-3440. 47. Zhu, Y. G.; Zhao, Y.; Li, B.; Huang, C. L.; Zhang, S. Y.; Yu, S.; Chen, Y. S.; Zhang, T.; Gillings, M. R.; Su, J. Q., Continental-scale pollution of estuaries with antibiotic resistance genes. Nat Microbiol 2017, 2, 16270. 48. Su, J.-Q.; Wei, B.; Ou-Yang, W.-Y.; Huang, F.-Y.; Zhao, Y.; Xu, H.-J.; Zhu, Y.-G., Antibiotic resistome and its association with bacterial communities during sewage sludge composting. Environmental science & technology 2015, 49, (12), 7356-7363. 49. Chen, Q.; An, X.; Li, H.; Su, J.; Ma, Y.; Zhu, Y. G., Long-term field application of sewage sludge increases the abundance of antibiotic resistance genes in soil. Environment international 2016, 92-93, 1-10. 50. Cole, J. R.; Wang, Q.; Cardenas, E.; Fish, J.; Chai, B.; Farris, R. J.; Kulam-Syed- Mohideen, A. S.; McGarrell, D. M.; Marsh, T.; Garrity, G. M.; Tiedje, J. M., The Ribosomal Database Project: improved alignments and new tools for rRNA analysis. Nucleic Acids Res 2009, 37, D141-D145. 51. Fierer, N.; Ladau, J.; Clemente, J. C.; Leff, J. W.; Owens, S. M.; Pollard, K. S.; Knight, R.; Gilbert, J. A.; McCulley, R. L., Reconstructing the microbial diversity and function of pre-agricultural tallgrass prairie soils in the United States. Science 2013, 342, (6158), 621-4. 52. Caporaso, J. G.; Kuczynski, J.; Stombaugh, J.; Bittinger, K.; Bushman, F. D.; Costello, E. K.; Fierer, N.; Pena, A. G.; Goodrich, J. K.; Gordon, J. I.; Huttley, G. A.; Kelley, S. T.; Knights, D.; Koenig, J. E.; Ley, R. E.; Lozupone, C. A.; McDonald, D.; Muegge, B. D.; Pirrung, M.; Reeder, J.; Sevinsky, J. R.; Turnbaugh, P. J.; Walters, W. A.; Widmann, J.; Yatsunenko, T.; Zaneveld, J.; Knight, R., QIIME allows analysis of high-throughput community sequencing data. Nature methods 2010, 7, (5), 335-6. 53. Edgar, R. C., Search and clustering orders of magnitude faster than BLAST. Bioinformatics 2010, 26, (19), 2460-1. 54. Hammesfahr, U.; Heuer, H.; Manzke, B.; Smalla, K.; Thiele-Bruhn, S., Impact of the antibiotic sulfadiazine and pig manure on the microbial community structure in agricultural soils. Soil Biol Biochem 2008, 40, (7), 1583-1591. 55. Chee-Sanford, J. C.; Mackie, R. I.; Koike, S.; Krapac, I. G.; Lin, Y. F.; Yannarell, A. C.; Maxwell, S.; Aminov, R. I., Fate and transport of antibiotic residues and antibiotic resistance genes following land application of manure waste. Journal of environmental quality 2009, 38, (3), 1086-108. 56. Heuer, H.; Focks, A.; Lamshoeft, M.; Smalla, K.; Matthies, M.; Spiteller, M., Fate of sulfadiazine administered to pigs and its quantitative effect on the dynamics of bacterial resistance genes in manure and manured soil. Soil Biol Biochem 2008, 40, (7), 1892-1900. 57. Fargione, J. E.; Tilman, D., Diversity decreases invasion via both sampling and complementarity effects. Ecol Lett 2005, 8, (6), 604-611. 58. Tilman, D., The ecological consequences of changes in biodiversity: A search for general principles. Ecology 1999, 80, (5), 1455-1474. 59. Goodman, A. L.; Kallstrom, G.; Faith, J. J.; Reyes, A.; Moore, A.; Dantas, G.; Gordon, J. I., Extensive personal human gut microbiota culture collections characterized and manipulated in gnotobiotic mice. Proc Natl Acad Sci U S A 2011, 108, (15), 6252-7. 60. Case, T. J., Invasion resistance arises in strongly interacting species-rich model competition communities. Proc Natl Acad Sci U S A 1990, 87, (24), 9610-4. 61. Hodgson, D. J.; Rainey, P. B.; Buckling, A., Mechanisms linking diversity, productivity and invasibility in experimental bacterial communities. Proc Biol Sci 2002, 269, (1506), 2277-83. 62. Eisenhauer, N.; Bowker, M. A.; Grace, J. B.; Powell, J. R., From patterns to causal understanding: structural equation modeling (SEM) in soil ecology. Pedobiologia 2015, 58, (2), 65-72. 63. Hu, H. W.; Wang, J. T.; Li, J.; Shi, X. Z.; Ma, Y. B.; Chen, D.; He, J. Z., Long-Term Nickel Contamination Increases the Occurrence of Antibiotic Resistance Genes in Agricultural Soils. Environ Sci Technol 2017, 51, (2), 790-800.