Burned area nitrogen cycle multi-factor correlation analysis method

Abstract

The invention discloses a burned area nitrogen cycle multi-factor correlation analysis method, and belongs to the technical field of burned area nitrogen cycle analys.The burned area nitrogen cycle multi-factor correlation analysis method comprises the steps that updated larch forests of different types of burned areas are selected as research objects, and forest stand characteristics, soil physicochemical properties and soil microbial community characteristics are investigated and analyzed; the influence of soil physicochemical properties on nitrogen cycle related functional microorganisms in a burned area recovery process is explored, a microbial mechanism of vegetation recovery influencing a soil nitrogen cycle key process is revealed, and then a theoretical basis is provided for optimal management of burned area forest ecosystem vegetation recovery.

Description

TECHNICAL FIELD

[0001] The application belongs to the technical field of fire-burned site nitrogen cycle analysis, and particularly relates to a fire-burned site nitrogen cycle multi-factor correlation analysis method. BACKGROUND

[0002] Nitrogen is one of the most important macro nutrients limiting plant growth, and is also the soil nutrient most significantly affected by forest fire. Nitrification and denitrification are key processes for regulating soil nitrogen loss and accumulation, and play a key role in regulating the dynamic changes of nitrogen accumulation and loss, and the occurrence of the two processes is executed by corresponding microorganisms. In the succession process of fire-burned site community, there is uncertainty between soil nitrification and denitrification and soil nitrogen accumulation and loss in view of the changes of soil resource availability and vegetation type, and it is urgent to analyze the microbe-driven mechanism of soil nitrification and denitrification and soil nitrogen, so as to provide a theoretical basis for the optimized management of vegetation restoration of fire-burned site forest ecosystem. SUMMARY

[0003] The purpose of the application is to provide a fire-burned site nitrogen cycle multi-factor correlation analysis method, so as to provide a theoretical basis for the optimized management of vegetation restoration of fire-burned site forest ecosystem.

[0004] The technical scheme adopted by the application is as follows:

[0005] A fire-burned site nitrogen cycle multi-factor correlation analysis method comprises the following steps:

[0006] (1) Select three different fire-burned times of light, moderate and severe fire-burned larch forests as research areas, and set up control plots in the nearby non-fire-burned areas, the slope of the plots is controlled to be 5-20 degrees, the altitude range is 800-950 m, the area of each plot is 20*20 m, and 3 repetitions are added to each plot under each treatment, and a total of 36 plots are set up;

[0007] (2) In each plot, 3 5*5 m small sample plots and 3 1*1 m small sample plots are set up by using PVC pipes, and the species, abundance and coverage of understory shrubs and herbs are investigated in the vegetation growth flourishing period, and the species and growth forms of all plants in the sample plots are identified;

[0008] Among them, the indexes for measuring the diversity of vegetation composition are Margalef richness index D, Shannon index H ′ , Simpson index H and Pielou evenness index J;

[0009] (3) In the vegetation growth flourishing period, soil profiles are dug in each plot, and soil is taken by using a cutting ring in the 0-15 cm, 10-20 cm and 20-40 cm soil layers, and the soil physical properties are tested by repeating sampling 3 times in each layer;

[0010] Use a soil drill to collect soil in layers along the median line of each plot. Take 10 soil cores from each plot, and mix the soil from the same layers together. Spread the soil onto clean paper and place it in a cool, ventilated place indoors to air-dry. Remove stones and undecomposed organic matter, grind the soil, and pass it through a 2mm sieve. Use the quartering method to select a 1kg soil sample from the mixed soil. After thorough mixing, take a portion of the sample and place it in a sealed bag with the sampling information noted. Grind the remaining sample that has passed the 2mm sieve until it passes the 0.149mm sieve, place it in a sealed bag with the sampling information noted. Use the samples that have passed the 2mm sieve and the 0.149mm sieve to test the soil chemical properties.

[0011] Among them, exchangeable calcium (Ca 2+ ) and exchangeable magnesium (Mg 2+ ) was determined by ammonium acetate exchange followed by ICP-OES analysis; microbial carbon and microbial nitrogen were determined by chloroform fumigation extraction;

[0012] Using an ethanol-burned shovel, the soil layer from 0 to 15 cm below the surface was dug. After removing visible roots, the soil was passed through a 2 mm sieve. Each sample was collected from three sampling points and mixed. A portion of the soil was placed in a 4°C cold storage as the first soil sample for isotope labeling culture experiments and microbial biomass determination; a portion of the soil was immediately wrapped in tin foil and placed in a liquid nitrogen tank and placed in a -80°C freezer as the second soil sample for DNA extraction and high-throughput sequencing analysis.

