A sandy loam soil sample was collected from the 0–10 cm soil layer with a previous history of oil contamination in an oil field of Shanxi, stored at 4°C and transferred to laboratory. The main characteristics of soil sample was provided in Supplementary Table S1.

The DBT degrader was isolated from sandy loam soil samples contaminated with oil and purified using liquid and solid mineral salt medium with increasing DBT concentration (100, 200, 300, and 400 mg/L) according to the method of Jiang et al. (2022). The isolated cultivable DBT degrader was identified as Rhodococcus sp. CGMCC 24903 and preserved in the ecological engineering laboratory of China Agricultural University and China General Microbiological Culture Collection Center. Prior to experimentation, the strain was cultured in liquid LB medium at 180 rpm for 3 days to 10⁹ cells/mL determined by the microplate method (Estrada-Bonilla et al., 2017), and subsequently centrifuged bacterial cell was used as inoculum for further study on DBT degradation microcosms.

Three treatments with 1,000 g soil were carried out to study the effect of Rhodococcus bioaugmentation and biostimulation on DBT biodegradation performance in soils and experiment lasted 40 days (Liu et al., 2023). The natural degradation treatment (NAT) consisted of original soil. The biological augmentation treatment (BIOE) involved soil combined with Rhodococcus sp. CGMCC 24903 inoculation at the total populations of 10⁹ cells/g soil. The biostimulation treatment (BIOS) included soil with optimizing nutrient levels by adding nutrient solution, e.g., CH₄N₂O and NaH₂PO₃ to a final C:N:P ratio of 100:10:1, with concentrations of C at 24.36 mg/g, N at 2.427 mg/g, and P at 0.24 mg/g. (Xu and Lu, 2010). The soil used in Microcosms was air-dried, and particles larger than 2 mm were removed to achieve soil homogeneity. DBT in the three treatments was set at a final concentration of 100 mg/kg. The microcosms were incubated at 25°C ± 1°C for 40 days in the dark in triplicates, with soil moisture maintained at 20% through periodic addition of distilled water. 50 g of soils were sampled on day 0, 10, 20, 30, and 40. All samples were stored at −20°C until further physicochemical and DNA analysis.

All solvents used were laboratory grade chemicals obtained from Beijing Chemical Company (China). Dibenzothiophene (CAS 132-65-0, purity >99%) were obtained from Sigma Aldrich. Deionized water was used to prepare all media and solutions.

DBT analysis was performed by extracting 5 g of soil using ethyl acetate at a ratio of 1 mL per Gram of soil. Each soil-solvent mixture was sonicated for 1 h, followed by nitrogen drying and resuspension of the residue in 4 mL of the same solvent. Prior to GC-MS analyses (Agilent 7890A-5975C) using a TG-5SILMS column (30 × 0.25 mm, 0.25 μm), samples were filtered (0.22 μm filter, Pall). The GC oven temperature was programmed to increase from 50°C to 250°C at a rate of 15°C per minute, and held at 250°C for 2 min. The injector temperature was set to 250°C. The DBT determination were conducted in triplicate, and the mean values were taken as the results.

Total DNA of soil samples was extracted using the FastDNA spin Kit for soil (MP Biomedicals, Santa Ana, Carlsbad, CA, United States). PCR amplification of 16S rRNA gene fragments was performed using the primers 338F (5′-ACT​CCT​ACG​GGA​GGC​AGC​AG-3′) and 806R (5′-GGACTACHVGGGTWTCTAAT-3′). The PCR products were sequenced in a 2 × 300-bp paired-end format on the Illumina MiSeq platform, which was conducted by Majorbio BioPharm Technology Co. Ltd. (Shanghai, China). The sequences were deposited in the NCBI database under the project accession number PRJNA914518.

The raw 16S rRNA gene sequencing reads were demultiplexed and quality-filtered using fastp version 0.20.0 (Chen et al., 2018), and then merged using FLASH version 1.2.7 (Magoč and Salzberg, 2011). Operational taxonomic units (OTUs) with a 97% similarity cutoff were clustered using UPARSE version 7.1 (Stackebrandt and Goebel, 1994; Edgar, 2013), and chimeric sequences were identified and removed. The taxonomic classification of each OTU representative sequence was determined using RDP Classifier version 2.2 (Wang et al., 2007) with a confidence threshold of 0.7, based on the silva reference database (Quast et al., 2012).

