Impact of Dietary Protein Content on Soil Bacterial and Fungal Communities in a Rice–Crab Co-culture System

Although co-culture of paddy fields with aquatic animals is lucrative, the ecological impacts of high-protein content entering the agricultural soil via animal pellet feed and feces have not been well studied. Moreover, the effects of dietary protein on soils and soil microbes remain unclear. To elucidate this, we examined soil bacterial and fungal community composition and temporal changes in paddy fields subjected to different protein-content diets via 16S/18S rRNA gene amplicon sequencing analysis with a high-throughput next-generation sequencer. MiSeq sequencing revealed that protein content significantly impacted fungal community structure. High-protein diets reduced bacterial community diversity and relative abundance in both July and October. The phylum-level bacterial taxonomic composition was not affected by diet treatment, while in fungi, a major phylum-level shift was evident. Hierarchically clustered analysis showed that high-protein diets significantly reduced the relative abundance of Brevundimonas in both July and October. Saprotrophic macrofungal diversity was negatively related to dietary protein content. Considering microbial community structure and environmental factors, ca. 15% protein content is appropriate for the rice-crab co-culture system that we studied.


INTRODUCTION
Rice is a staple food for over half of the global population, and its production has increased in the last two decades due to technological advances, including the use of pesticides and chemical fertilizers. Although improving grain yield and rice quality are important goals for many governments, minimizing chemical fertilizer and pesticide inputs is important for longterm agricultural and environmental sustainability (Yao et al., 2017;Stuart et al., 2018). Rice fields can provide a habitat mosaic of temporary and more permanent waters and sustain populations of various aquatic species, including carp, crabs, crayfish, soft-shelled turtles, and frogs (Lawler, 2001). Rice-aquatic-animal co-culture has been widely studied; such systems can contribute to the ecological intensification of agriculture by providing multiple ecosystem services, promoting biological pest control, reducing the use of pesticides, improving soil quality, and enhancing crop yields (Xie et al., 2011;Hu et al., 2016;Sha et al., 2017;Wan et al., 2019). Moreover, such systems increase protein production for human consumption and reduce water use and environmental pollution (Bashir et al., 2020). Thus, co-culture systems have gained increasing attention in recent years because of their benefits.
In China, the government has established many projects to help farmers develop diverse paddy field co-culture models, including rice-duck, rice-fish, rice-crayfish, and rice-crab. However, high-protein diets widely used in pond aquaculture are also used in rice-aquatic animal co-culture. Undigested protein-rich feed and feces are excreted into the bottom soil of paddy fields, potentially causing severe nitrogen pollution of the soil and irrigation water. However, studies have shown that rice-animal co-culture can enhance nutrient use efficiency and reduce nutrient loss to the environment via the complementary use of nutrients by aquatic animals and rice (Hu et al., 2020). Nonetheless, the relationships between high-protein diets and paddy field soil conditions in co-culture systems remain unclarified.
Rice-crab co-culture is the characteristic rice-based ecological aquaculture system used in northern China (Song et al., 2019;Wang et al., 2019). It improves soil nutrient levels and increases soil-to-rice nutrient translocation capacity (Bashir et al., 2021). Although high-throughput sequencing (e.g., Illumina) is increasingly used to analyze rice soil microbial diversity and community structure, the microbial response to different rice-crab co-culture systems is yet to be adequately studied. Determining soil bacterial community structure in rice-crab co-culture systems can provide key initiatives for improving aquaculture. Moreover, the soil bacterial community structure has a strong influence on the bacterial microbiota of crabs. The relationships between rice-crab co-culture soil microbiota communities and dietary protein content have not previously been studied. Therefore, we examined these dynamics and aimed to provide a theoretical reference for these relationships.

Ethics Statement
Our study did not involve endangered or protected species. In China, breeding and catching Chinese mitten crabs, Eriocheir sinensis, in rice fields does not require specific permits. All efforts were made to minimize animal suffering and discomfort. The experimental protocol was approved by the Animal Ethics Committee of Shenyang Agriculture University.

Field Description and Experimental Design
The experiments were performed in 12 paddy fields (each 6 × 7 m 2 ) at Panjin Guanghe Crab Industry Co., Ltd., Panjin, Liaoning Province, China, during the 2020 rice-growing season (over 5 months). In total, 30,000 late megalopa-stage Chinese mitten crabs (average weight was 6.25 ± 1.21 mg) were obtained from Panjin Guanghe Crab Industry and were randomly distributed in the 12 fields (2,500 crabs/field). Three different protein-content diets were formulated (15, 30, and 45% crude protein). As a control, no supplemental food was provided ( Table 1). During the experiment, crabs were fed at 10% of body weight per day (at 08:00, 12:00, and 18:00) for 5 months.