[0013] (4) Adoption 15 The soil nitrification rate was determined by the N stable isotope tracer method. Before incubation, 2 mol L -1 25 g of the first soil sample was leached with KCl solution, and then the sample was placed in a 25 ° C constant temperature shaker at 180 rpm for 1 hour, and the filtrate was filtered and the NH4 + -N and NO3 - -N content and 15 N abundance, take 40g of the first soil sample and adjust the soil moisture to 50% of the field water holding capacity, culture it in a constant temperature incubator at 25℃ for one week, and add water on the 2nd, 4th and 6th days to maintain the soil moisture content. After taking it out on the 7th day, use K 15 The sample was sprayed with NO3 solution to obtain 1.1 mg 15 The labeled culture sample of N per gram of soil was divided into 4 equal parts and placed in a 25℃ constant temperature and humidity incubator for culture. Two parts were taken out for parallel extraction and filtration after 15 minutes and 24 hours respectively. The NO3 - -N content and 15 Determination of nitrogen abundance, including NH4 + -N and NO3 --N content was determined using a flow analyzer. 15 N uses diffusion method to convert NO3 - -N is converted into N2O and then measured using a stable isotope mass spectrometer;

[0014] (5) Adoption 15 The soil denitrification rate was determined by the N stable isotope tracer method. First, 2 g of the first soil sample was placed in a 12 ml headspace culture bottle, and three parallel samples were set up. The bottles were sealed with butyl rubber stoppers and the lids were tightened. The entire operation was completed in an anaerobic glove box filled with nitrogen. Then, all the headspace culture bottles containing samples were vacuumed and flushed with high-purity helium until the pressure was balanced. After repeating the operation three times, the headspace culture bottles were placed in a 25 ° C constant temperature incubator for pre-incubation to eliminate the residual NO3 in the samples. - and oxygen, and on the 2nd, 4th, and 6th days, the bottle was removed and vacuumed and flushed with helium three times before continuing the culture. After the 7th day, the bottle was removed and the pre-culture was terminated. One of the samples was taken to determine the NO3 - The remaining two were used as parallel samples to add K 15 NO3 and use a vortexer to process it so that the marker and the sample are evenly and fully contacted, and then continue to culture in a constant temperature incubator under dark conditions. After 24 hours, add 200 μL of 7 mol L -1 The incubation was terminated by adding ZnCl2 solution and vortexing. Two headspace glass bottles filled with N2 were used as blanks. The contents of all headspace bottles were immediately measured by stable isotope mass spectrometry. 28 N2, 29 N2, 30 The charge-to-mass ratio of N2 is used, and the production of N2 is used to characterize the denitrification rate;

[0015] (6) Microbial genomic DNA from the second soil sample was extracted using a magnetic bead-based nucleic acid extraction kit, and the DNA quality was tested using a Qubit fluorescence quantifier and a Nanodrop micro-UV-visible spectrophotometer;

[0016] (7) The genomic DNA that passed the quality inspection was fragmented into 350 bp, then end-repaired, A-tailed, and Illumina sequencing adapters were added. The 300-400 bp DNA fragments were amplified and enriched by PCR. Finally, the PCR products were purified using the AMPure XP system, the sequencing library was detected using an Agilent 2100 bioanalyzer, and the library was quantified using real-time PCR. Sequencing was performed on a NovaSeq X Plus sequencer using a PE150 sequencing strategy. The raw data from the Illumina platform was filtered using FASTP, and the filtered clean reads were used for assembly analysis.

[0017] (8) Using MEGAHIT, assemble Clean reads of each sample by k-mer 21, 41, 61, 81, 101, 121, 141 respectively, retain contig sequences ≥500bp, count the assembly of each sample; use MetaGeneMark to predict genes of the assembled sequences of each sample, use CD-HIT to combine sequences with sequence similarity ≥95% and reads coverage greater than 90% into a cluster, select the longest sequence in each cluster as the representative sequence to form a sequence set; use Bowtie2 to re-align the reads to the sequence set and count the number of reads, use the software Pathoscope to re-allocate the reads to the "most likely source" gene, and finally filter out genes with less than 2 reads to obtain the final non-redundant gene set, based on the Pathoscope software, use the Bayesian framework to test the sequence and alignment quality of each read to re-allocate, and then calculate the gene abundance for each sample S;

[0018] (9) Use the non-redundant gene set to align the NCycDB database using the DIAMOND software, select the highest similarity protein sequence through alignment and finally obtain the corresponding protein function annotation; use Kaiju to align the NR library to perform species annotation at each taxonomic level; select gene, species, and function richness and gene, species, and function evenness diversity index to evaluate Alpha diversity;

[0019] Use the R language gmodels package for principal component analysis, use the R language Vegan package to calculate the Bray-curtis distance matrix between samples based on gene, species, and function abundance for multivariate statistical analysis, use the R language Vegan package for permutational multivariate analysis of variance, similarity analysis to analyze the statistical significance of species and function grouping; use Metastats for species difference analysis between two comparison groups;

[0020] Call the lmer package in R language to build direct or indirect influence relationships between variables and calculate the standard coefficients between variables, call the piecewiseSEM package in R language, use the lm() function to build a regression model, call the summary function to view the model results, delete indicators with high collinearity between variables to ensure that the p-value is greater than 0.05; use the screened indicators to build a structural equation model to analyze the relationship between soil physical and chemical properties, nitrification or denitrification microorganisms, and nitrification or denitrification rates under the influence of vegetation restoration.

[0021] Further, in the step (2), the Margalef richness index D, the Shannon index H ′, Simpson index H and Pielou uniformity index J are calculated as follows:

[0022]

[0023]