The basic indicators were visualized using Origin Pro 2022. Heatmaps were generated using the pheatmap package in R 4.1.3. Bacterial diversity metrics, including Shannon diversity, Pielou’s evenness, Chao1 index, and Faith’s phylogenetic diversity, were estimated using the Picante and vegan packages. Principal coordinates analysis (PCoA) was conducted based on Bray-Curtis dissimilarity index matrices (Sheng et al., 2023).

Correlation network analysis was conducted to ascertain the role of core microorganisms in DBT degradation. Based on Pearson correlation coefficients, bacterial taxonomic groups at the OTU level that exhibited significant correlations with DBT degradation (p < 0.05, r > 0.97) were selected as core microorganisms, and the visualization was carried out using Cytoscape 3.9.1. Co-occurrence networks at the OTU level were conducted using Molecular Ecology Network Analysis with a cutoff value of 0.85. The networks were visualized using Gephi 0.9.2. Keystones, including connectors, network hubs, and module hubs, were categorized based on within-module connectivity (Zi) and among-module connectivity (Pi) in the network, following the classification proposed by (Deng et al., 2012).

Gas chromatographic analysis revealed that the degradation efficiency of DBT in BIOE was significantly higher than that in NAT and BIOS over a period of 40 days (Figure 1). Within the first 20 days of soil incubation, the DBT degradation efficiency in BIOE and BIOS was significantly higher compared to natural degradation, suggesting that both Rhodococcus bioaugmentation and biostimulation could enhance the initial biodegradation ability of soil microbes. In the subsequent 20 days of incubation, though bioaugmentation group still had the highest degradation rate and increased to ∼60% at day 40, the DBT degradation performance in BIOS was lower than that in NAT. After day 30, the degradation level of the three treatments was gradually flatten, suggesting that the majority of DBT degradation occurred within the first 30 days. Biostimulation had only a positive impact at the early stage of soil incubation, which could be attributed to the proliferation of microorganisms capable of degrading DBT that were promoted by the additional nutrients. However, during the later stage of soil incubation, the nutrients were depleted and the dominant species might be changed to unfavorable conditions for DBT-degrading microorganisms, which reflected the complex interaction between microbial populations and nutrient availability (Fuhrman et al., 2015). The above results suggested that the addition of inoculum (Rhodococcus) for bioaugmentation significantly stimulated the biodegradation of DBT and improved the bioremediation performance in petroleum-contaminated soils.

After the optimization of sequencing, a total of 1,125,468 bacterial reads of the 16S rRNA gene were obtained, which were subsequently clustered into 1,362 bacterial operational taxonomic units (OTUs). Among these OTUs, the NAT, BIOE, and BIOS treatments had 119, 6 and 42 unique OTUs, respectively, and had 344 OTUs shared among all treatments (Figure 2A). Firmicutes and Actinobacteria were identified as the dominant phyla in all treatments (Figure 2B). Notably, the addition of nutrients increased the relative abundance of Firmicutes, while reduced the levels of Actinobacteria. To gain deeper insights into the dynamic changes within the bacterial community, a heatmap analysis was performed on the 20 most abundant genera (Figure 2C). Based on the relative abundance, the genera were classified into two clusters. Bacillus, Micromonospora, Streptomyces, Bhargavaea, and Chungangia formed one cluster, while Sporosarcina, Sphingobacterium, Paracoccus, and others clustered together, indicating potential functional similarities or ecological interactions among these taxa (Bao et al., 2023). Bacillus and Streptomyces were dominant genera during soil incubation, which had the ability to degrade dibenzothiophene through their participation in the 4S pathway (Buzanello et al., 2014). This finding highlights their potential as key DBT-degrading bacteria in the remediation of petroleum-contaminated soil. Surprisingly, Rhodococcus genus was not shown in the heatmap even though it was inoculated in the BIOE, which might be due to its low relative abundance in soil and its colonization related to the concentration of PAHs (Cunliffe and Kertesz, 2006).

The alpha diversity of soil bacterial community during the incubation periods was shown in Figure 3A. The Shannon diversity of the NAT treatment was significantly higher than that of the BIOE and BIOS treatments on the 10th day of incubation. However, after the 40 days of incubation, there were no significant differences in Shannon diversity among the three treatments. During the incubation period, the BIOE treatment exhibited significantly lower Chao1 richness and phylogenetic diversity compared to the other treatments. However, the Simpson evenness of the BIOE treatment was significantly higher than that of the other treatments. The BIOS treatment showed higher Simpson evenness but lower Chao1 richness and phylogenetic diversity compared to the NAT treatment during the incubation period. The above results indicated that Rhodococcus bioaugmentation increased the dominance of species in petroleum-contaminated soil while reducing the richness of total bacterial community. principal coordinates analysis (PCoA) based on Bray-Curtis distance matrix was conducted to visualize the succession patterns of soil bacterial communities in diverse groups (Figure 3B). The first two principal components accounted for 74.7% of the community variance. Notably, the bacterial communities during incubation in all samples were divided into three classes, indicating dissimilar bacterial composition in the microcosms of NAT, BIOE and BIOS. The results highlighted the significant impact of different treatments on the beta diversity of bacterial communities, especially the critical role of nutrient addition and inoculation on bacterial community composition in the degradation process of DBT during soil incubation (Strickland et al., 2009; Wang et al., 2018).