Sample Collection
Soil samples were collected on July 15, 2020, during the tillering stage (15 days after crab introduction), and October 1, during the ripening stage (100 days after seed sowing). For each plot, five soil samples (diameter 2.5 cm and depth 0-20 cm) were pooled to obtain one biological replicate. After removing rice roots and stones, the pooled soil sample was placed in a sterile plastic bag in an icebox until transport to a laboratory for further sample storage and pH analysis. pH was quantified using a PHS-3C pH meter (Shanghai, China) with 1:2.5 vol soil/H 2 O solutions.

DNA Extraction and Purification
Soil samples were homogenized, and a 0.8-g subsample of soil was used for total genomic DNA extraction. DNA was extracted using the OMEGA Soil DNA kit (Omega Bio-Tek, Norcross, GA) according to the manufacturer's instructions and stored at -20 • C before further analysis. The quantity and quality of the extracted DNA were evaluated using a NanoDrop ND-1000 spectrophotometer (Thermo Fisher Scientific, Waltham, MA) and agarose gel electrophoresis, respectively.
Thermal cycling consisted of initial denaturation at 98 • C for 5 min, followed by 28 cycles of denaturation at 98 • C for 30 s, annealing at 55 • C for 30 s, and extension at 72 • C for 45 s, with a final extension of 5 min at 72 • C. Amplicons were extracted from 2% agarose gels, purified using Vazyme VAHTS DNA Clean Beads (Vazyme, Nanjing, China) according to the manufacturer's instructions, and quantified using the Quant-iT PicoGreen dsDNA Assay Kit (Invitrogen, Carlsbad, CA). The soil sample bacterial community diversity and composition were analyzed based on the raw sequencing data obtained using the Illumina MiSeq platform at Shanghai Personal Biotechnology Co., Ltd. (Shanghai, China), according to standard protocols.