[0024] Among them, N i is the importance value of the i-th species, N is the sum of the importance values ​​of all species, and S is the number of species in each plot.

[0025] Furthermore, in step (4), the calculation formula for the total nitrification rate is as follows:

[0026]

[0027] Where GNR is the total nitrification rate (μg N g -1 soil d -1 ); C0 and C1 represent the NO3 measured after 15 min and 24 h of labeling culture, respectively - -N content (mg kg -1 ); H0 and H1 represent the NO3 measured after 15 min and 24 h of labeling culture, respectively - -N 15 The percentage of N atoms exceeds 100%. T1 and T0 represent 15 min and 24 h, respectively.

[0028] Furthermore, in step (5), the calculation formula for the denitrification rate is:

[0029]

[0030] Among them, f 30 for 30 Molar fraction of N2; 28 N2, 29 N2, 30 N2 is the measured charge-to-mass ratio (m / z);

[0031]

[0032] in, is the estimated value of the N2 mass in the headspace bottle; is the average density of N2, which is 1.25 gL -1 ; V is the N2 space at the top of the headspace glass bottle;

[0033]

[0034] Among them, P 30 is the production rate of N2; M soil is the sample weight; t is the incubation time;

[0035]

[0036] wherein, F N is the total 15 N abundance in the soil after adding the labeled solution; C1 is the NO3 - concentration in the added label, C0 is the NO3 - concentration in the soil sample; H1 is the N abundance in the added label, H0 is the natural abundance of N in the soil sample. 15 15

[0037]

[0038] wherein, DR is the N2 production rate, indicating the denitrification rate (μg N g -1 soil h -1 ).

[0039] Further, in the step (7), the filtering standard for filtering the raw data of the Illumina platform by using FASTP is: removing reads containing unknown nucleotides (N) ≥10%, removing reads containing bases with phred quality score ≤20 ≥50%, removing reads containing adapters, and removing reads with read length <50.

[0040] Further, in the step (8), the calculation formula of gene abundance is:

[0041]

[0042] wherein, TPM represents the relative abundance of gene i in sample S, Li represents the length of gene i, Ni represents the number of times that gene i is detected in sample S, and j represents the number of genes.

[0043] To sum up, due to the adoption of the above technical solutions, the beneficial effects of the present application are:

[0044] 1、In the present application, based on stable isotope labeling technology, the soil nitrogen transformation process rate is quantified through indoor culture experiment, including soil total nitrification rate and denitrification rate, and the influence of vegetation restoration on soil nitrogen cycle in fire-burned sites is evaluated. Combined with the restoration characteristics of different types of fire-burned sites, the variation characteristics of soil nitrification and denitrification rates in different types of fire-burned sites are clarified.

[0045] ​​2、In the application, by high-throughput sequencing method, the abundance, diversity and community composition of soil nitrification and denitrification microorganisms are analyzed, and the distribution characteristics of nitrification and denitrification microorganisms in different types of fire-burned land are clarified. Combined with soil properties, the key influencing factors of soil nitrification and denitrification microorganisms are determined, and the differences of soil nitrification and denitrification functional microorganisms in different types of fire-burned land are revealed.

[0046] 3、In the application, based on the data of soil total nitrification rate and denitrification rate in different types of fire-burned land vegetation restoration, combined with regression analysis and structural equation model, the internal relationship between nitrification and denitrification rate, microbial community characteristics and soil nitrogen is comprehensively analyzed, and the microbial mechanism of vegetation restoration regulating soil nitrification and denitrification rate under the background of different intensity and different recovery time of forest fire is revealed, which provides a theoretical basis for the optimization management of vegetation restoration in fire-burned land forest ecosystem. BRIEF DESCRIPTION OF DRAWINGS

[0047] In order to more clearly illustrate the technical solutions of the embodiments of the application, the following will briefly introduce the drawings needed to be used in the embodiments, and it should be understood that the following drawings only show some embodiments of the application, and therefore should not be regarded as a limitation to the scope, and for those skilled in the art, other related drawings can also be obtained without creative labor on the basis of these drawings, wherein:

[0048] Figure 1 It is an analysis flowchart of the application. DETAILED DESCRIPTION

[0049] In order to make the purpose, technical scheme and advantages of the embodiments of the application more clear, the technical scheme in the embodiments of the application will be described clearly and completely in combination with the drawings in the embodiments of the application, and obviously, the described embodiments are part of the embodiments of the application, rather than all the embodiments. The components of the embodiments of the application described and shown in the drawings here can be arranged and designed in various different configurations.

[0050] The application is realized by the following technical solutions:

[0051] According to the fire records of Inner Mongolia Genhe Forest Industry Co., Ltd., the test was carried out in Hulun Buir City in the northern section of Daxing'anling (50°30'-52°30'N, 120°30'-122°30'E). Hulun Buir City is located in the northern part of Inner Mongolia Autonomous Region, belonging to the humid monsoon climate of cold temperate zone, dry in spring, easy to cause forest fire caused by lightning. The forest type in the study area is cold temperate coniferous forest with Larix gmelinii as the dominant tree species, accompanied by tree species such as Betula platyphylla, Pinus sylvestris, Populus davidiana, and shrubs such as Rhododendron dauricum and Rhododendron tomentosum, with typicality and representativeness, and obvious regional advantage.