The degradation of DBT was inevitably influenced by the dynamics of the bacterial community and their interaction during soil incubation process. Here, network analysis was employed to identify the key microbial taxa at the OTU level that were associated with DBT degradation, which were screened based on significant correlations with DBT degradation (p < 0.05). In NAT, BIOE, and BIOS groups, there were 15, 1, and 7 nodes, respectively, significantly participating in the DBT degradation process (Figure 4). The OTU1291 (Shimazuella), OTU1532 (Paenibacillus) in NAT, OTU 1851 (Bacillus) in BIOE and OTU1292 (Bacillus), OTU1322 (Bacillus) in BIOS were not only the core OTUs involved in DBT degradation but also belonged to the top 20 bacterial genera during soil incubation (Figure 2C; Supplementary Table S2). Previous studies also reported that Shimazuella and Bacillus had strong capabilities in DBT degradation (Buzanello et al., 2014).

The bacterial interactions under different treatments were investigated by constructing co-occurrence networks based on Pearson correlations among OTUs (Figure 5A). All co-occurrence networks demonstrated general network properties of scale-free, small-world, and modularity (Supplementary Table S3), with the network topology following a power-law distribution (R ² = 0.886, 0.756, 0.804 in NAT, BIOE and BIOS, respectively). Compared to other treatments, BIOE network showed fewer nodes (163) and links (606), but displayed higher average degree (avgK) and average clustering coefficient clustering coefficient (avgCC), as well as lower average path distance (GD) and average clustering coefficient (HD). These findings suggested that the BIOE network was more tightly connected, and external disturbances could quickly propagate throughout the entire network. Compared to NAT, the BIOS network also exhibits tight connectivity based on its topological characteristics, with the highest number of links (2467). The correlations between bacterial species in the network varied, with NAT showing a higher percentage of negative correlations (45.85%) compared to BIOE (16.17%) and BIOS (35.63%). These findings collectively indicated significant differences in the networks among the three treatments, and more bacteria in BIOE and BIOS groups showed strongly cooperative symbiosis to the environmental pressure (DBT pollution) (Zhang et al., 2020).

Based on the analysis of within-module connectivity (Zi) and among-module connectivity (Pi) (Olesen et al., 2007), a total of 14, 2, and 4 module hubs were observed in the NAT, BIOE, and BIOS networks, respectively, along with 6, 4, and 8 connectors (Figure 5B; Supplementary Table S4). However, no network hubs were identified in these groups. More key nodes in the NAT network indicated its higher complexity. All these nodes represented keystones that played crucial roles in the network. The 20 keystones in the NAT network mainly belonged to Firmicutes, Proteobacteria, and Actinobacteriota, the 6 keystones in the BIOE network belonged to Firmicutes, Actinobacteriota, and Proteobacteria, and the 12 keystones in the BIOS network were predominantly from Myxococcota, Firmicutes, and Actinobacteriota. Members of Firmicutes and Actinobacteriota were the predominant keystone taxa during the incubation process of petroleum-contaminated soil, which accounted for approximately 70% of all module hubs and connectors. Furthermore, key microorganisms belonging to the genera Paenibacillus (OTU1532, OTU861, OTU 2001, OTU1317, OTU1325, OTU206, OTU1390) and Bacillus (OTU1020, OTU1696) related to DBT degradation were also identified as keystones in bacterial interaction network, suggesting that these microorganisms might be indigenous bacteria with a potential function of DBT degradation by regulating core bacterial collaborative symbiosis in the co-occurrence network. The above results indicated that Rhodococcus bioaugmentation (BIOE) was the simplest and effective method to enhance the clustering degree of bacterial community network compared to natural degradation for enhancing DBT degradation in petroleum-contaminated soil. Additionally, nutrient addition (BIOS) reduces the complexity of bacterial community in petroleum-contaminated soil without obviously enhancing DBT degradation capability.