Bioinformatics and Statistical Analysis
After sequencing, raw FASTQ data were processed using QIIME2 and R v. 3.2.0, with slight modifications according to the official tutorials 1 . Briefly, the raw sequence data were demultiplexed, quality filtered, denoised, and merged, and the chimeras were removed using the DADA2 plugin. Amplicon sequencing variant (ASV)-level alpha diversity indices (Chao1 richness estimator, observed species, Shannon diversity index, Simpson index, Faith's phylogenetic diversity index, Pielou's evenness, and Good's coverage) were calculated using the ASV table in QIIME2 and were visualized using box plots. Taxonomic composition and abundance were investigated using MEGAN and GraPhlAn, respectively. A Venn diagram was generated using the R package "VennDiagram" to visualize the shared and unique ASVs among samples or groups, based on the occurrence of ASVs across samples and groups, regardless of their relative abundance. The raw reads have been deposited with the NCBI (BioProject number PRJNA693650).

Analysis of Pyrosequencing Data
In total, 55,522 bacterial and 75,786 fungal non-singleton sequences were obtained from all 24 samples (Supplementary  Figures 1, 2). After sequence processing and quality filtering, 88,281 bacterial OTUs and 4,356 fungal OTUs were identified ( Supplementary  Figures 3, 4). The observed and estimated community structure indices were higher for bacterial than for fungal OTUs (Supplementary Table 2). The pH data of those soil samples are shown in Supplementary Table 3.

Bacterial and Fungal Community Structure
The ten most abundant soil bacterial phyla were Proteobacteria, Chloroflexi, Bacteroidetes, Acidobacteria, Actinobacteria, Firmicutes, Gemmatimonadetes, Nitrospirae, Patescibacteria, and Verrucomicrobia (Supplementary Figure 5A). The phylumlevel bacterial taxonomic composition was not affected by diet treatment, while in fungi, a major phylum level shift was evident ( Figure 1A). Five bacterial phyla (Proteobacteria, Chloroflexi, Bacteroidetes, Acidobacteria, and Actinobacteria) accounted for 90% of the sequences. Regarding the fungal communities, 10 phyla were observed in the soil samples: Ascomycota, Basidiomycota, Mortierellomycota, Rozellomycota, Chytridiomycota, Basidiobolomycota, Zoopagomycota, Blastocladiomycota, Mucoromycota, and Aphelidiomycota (Supplementary Figure 5B). Fungal community composition underwent major phylum-level shifts between July and October ( Figure 1B): Ascomycota dominated in October (especially in the control group), and Basidiomycota comprised a large proportion of the fungal community in the OTD_15, _30, and _45 groups. Table 2 shows the bacterial and fungal within-habitat (alpha) diversity. According to the Chao1 analysis. The bacterial richness was highest in the JTD_0 group and lowest in the OTD_45 group. The rarefaction curves show that bacterial community diversity was lower in the diet treatment groups than in the control group (Supplementary Figure 6). However, soil bacterial diversity and evenness were not different between samples obtained under the different treatments in July and October (Supplementary Figure 7). Regarding fungal communities, the Chao1, Shannon, and Simpson indices were higher in the JTD_30 and OTD_30 groups than in the other groups (Supplementary Figure 8).

Microbial Alpha and Beta Diversity
Between-habitat (beta) diversity was evaluated via Principal coordinate analysis (PCoA), Non-metric Multidimensional scaling (NMDS), and cluster analysis. The PCoA based on Bray-Curtis distances accounted for 30.7% of the bacterial community variability and 47.9% of the fungal community variability (Figure 2A,B). For the fungal community, the July samples were separated from October samples in the PCoA (Figure 2B). According to the NMDS analysis, the bacterial and fungal communities had the same stress values (0.126 and 0.157, respectively); only fungal communities were different between July and October soils (Figure 2C,D).

Soil Microbial Community Taxonomic Differences and Biomarkers
Hierarchically clustered genus-level heat-map analysis was used to identify high-dimensional biomarker taxa with significantly different abundances among the four diet treatments, comparing the samples collected in July and October. The clustering of these groups into two groups ("JTD_0 and JTD_15 ' and "JTD_30 and JTD_45 ') revealed that the bacterial and fungal taxa showed similar patterns (Figure 3). Compared to other groups, Bacteroidetes sp. and SJA-15 sp. were more abundant in the JTD_45 group in both the July and October samples (Figure 3A,B). However, although A4b, SB-5, and Brevundimonas were less abundant in JTD_30 and JTD_45 in the July samples, they were dominant in OTD_30 and OTD_45 (Figure 3A,B). In terms of fungal genera, the relative abundance of Phialemoniopsis, Pseudaleuria, and Echria decreased with increased protein content in both July and October (Figure 4A,B), whereas that of Pichia, Schizothecium, and Pseudeurotium was higher in October than in July. In JTD_45, although Malassezia and Mortierella had lower   Frontiers in Microbiology | www.frontiersin.org relative abundance in July, they were dominant in October (Figure 4A,B).

DISCUSSION
As an important part of the paddy ecosystem, the soil microbial community enhances paddy nutrient cycling and soil health and protects the environment (Bruggen et al., 2015). However, it is highly sensitive to the effects of external agricultural chemicals such as fertilizers and pesticides (Itoh et al., 2003;Zhong and Cai, 2007;Tang et al., 2020). Rice-fish co-culture can reduce chemical application by improving land productivity and soil fertility; this is because the excreta of aquatic animals provide sufficient nutrients for the growth of soil microorganisms, which, in turn, alter the oxygen conditions of paddy surface soil (Xie et al., 2010;Wang et al., 2019). To examine these dynamics further, we compared the soil microbial community in a rice-crab co-culture paddy under different dietary protein treatments in July and October.
Our findings indicate that bacterial communities are more diverse than fungal communities. Similar results have been reported in rice-rotation systems (Ma et al., 2020;Maguire et al., 2020) and other rice co-culture systems (Li et al., 2018;Hou et al., 2021). Compared to fungi, bacteria may be more abundant in these systems because their range of metabolic and nutritional strategies enables them to adapt to more complex ecological niches and habitats. This explains why bacteria are more abundant than fungi in soils with complex environments. Although the rice-crab co-culture system can increase aquatic food supply and farmer income and help sustain rice production, the challenge of sustaining rice yield without increasing the use of fertilizer-N, pesticides, and other costly inputs is on the increase. Our results may provide referential experiences for several other countries (e.g., Egypt, India, Indonesia, Thailand, Vietnam, the Philippines, Bangladesh, Myanmar, and Malaysia) that are practicing rice co-culture.
Co-culture systems involving rice and aquatic animals cause nutrients such as carbon, nitrogen, and phosphorus to accumulate in the paddy soil, increasing its nutrient load (Wan et al., 2019;Bashir et al., 2020). However, in our study, phylum-level bacterial community composition did not change significantly with culture time, and soil bacterial diversity and evenness were not different between the July and October samples. Xu (2020) found that soil microbial biomass does not change significantly over years of co-culture (Xu, 2020). Our study site at Panjin city, which has over 30 years of co-culture history, had a higher comprehensive index of soil nutrients and decomposing bacteria (Proteobacteria and Chloroflexi) than other locations in China (Bashir et al., 2021), suggesting that Proteobacteria, Chloroflexi, and Bacteroidetes play important roles in structuring the ricecrab co-culture soil flora.
However, our beta diversity analysis revealed that the fungal community differed between the July and October samples. To our surprise, the soil mycobiota were dominated by relatively few saprotrophic fungal species instead of hundreds of bacterial taxa. Saprotrophic macrofungi such as Ascomycota and Basidiomycota are key regulators of soil nutrient cycling and decomposition. In the present study, Ascomycota was dominant in the control group, whereas Basidiomycota was the most abundant in the protein-diet treatments. Although soil macrofungi are generalists and occur in heterogeneous environments, their diversity and richness can be impacted by indiscriminate feeding during aquaculture (Crowther et al., 2012). Additionally, fungal alpha diversity indices (Chao1, Shannon, and Simpson) were higher in the JTD_30 and OTD_30 groups than in the low-protein treatment groups, indicating that saprotrophic macrofungal diversity was negatively related to dietary protein content. Although a high protein diet can effectively promote the growth performance and digestive physiology level of mitten crab, it can decrease soil pH from the rice-crab coculture system and negatively impact soil bacteria community diversity. Moreover, formulated feed does not make up a larger portion of the diet of E. sinensis, and excess formulated feed added to the rice-crab system sinks to the bottom of field, which may lead to water pollution (Guo et al., 2015). Therefore, a lower dietary protein level such as ca. 15% is appropriate for the rice-crab co-culture system both environmentally and economically.
In both July and October samples, the high-protein groups (JTD_30, JTD_45, OTD_30, and OTD_45) had greater relative abundances of Bacteroidetes and Thiobacillus than the control and low-protein groups (JTD_0, JTD_15, OTD_0, and OTD_15). We hypothesize that high-protein feeding might gradually reduce soil pH while promoting organic matter, nitrate nitrogen, phosphate, and other nutrient accumulation. This hypothesis is supported by previous findings that the relative abundance of Bacteroidetes and Thiobacillus increases during biogas irrigation (Yuge et al., 2014). The occurrence of these bacteria in the rice soil microbiota may reflect its physiological and biochemical conditions, such as its high organic matter, protein content, and pH. Therefore, soil pH may play a critical role in driving the topological structure of the co-occurrence network for soil microbiota, and the network was the most stable at neutral pH (Chen et al., 2021).
Brevundimonas, as a plant growth-promoting rhizobacterial genus, enhances the yield and growth of rice (Singh et al., 2016), potato (Naqqash et al., 2020), and wheat (Anuj et al., 2012), thereby allowing reduced application of chemical fertilizers and causing minimal environmental impact. We found that highprotein diets (JTD_30 and JTD_45) significantly reduced the relative abundance of Brevundimonas in both the July and October samples.

CONCLUSION
In conclusion, we evaluated how the soil microbial community of a rice-crab co-culture system responded to different dietary protein levels using partial 16S/18S rRNA gene analysis. A number of dominant and rare bacterial phyla (groups) were detected. Our results indicated that bacterial community structure was not affected by diet treatment, while in fungi, significant community structure change was evident. Moreover, high-protein diets reduced Chao 1 index and the Observed_species index in both July and October. Low-protein diet soils displayed greater bacterial diversity than high-protein diets soils, which reveal that a dietary protein level of ca. 15% is appropriate for this rice-crab co-culture system. Owing to the limited sampling and restricted sampling area (Panjin City, northern China), these conclusions cannot be generalized to other locations. However, our findings represent an important step toward understanding how diet alters soil microbial communities in rice-crab co-culture systems. In the future, we will conduct a more comprehensive study of the interaction between feeding patterns and soil microbial communities in rice-crab co-culture systems to formulate more accurate recommendations for this co-culture model.

DATA AVAILABILITY STATEMENT
The data presented in the study are deposited in the NCBI Sequence Read Archive (SRA) database and Biosample accession number was SAMN17393294.

AUTHOR CONTRIBUTIONS
YL and XL designed the experiments, analyzed the data, contributed reagents, materials, and analysis tools, and wrote the manuscript. YY and LL generated biological samples. YY, LL, and QH performed the experiments. YL and LL performed statistical analysis. All authors read and approved the final manuscript.

ACKNOWLEDGMENTS
We thank Panjin Guanghe Fisheries Co., Ltd. for providing the crabs.