[0052] The light, moderate and severe fire-burned larch forests in 1987, 2003 and 2015 were selected as the research area, and the control sample plots were set up in the nearby unburned area. The forest fire intensity classification was determined combined with the characteristics of surface combustible material and soil layer, herbaceous low shrub layer, tree layer and canopy layer, as shown in Table 1. The slope of the sample plot was controlled between 5-20 degrees, and the altitude range was 800-950m, to ensure that the site conditions were basically similar. The area of each sample plot was 20x20m, and 3 repetitions were set for each treatment, a total of 36 sample plots.

[0053]

[0054]

[0055] Table 1 Investigation content and evaluation standard of comprehensive fire-burned index in Daxing'anling forest area

[0056] In each sample plot, 3 5x5m small sample plots and 3 1x1m small sample plots were set up with PVC pipes, and the species, richness, coverage of understory shrubs and herbs were investigated in July respectively. The species and growth form of all plants in the sample plot were identified through Chinese Plant Species Information System (http: / / www.iplant.cn).

[0057] The most commonly used index to measure the diversity of vegetation composition is Margalef richness index (D), Shannon index (H ′ ), Simpson index (H) and Pielou evenness index (J), which are calculated as follows:

[0058]

[0059] In the formula, N iImportance value for the ith species, N is the sum of importance values for all species, and S is the number of species within each plot.

[0060] Soil profiles were dug in the experimental plots in July, and soil was sampled with a ring knife at 0-15 cm, 10-20 cm, 20-40 cm, with 3 replicates per layer. Soil physical properties were tested. Soil was sampled with an auger along the midline of each plot, with 10 cores per plot. The same layers of soil were combined to form a mixed sample, which was spread onto a clean paper and air-dried in a cool, well-ventilated room. Stones and undecomposed organic matter were picked out, ground, and passed through a 2 mm sieve. The mixed soil was quartered to obtain a 1 kg sample, which was mixed thoroughly and then divided into two parts. One part was placed in a sealed bag with the sampling information noted. The remaining sample was further ground to pass through a 0.149 mm sieve and placed in a sealed bag with the sampling information noted. Soil chemical and physical properties were analyzed according to the relevant standards in the Forestry Industry Standards of the People’s Republic of China. Exchangeable calcium (Ca 2+ ) and exchangeable magnesium (Mg 2 + ) were analyzed by ICP-OES (Agilent, Santa Clara, CA, USA) after exchange with ammonium acetate; microbial carbon (MBC) and microbial nitrogen (MBN) were determined by the chloroform fumigation extraction method.

[0061] The soil physical properties include soil density, maximum moisture capacity, minimum moisture capacity, and capillary moisture capacity; the soil chemical properties include pH value, soil organic carbon, total nitrogen, total phosphorus, total potassium, available nitrogen, available phosphorus, available potassium, invertase activity, urease activity, phosphatase activity, and catalase activity.

[0062] An ethanol-burned shovel was used to dig the soil layer from 0 to 15 cm underground. After removing visible roots, the soil was passed through a 2 mm sieve. Each sample was collected from three sampling points and mixed. The sample size was 50 g. One portion was placed in a 4°C cold storage as the first soil sample for isotope labeling culture experiments and microbial biomass determination; another portion of the soil (20 g) was wrapped in tin foil and immediately placed in a liquid nitrogen tank back to the laboratory, placed in a -80°C refrigerator, and used as the second soil sample for DNA extraction and high-throughput sequencing analysis.

[0063] use 15 The soil nitrification rate was determined by the N stable isotope tracing method. Before the incubation, 2 mol L -1 25 g of the first soil sample was extracted with KCl solution (soil-liquid ratio 1:5), and then the sample was placed in a 25 ° C constant temperature shaker at 180 rpm for 1 hour, and the filtrate was filtered and the NH4 + -N and NO3 - -N content and 15 N abundance. Take 40g of the first soil sample and adjust the soil moisture to 50% of the field water holding capacity. Incubate in a 25℃ constant temperature incubator for one week and add water on the 2nd, 4th and 6th days to maintain the soil moisture content. Take it out on the 7th day and use K 15 NO3 solution ( 15 The abundance of N was 10%), and 1.1 mg 15 The labeled culture sample of each gram of soil was divided into 4 equal parts (10 g each) and placed in a 25°C constant temperature and humidity incubator for culture. Two parts were taken out for parallel extraction and filtration after 15 minutes and 24 hours respectively. The NO3 - -N content and 15 Determination of N abundance. NH4 + -N and NO3 - The -N content was determined using a flow analyzer (FIAstar5000, FOSS, Sweden). 15 N is converted into NO3 by diffusion method - -N was converted into N2O and then determined using a stable isotope mass spectrometer (ThermoFinnigan MAT, Bremen, Germany).

[0064] The calculation formula for the total nitrification rate is as follows:

[0065]

[0066] Where GNR is the total nitrification rate (μg N g -1 soil d -1 ); C0 and C1 represent the NO3 measured after 15 min and 24 h of labeling culture, respectively- -N content (mg kg -1 ); H0 and H1 represent the NO3 measured after 15 min and 24 h of labeling culture, respectively - -N 15 The percentage of N atoms exceeds 100%. T1 and T0 represent 15 min and 24 h, respectively.

[0067] Soil denitrification rate was also used 15 N stable isotope tracer method. 15 The N labeling method has been widely used to determine the loss pathways of soil N2, but all N2 fluxes in this experiment are potential rates and do not represent the actual situation of soil denitrification. First, take 2g of the first soil sample into a 12ml headspace glass bottle, and set up 3 parallel samples, seal it with a butyl rubber stopper and tighten the lid. The entire operation is completed in an anaerobic glove box filled with nitrogen. Then vacuum all sample bottles and flush them with high-purity helium (purity 99.99%) until equilibrium is reached at normal pressure. After repeating the operation 3 times, place the headspace culture bottle in a 25℃ constant temperature incubator for pre-culture to eliminate residual NO3 in the sample. - and oxygen. On the 2nd, 4th and 6th day, remove the bottle and repeat the vacuum and helium flushing 3 times before continuing the culture. After the 7th day, remove the bottle and end the pre-culture. Take one of the samples to determine the NO3 - The remaining two samples were added with K 15 NO3( 15 N abundance is 99.99%) and vortexed to ensure that the marker is evenly and fully exposed to the sample, and then placed in a constant temperature incubator under dark conditions for incubation. After 24 hours, 200 μL of 7 mol L -1 The incubation was terminated with ZnCl2 solution and vortexed, and two headspace glass bottles filled with N2 were used as blanks. The contents of all headspace bottles were immediately measured using a stable isotope mass spectrometer (Thermo Finnigan MAT, Bremen, Germany). 28 N2, 29 N2, 30 The charge-to-mass ratio of N2. The denitrification rate is characterized by the production of N2, which is calculated as follows:

[0068]

[0069] Where, f 30 for 30 Molar fraction of N2; 28 N2, 29 N2, 30 N2 is the measured mass-to-charge ratio (m / z).

[0070]

[0071] where, is the estimated mass of N2 in the headspace vial; is the average density of N2, which is 1.25 gL -1 ; V is the volume of the headspace of the glass vial.

[0072]

[0073] where, P 30 is the production rate of N2; M soil is the mass of the sample; t is the incubation time;

[0074]

[0075] where, F N is the total N abundance in the soil after adding the labeled solution; C1 is the concentration of NO3 15 in the added label; C0 is the concentration of NO3 - in the soil sample; H1 is the abundance of N in the added label; H0 is the natural abundance of N in the soil sample. - 15 15

[0076]

[0077] where, DR is the N2 production rate, which represents the denitrification rate (pg N g -1 soil h -1 ).

[0078] DNA extraction:

[0079] Microbial genomic DNA was extracted from the second soil sample using Magen Pure Soil DNA Kit or Magen Pure Stool DNA Kit (Guangzhou, China) according to the manufacturer’s instructions. DNA quality was checked using Qubit (Thermo Fisher Scientific, Waltham, MA) and Nanodrop (Thermo Fisher Scientific, Waltham, MA).

[0080] Illumina sequencing:

[0081] Genomic DNA that passed quality control was fragmented to 350 bp size, then end repaired, A-tailed, and adapters were added using MLtra TM ​​​DNA Library Prep Kit (NEB, USA) with Illumina sequencing adapters. DNA fragments of 300-400 bp were PCR amplified and enriched, and finally, the PCR products were purified using the AMPure XP system (Beckman Coulter, Brea, CA, USA), the sequencing library was detected using an Agilent 2100 Bioanalyzer (Agilent Technologies, Germany), and the library was quantified using real-time PCR. Sequencing was performed on a NovaSeq X Plus sequencer using a PE150 sequencing strategy.

[0082] The raw data of the Illumina platform was filtered using FASTP (version 0.20.0), and the filtering criteria were as follows: 1) remove reads containing unknown nucleotides (N) ≥10%; 2) remove reads with base phred quality score ≤20 ≥50%; 3) remove reads containing adapters; 4) remove reads with read length <50, and the filtered clean reads were used for assembly analysis.

[0083] Data analysis:

[0084] MEGAHIT (version 1.2.9) was used to assemble Clean reads of each sample according to k-mer 21, 41, 61, 81, 101, 121, 141, respectively, and sequences with contig ≥500 bp were retained. The assembly of each sample was statistically analyzed (including contig number, sequence length, N50, etc.). MetaGeneMark (version 3.38) was used to predict genes for each sample's assembled sequence, and the resulting gene sequences were merged into a cluster using CD-HIT (version 4.6) with sequence similarity ≥95% and read coverage greater than 90%. The longest sequence in each cluster was selected as the representative sequence to form a sequence set. Bowtie2 (version 2.3.5.1) was used to re-align the reads to the sequence set and count the number of reads. The reads were re-assigned to the "most likely source" gene using the software Pathoscope (v2.0.7), and genes with read support less than 2 were filtered out to obtain the final non-redundant gene set. Pathoscope software uses a Bayesian framework to test the sequence and alignment quality of each read to reassign it. Then, for each sample S, the gene abundance was calculated according to the following formula:

[0085]

[0086] TPM represents the relative abundance of gene i in sample S, Li represents the length of gene i, Ni represents the number of times gene i is detected in sample S (i.e. the number of reads on the alignment), and j represents the number of genes. After calculating the gene abundance results, the R language ggplot2 package is used to draw a violin plot of the number of genes in the comparison groups.

[0087] The functions and key mechanisms of the microbial community were mined using functional annotation with multiple databases. The non-redundant gene set was aligned with the NCycDB database using DIAMOND software (version 0.9.25, parameter e-value <0.00001), the most similar protein sequence was selected by alignment, and the corresponding protein functional annotation was finally obtained. The R language pheatmap package is used to draw a functional abundance heatmap, and the R language ggplot2 package is used to draw a functional abundance column chart.

[0088] Kaiju (version 1.10.1) is used to align the NR library for species annotation at each taxonomic level. Basic process: translate reads into amino acid sequences, align with Nr microbial library, and obtain species annotation of reads. Kaiju software: based on NR data sorting to construct sequence index of microorganisms (bacteria, archaea, viruses, fungi, small animals and plants), translate sequencing reads into amino acid sequences, break from the terminator place, then use Greedy mode, based on score screening target amino acid sequence, obtain species annotation of reads. The relative abundance of the species is the ratio of the number of reads of the species to the total number of sample reads.

[0089] Alpha diversity is evaluated based on gene / species / function richness (Chao, ACE) and gene / species / function evenness (Shannon, Simpon) diversity indices. Chao1 and ACE indices only consider the richness of species or, the larger the value, the higher the diversity. Shannon and Simpson comprehensively reflect the richness and evenness of species, the larger the value, the higher the balance of the sample. According to the formula of Alpha diversity index, it is calculated in R language. The stats package of R language is used to perform Welch's t test and Wilcox rank sum test between groups, Kruskal-Wallis rank sum test between three groups and more than three comparison groups, Turkey HSD multiple comparison after multiple group difference test, to further evaluate the significance of differences in different Alpha diversity indices between comparison groups.

[0090] Principal Components Analysis (PCA) was performed using R package gmodels, Bray-curtis distance matrix between samples was calculated based on gene / species / function abundance using R package Vegan for multivariate statistical analysis including principal co-ordinates analysis (PCoA), nonmetric multidimensional scaling (NMDS) and UPGMA (Unweighted Pair Group Method with Arithmetic Mean) clustering tree, and multivariate statistical analysis scatter plots were plotted using R package ggplot2. Permutation multivariate analysis of variance (Adonis) and Analysis of similarities (Anosim) were performed using R package Vegan to analyze the significance of species / function grouping, and boxplots were plotted using R package ggplot2.

[0091] Species difference analysis between two comparison groups was performed using Metastats. This analysis is a combination of different methods: first, T-test calculation, if the number of species in the group is less than the number of sample replicates, then calculate the P value based on Fisher's exact test; if the number of species in the group is greater than the number of sample replicates, and the number of replicates is greater than or equal to 8, then perform single-species Permutation test to calculate the P value; if the number of species in the group is greater than the number of sample replicates, and the number of replicates is less than 8, then mix the entire sample based on Permutation test to calculate the P value. Finally, multiple test correction is performed to calculate the q value. Boxplots were plotted using R package ggplot2.

[0092] Only for the presentation of significantly different species among the three groups, first, difference analysis among the three groups was performed to find the species significantly different among the three groups, then multiple comparisons among the three comparison groups were performed for the significantly different species to further determine whether the species had significant difference between the two comparison groups. Only when species S1 has the highest abundance in group A, and has significant difference between A-vs-B and A-vs-C, it is considered that species S1 is significantly enriched in group A. The species enrichment ternary plot shows the proportion of the abundance of the three biological groups, which can further view the changes or distribution of the three biological groups along some environmental gradients, and is expected to project the environmental information in the ternary plot. Species enrichment ternary plot was plotted using R package ggtern.

[0093] The lmer() package is called in R language to build direct or indirect influence relationship between each variable and to calculate standard coefficient between each relationship. The piecewiseSEM package is called in R language to build regression model by using lm() function, to view model result by using summary() function, to delete indexes with high collinearity between variables and to ensure that p value is greater than 0.05; the screened indexes are used to build structural equation model to analyze relationship between soil physical and chemical properties, nitrification or denitrification microorganism and nitrification or denitrification rate under influence of vegetation restoration.

[0094] The above is the embodiment of the present application. The foregoing is each preferred embodiment of the present application, and each preferred embodiment can be arbitrarily combined if not obviously contradictory or with a certain preferred embodiment as a prerequisite. The embodiment and the specific parameters in the embodiment are only for clearly describing the verification process of the application and are not used to limit the patent protection scope of the application. The patent protection scope of the application is still subject to the claims, and any equivalent structural changes made by using the content of the specification and the drawings should also be included in the protection scope of the application.

Claims

1. A multi-factor correlation analysis method for nitrogen cycle in burned areas, characterized in that: The following steps are involved: (1) Three lightly, moderately, and severely burned larch forests with different burn times were selected as the study area. Control plots were set up in nearby unburned areas. The slope of the plots was controlled between 5 and 20 degrees, and the altitude range was 800-950 m. Each plot was 20 × 20 m in size. Three replicates were added to each treatment, for a total of 36 plots. (2) In each plot, three 5 × 5 m and three 1 × 1 m plots were set up using PVC pipes. The species, abundance, and cover of understory shrubs and herbs were investigated during the period of lush vegetation growth. The species and growth forms of all plants in the plots were identified. Among them, the indicators for measuring vegetation composition diversity are Margalef richness index D and Shannon index H ′ , Simpson index H and Pielou evenness index J; (3) During the period of vigorous vegetation growth, soil profiles were dug from each plot. Soil was sampled at the 0-15 cm, 10-20 cm, and 20-40 cm layers using a circular cutter. Each layer was sampled three times to test soil physical properties. Use a soil drill to collect soil in layers along the median line of each plot. Take 10 soil cores from each plot, and mix the soil from the same layers together. Spread the soil onto clean paper and place it in a cool, ventilated place indoors to air-dry. Remove stones and undecomposed organic matter, grind the soil, and pass it through a 2mm sieve. Use the quartering method to select a 1kg soil sample from the mixed soil. After thorough mixing, take a portion of the sample and place it in a sealed bag with the sampling information noted. Grind the remaining sample that has passed the 2mm sieve until it passes the 0.149mm sieve, place it in a sealed bag with the sampling information noted. Use the samples that have passed the 2mm sieve and the 0.149mm sieve to test the soil chemical properties. Among them, exchangeable calcium (Ca 2+ ) and exchangeable magnesium (Mg 2+ ) was determined by ammonium acetate exchange followed by ICP-OES analysis; microbial carbon and microbial nitrogen were determined by chloroform fumigation extraction; Using an ethanol-burned shovel, the soil layer from 0 to 15 cm below the surface was dug. After removing visible roots, the soil was passed through a 2 mm sieve. Each sample was collected from three sampling points and mixed. A portion of the soil was placed in a 4°C cold storage as the first soil sample for isotope labeling culture experiments and microbial biomass determination; a portion of the soil was immediately wrapped in tin foil and placed in a liquid nitrogen tank and placed in a -80°C freezer as the second soil sample for DNA extraction and high-throughput sequencing analysis. (4) Adoption 15 The soil nitrification rate was determined by the N stable isotope tracer method. Before incubation, 2 mol L - 1 25 g of the first soil sample was leached with KCl solution, and then the sample was placed in a 25 ° C constant temperature shaker at 180 rpm for 1 hour, and the filtrate was filtered and the NH4 + -N and NO3 - -N content and 15 N abundance, take 40g of the first soil sample and adjust the soil moisture to 50% of the field water holding capacity, culture it in a constant temperature incubator at 25℃ for one week, and add water on the 2nd, 4th and 6th days to maintain the soil moisture content. After taking it out on the 7th day, use K 15 The sample was sprayed with NO3 solution to obtain 1.1 mg 15 The labeled culture sample of each gram of soil was divided into 4 equal parts and placed in a 25℃ constant temperature and humidity incubator for culture. Two parts were taken out for parallel extraction and filtration after 15 minutes and 24 hours respectively. The NO3 - -N content and 15 Determination of nitrogen abundance, including NH4 + -N and NO3 - -N content was determined using a flow analyzer. 15 N uses diffusion method to convert NO3 - -N is converted into N2O and then measured using a stable isotope mass spectrometer; (5) Adoption 15 The soil denitrification rate was determined by the N stable isotope tracer method. First, 2 g of the first soil sample was placed in a 12 ml headspace culture bottle, and three parallel samples were set up. The bottles were sealed with butyl rubber stoppers and the lids were tightened. The entire operation was completed in an anaerobic glove box filled with nitrogen. Then, all the headspace culture bottles containing samples were vacuumed and flushed with high-purity helium until the pressure was balanced. After repeating the operation three times, the headspace culture bottles were placed in a 25 ° C constant temperature incubator for pre-incubation to eliminate the residual NO3 in the samples. - and oxygen, and on the 2nd, 4th, and 6th days, the bottle was removed and vacuumed and flushed with helium three times before continuing the culture. After the 7th day, the bottle was removed and the pre-culture was terminated. One of the samples was taken to determine the NO3 - The remaining two were used as parallel samples to add K 15 NO3 and use a vortexer to process it so that the marker and the sample are evenly and fully contacted, and then continue to culture in a constant temperature incubator under dark conditions. After 24 hours, add 200 μL of 7 mol L -1 The incubation was terminated by adding ZnCl2 solution and vortexing. Two headspace glass bottles filled with N2 were used as blanks. The contents of all headspace bottles were immediately measured by stable isotope mass spectrometry. 28 N2, 29 N2, 30 The charge-to-mass ratio of N2 is used, and the production of N2 is used to characterize the denitrification rate; (6) Microbial genomic DNA from the second soil sample was extracted using a magnetic bead-based nucleic acid extraction kit, and the DNA quality was tested using a Qubit fluorescence quantifier and a Nanodrop micro-UV-visible spectrophotometer; (7) The genomic DNA that passed the quality inspection was fragmented into 350 bp, then end-repaired, A-tailed, and Illumina sequencing adapters were added. The 300-400 bp DNA fragments were amplified and enriched by PCR. Finally, the PCR products were purified using the AMPure XP system, the sequencing library was detected using an Agilent 2100 bioanalyzer, and the library was quantified using real-time PCR. Sequencing was performed on a NovaSeq X Plus sequencer using a PE150 sequencing strategy. The raw data from the Illumina platform was filtered using FASTP, and the filtered clean reads were used for assembly analysis. (8) MEGAHIT was used to assemble the clean reads of each sample according to k-mer 21, 41, 61, 81, 101, 121, and 141, and sequences with contig ≥ 500 bp were retained. The assembly status of each sample was statistically analyzed. MetaGeneMark was used to predict genes for the assembled sequences of each sample. The obtained gene sequences were merged into a cluster using CD-HIT with sequence similarity ≥ 95% and read coverage greater than 90%. The longest sequence in each cluster was selected as the representative sequence to form a sequence set. Bowtie2 was used to realign the reads to the sequence set and count the number of reads. The software Pathoscope was used to reallocate the reads to the genes of the "most likely source". Finally, genes with read support numbers less than 2 were filtered out to obtain the final non-redundant gene set. The sequence and alignment quality of each read were tested using the Bayesian framework based on the Pathoscope software to perform reallocation. Then, for each sample S, the gene abundance was calculated. (9) The non-redundant gene set was compared with the NCycDB database using DIAMOND software. The protein sequences with the highest similarity were selected through comparison and the corresponding protein function annotations were finally obtained. The species annotations at each classification level were performed by comparing the NR database using Kaiju. The alpha diversity was evaluated by selecting the diversity index based on gene, species, and function richness and gene, species, and function evenness. Principal component analysis was performed using the R language gmodels package. The R language Vegan package was used to calculate the Bray-Curtis distance matrix between samples based on gene, species, and function abundance for multivariate statistical analysis. The R language Vegan package was also used to perform permutational multivariate analysis of variance and similarity analysis to analyze the statistical significance of species and function groupings. Metastats was used to analyze species differences between the two comparison groups. The lmer package was called in R language to construct the direct or indirect influence relationship between each variable and calculate the standard coefficient between each relationship. The piecewiseSEM package was called in R language, and the lm() function was used to build a regression model. The summary function was called to view the model results. Indicators with high collinearity between variables were deleted to ensure that the p-value was greater than 0.

05. The screened indicators were used to construct a structural equation model to analyze the relationship between soil physical and chemical properties, nitrifying or denitrifying microorganisms, and nitrification or denitrification rates under the influence of vegetation restoration.

2. The multi-factor correlation analysis method for nitrogen cycle in burned areas according to claim 1, characterized in that: In the step (2), the Margalef richness index D and the Shannon index H ′ , Simpson index H and Pielou uniformity index J are calculated as follows: Among them, N i is the importance value of the i-th species, N is the sum of the importance values ​​of all species, and S is the number of species in each plot.

3. The multi-factor correlation analysis method for nitrogen cycle in burned areas according to claim 1, characterized in that: In the step (4), the calculation formula of the total nitrification rate is as follows: Where GNR is the total nitrification rate (μg N g -1 soil d -1 ); C0 and C1 represent the NO3 measured after 15 min and 24 h of labeling culture, respectively - -N content (mg kg -1 ); H0 and H1 represent the NO3 measured after 15 min and 24 h of labeling culture, respectively - -N 15 The percentage of N atoms exceeds 100%. T1 and T0 represent 15 min and 24 h, respectively.

4. The multi-factor correlation analysis method for nitrogen cycle in burned areas according to claim 1, characterized in that: In the step (5), the calculation formula for the denitrification rate is: Among them, f 30 for 30 Molar fraction of N2; 28 N2, 29 N2, 30 N2 is the measured charge-to-mass ratio (m / z); in, is the estimated value of the N2 mass in the headspace bottle; is the average density of N2, which is 1.25 gL -1 ; V is the N2 space at the top of the headspace glass bottle; Among them, P 30 is the production rate of N2; M soil is the sample weight; t is the incubation time; Among them, F N is the total amount of soil after adding the labeled solution 15 N abundance; C1 is NO3 in the added marker - concentration, C0 is the NO3 in the soil sample - Concentration; H1 is the concentration of added marker 15 N abundance, H0 is the abundance of N in the soil sample. 15 The natural abundance of N. Where DR is the N2 production rate, which represents the denitrification rate (μg N g -1 soil h -1 ).

5. The multi-factor correlation analysis method for nitrogen cycle in burned areas according to claim 1, characterized in that: In step (7), the filtering criteria for filtering the raw data of the Illumina platform using FASTP are as follows: removing reads containing unknown nucleotides (N) ≥ 10%, removing reads with bases ≥ 50% with a phred quality score ≤ 20, removing reads containing adapters, and removing reads with a read length < 50.

6. The multi-factor correlation analysis method for nitrogen cycle in burned areas according to claim 1, characterized in that: In step (8), the calculation formula for gene abundance is: Where TPM represents the relative abundance of gene i in sample S, Li represents the length of gene i, Ni represents the number of times gene i is detected in sample S, and j represents the number of genes.

Images