Method and system for supervising recognition of soil nitrous oxide Birch effect and pulse emission calculation

Abstract

The invention relates to a method and a system for supervised identification of soil nitrous oxide Birch effect and pulse emission calculation, and the method comprises the steps: obtaining continuous N2O flux observation data and related covariables, and organizing the data to a time axis to form a time sequence; rolling window judgment and event judgment scoring are carried out on the time sequence, an effective pulse event segment is generated, and Birch pulse form constraint is applied; determining event start-stop boundaries in the morphologically constrained pulse event segment, and estimating a background baseline level of the flux for fine tuning the event start-stop boundaries, generating a time interval of the pulse from start to fall back to the baseline; unequal interval time integration is carried out on the flux exceeding the base line in the time interval, the pulse emission amount of the pulse event segment and the contribution rate to annual or seasonal total emission are obtained, uncertainty propagation is carried out through Monte Carlo sampling, and the confidence interval between the event level and the annual level is output.

Description

Technical Field

[0001] This invention relates to the field of ecological environment monitoring and greenhouse gas emission assessment technology, and in particular to a supervised method and system for identifying the Birch effect and pulse emission of soil nitrous oxide. Background Technology

[0002] Soil nitrous oxide (N2O) often exhibits short-duration high-flux "hot periods / pulses" (Birch effect) after triggering events such as rainfall, fertilization, irrigation, and freeze-thaw cycles. These events are relatively rare in time, have amplitudes higher than background fluxes, and contribute disproportionately to annual emissions. High-frequency automatic observations have revealed a significant increase in such events, but there is still no consensus on the quantitative definition and unified identification process for "hot periods." Low-frequency manual observations are also prone to missing peak values, thus affecting the accuracy of regional to global scale budget assessments and inventory compilation.

[0003] Existing thresholding methods (mean ± standard deviation, 1.5×IQR, etc.) and unsupervised anomaly detection (such as IsolationForest, IF) exhibit significant threshold fluctuations with background variance and seasonal segmentation under different seasons and data distribution conditions. SD and IQR are susceptible to skewed heavy-tailed distributions in threshold setting, leading to missed detections, while IF may cause over-labeling due to hyperparameters and distribution assumptions. In some scenarios, it even attributes almost all seasonal cumulative flux to "hot periods." Comprehensive evaluation also indicates that a single statistical method is difficult to apply stably across locations and management scenarios.

[0004] Therefore, there is an urgent need to design an automated identification and pulse emission calculation scheme that combines statistical robustness with event pattern determination capabilities. This would help reduce the false negatives / false positives of the traditional fixed threshold method and the over-labeling problem of a single unsupervised method, and provide a foundation for consistent identification and accounting across seasons and regions. Summary of the Invention

[0005] To address the short-term, high-flux "Birch effect" of soil nitrous oxide (N2O) triggered by rainfall, fertilization, irrigation, and freeze-thaw cycles, existing methods, centered on fixed thresholds or single unsupervised anomaly detection, suffer from insufficient robustness across seasons and management scenarios, high false negative and false positive rates, difficulty in automatically refining event start and end boundaries, and a lack of a standardized accounting chain for the quantification of uncertainty throughout the entire process. This invention aims to propose a supervised method and system for identifying the Birch effect and calculating pulse emissions of soil nitrous oxide. This system automatically and robustly identifies Birch pulses without altering the original flux sequence sampling grid, accurately defines event start and end boundaries, calculates peak value, duration, pulse emissions, and contribution rate using a unified standard, and outputs event-level and year-level confidence intervals.

[0006] To achieve the above objectives, the present invention provides the following solution: A supervised method for identifying the Birch effect and pulse emission calculation of soil nitrous oxide includes: Acquire continuous N2O flux observation data and related covariates, and organize them into a time series on the same time axis; The time series is subjected to rolling window judgment and event discrimination scoring to generate valid pulse event segments, and Birch pulse shape constraints are applied to exclude pulse event segments that do not conform to typical pulse characteristics. The start and end boundaries of the events in the morphologically constrained pulse event segment are determined, and the background baseline level of the flux is estimated to fine-tune the start and end boundaries of the events, generating the time interval from the start of the pulse to the fallback to the baseline. The flux exceeding the baseline is integraled over unequal time intervals over the time interval to obtain the pulse emission of the pulse event segment and its contribution to the total annual or seasonal emissions. Uncertainty propagation is then performed using Monte Carlo sampling to output the event-level and annual-level confidence intervals.

[0007] Optionally, generating the valid pulse event segment includes: Set a rolling window of target length, calculate the median and interquartile range of flux values ​​within each rolling window, and adaptively determine the threshold coefficient based on the target false alarm rate to generate a rolling threshold: T = A + k·IQR; where T is the rolling threshold, A is the median, IQR is the interquartile range, and k is the threshold coefficient; when the values ​​of multiple consecutive sampling points in the time series exceed the rolling threshold in the rolling window, they are marked as candidate impulse event segments; Multidimensional feature vectors characterizing the event amplitude, morphology, and temporal background are extracted from the candidate impulse event segments. The multidimensional feature vectors are then input into a supervised temporal event scoring model to obtain event scores. Candidate impulse event segments whose event scores exceed the screening threshold are considered valid impulse events. The supervised temporal event scoring model is obtained by training a gradient boosting decision tree algorithm or a one-dimensional temporal neural network model using a training set. During the training process, historical data from the training set is used for weakly supervised pre-training, and manually labeled pulse samples are used to fine-tune and correct the model. A focus loss function or a scoring threshold adjustment strategy is introduced to reduce the model's bias towards non-events.

[0008] Optionally, applying morphological constraints to the Birch pulse includes: The duration of the pulse event segment is limited to a target range, with the event flux curve exhibiting a single-peak shape as the target. If a secondary peak exists in the pulse event segment within the target range, the pulse event segment is divided into multiple single-peak sub-events according to the rules of local extrema and valley depth. The upward and downward slope characteristics and lag duration of the pulse event segment are constrained, i.e., the ratio of the maximum rise rate to the maximum fall rate is controlled within the target range, and the time lag between the pulse occurrence and the triggering event must not exceed the preset maximum duration.

[0009] Optionally, determining the start and end boundaries of the event includes: A parallel pruning precise linear time algorithm and a Bayesian online change point detection method are used to analyze the change points of the flux sequence in the morphologically constrained pulse event segment. If the difference between the start or end time of any pulse event segment given by the pruning precise linear time algorithm and the Bayesian online change point detection method does not exceed a preset tolerance, and the confidence level of the Bayesian online change point detection method at the corresponding time is not lower than the target threshold, then the time point is confirmed as the event start and end boundary of the pulse event segment. If the difference between the start or end time of any pulse event segment given by the pruning precise linear time algorithm and the Bayesian online change point detection method exceeds the time tolerance or the confidence level is insufficient, then local fine-tuning is triggered in the vicinity of the boundary, and secondary positioning is performed on the boundary positions with divergence to obtain the event start and end boundary of the pulse event segment.

[0010] Optionally, generating the time interval includes: The baseline of the soil N2O flux observation data is fitted using observation data during non-pulse periods. That is, the baseline is fitted after masking the currently identified pulses. Based on the baseline, the influence of unmasked residual pulses is reassessed. The mask is updated and fitted again until the baseline converges to a stable state. In the vicinity of the event's initiation and termination within the event's start and end boundaries, a local interpolation curve is constructed based on the original observation data points to determine the intersection point of the pulse segment flux curve and the baseline curve. This intersection point is then used as the time interval between the pulse's start and its return to the baseline.

[0011] Optionally, during the process of integrating fluxes exceeding the baseline over unequal intervals over the time interval, if there is a brief data gap within the event time window, a local linear state-space model is introduced to perform multiple interpolation on the data gap. That is, the data gap is predicted and interpolated using adjacent observations, and multiple interpolation sequences are generated by random sampling of the local linear state-space model to obtain multiple interpolation scenarios for the complete flux curve of the event. For each interpolation scenario, the pulse emission is calculated and multiple interpolation integration is performed.

[0012] Optionally, the output of the event-level and year-level confidence intervals includes: Joint Monte Carlo sampling of threshold coefficients, baselines, and event start and end boundaries is performed. Each simulation sampling yields a different set of threshold coefficients, baselines, and event start and end boundaries. The pulse emission and contribution rate are recalculated. The simulation is repeated a target number of times. All simulation results are summarized to obtain the probability distribution of event emission and contribution rate. Finally, the confidence range of emission and contribution rate for each pulse event is output in the form of confidence intervals.

[0013] Optionally, joint Monte Carlo sampling of the threshold coefficient, the baseline, and the event start and end boundaries includes: Based on the estimated error range of the threshold coefficient, sampling values ​​are taken, and several baseline curves are generated through the residual distribution of the baseline. These curves are used to recalculate the pulse flux for each simulation using a slightly deviated baseline. Based on the boundary posterior probability, the start and end times of each pulse event segment are sampled and perturbed in its neighborhood.

[0014] To achieve the above objectives, the present invention also provides a supervised identification system for the Birch effect and pulse emission calculation of soil nitrous oxide, comprising: The data management module is used to acquire continuous N2O flux observation data and related covariates, and organize them into a time series on the same time axis. The rolling window judgment and event scoring module is used to perform rolling window judgment and event discrimination scoring on the time series to generate valid pulse event segments; The shape constraint module is used to apply shape constraints to Birch pulses and to exclude pulse event segments that do not conform to typical pulse characteristics. The baseline estimation and integration module is used to determine the start and end boundaries of the event in the morphologically constrained pulse event segment and estimate the background baseline level of the flux. It is used to fine-tune the start and end boundaries of the event, generate the time interval from the start of the pulse to the fallback to the baseline, and perform unequal time integration on the flux exceeding the baseline in the time interval to obtain the pulse emission of the pulse event segment and its contribution rate to the annual or seasonal total emissions. The uncertainty assessment module is used to propagate uncertainty through Monte Carlo sampling and output confidence intervals at the event and year levels.

[0015] The beneficial effects of this invention are as follows: This invention combines adaptive thresholding with supervised scoring, taking into account both statistical robustness and event pattern discrimination, and significantly reduces missing peaks and false peaks.

[0016] This invention, based on the combination of PELT and BOCPD, enables the quantification and traceability of event boundaries, facilitating review and supervision.

[0017] This invention establishes a standardized integral link and uncertainty propagation mechanism, outputting event-level and year-level confidence intervals to improve the comparability and verifiability of accounting calculations.

[0018] The parameters of this invention are adaptive to the sampling frequency and context, making it suitable for multi-frequency data ranging from hourly to daily scales, and facilitating engineering and large-scale deployment. Attached Figure Description

[0019] To more clearly illustrate the technical solutions in the embodiments of the present invention or the prior art, the drawings used in the embodiments will be briefly introduced below. Obviously, the drawings described below are only some embodiments of the present invention. For those skilled in the art, other drawings can be obtained based on these drawings without creative effort.

[0020] Figure 1 This is a flowchart of a supervised method for identifying the Birch effect and pulse emission of nitrous oxide in soil according to an embodiment of the present invention. Figure 2 This is a schematic diagram of supervised event scoring and morphological constraints according to an embodiment of the present invention; Figure 3 This is a schematic diagram of N2O emission flux, rolling threshold, and pulse quantity for an hourly example of an embodiment of the present invention; Figure 4 This is a schematic diagram of N2O emission flux, rolling threshold, and pulse quantity as an example of daily scale in an embodiment of the present invention. Detailed Implementation

[0021] The technical solutions of the embodiments of the present invention will be clearly and completely described below with reference to the accompanying drawings. Obviously, the described embodiments are only some embodiments of the present invention, and not all embodiments. Based on the embodiments of the present invention, all other embodiments obtained by those skilled in the art without creative effort are within the scope of protection of the present invention.

[0022] To make the above-mentioned objects, features and advantages of the present invention more apparent and understandable, the present invention will be further described in detail below with reference to the accompanying drawings and specific embodiments.

[0023] like Figure 1As shown in the figure, this embodiment discloses a supervised method for identifying the Birch effect and calculating pulse emissions of soil nitrous oxide, including: acquiring continuous N2O flux observation data and related covariates, and organizing them uniformly onto the same time axis to form a time series; performing rolling window judgment and event discrimination scoring on the time series to generate effective pulse event segments, and applying morphological constraints on Birch pulses to exclude pulse event segments that do not conform to typical pulse characteristics; determining the start and end boundaries of events in the morphologically constrained pulse event segments, and estimating the background baseline level of flux to fine-tune the start and end boundaries of events, generating the time interval from the start of the pulse to its fall back to the baseline; performing unequal time interval integration on the flux exceeding the baseline within the time interval to obtain the pulse emission of the pulse event segment and its contribution rate to the annual or seasonal total emissions, and performing uncertainty propagation through Monte Carlo sampling to output event-level and annual-level confidence intervals.

[0024] Specifically, this embodiment discloses a supervised method for identifying the Birch effect and pulse emission calculation of soil nitrous oxide, including: (1) Detecting candidate segments using the rolling interquartile range adaptive method: Calculate the median and interquartile range (IQR) within a sliding window of length w. Set the dynamic threshold T = median + k·IQR, where k is determined by the target false alarm rate. α Adaptively determined, with a symmetric threshold set for negative flux; continuously not less than k min The out-of-threshold segments of each sampling point are recorded as candidate segments.

[0025] (2) Supervised time-series event scorer, such as Figure 2 As shown: Amplitude, morphology, and temporal context features (including but not limited to peak-to-baseline ratio, full width at half maximum, maximum up / down slope and its ratio, intra-segment interquartile range, skewness, kurtosis, unimodality statistics, time lag with the triggering event, and energy intensity in the frequency band corresponding to the typical duration of the pulse) are extracted from candidate segments. A gradient boosting model or a one-dimensional temporal neural network is used to output event scores s (between 0 and 1), and a cost-sensitive threshold is used for filtering. Weakly supervised pre-training and a small amount of manual fine-tuning are supported to alleviate the class imbalance problem.

[0026] (3) Morphological constraints: The candidate segments passed by the supervised time series event scorer are subject to constraints such as duration range, unimodality, slope ratio range, and maximum allowable lag with the triggering event; when there is a "shoulder peak" inside the candidate segment, it is split into multiple unimodal events according to the local extremum and valley depth rules to improve the physical interpretability and accounting consistency of the event.

[0027] (4) Joint verification of change points and boundary refinement: The start and end boundaries are jointly verified by combining the results of the "Pruned Exact Linear Time Algorithm" (PELT) and the "Bayesian Online Change Point Detection" (BOCPD). The former uses a piecewise linear mean model (including mean and slope terms) as the cost function for optimal segmentation, while the latter uses the Student's t-distribution likelihood to characterize the heavy-tailed residuals and calculate the posterior probability of the change point. When the boundary difference given by the two does not exceed the set tolerance and the posterior confidence is not lower than the threshold, the boundary is confirmed; otherwise, local fine-tuning is triggered in the boundary neighborhood.

[0028] (5) Baseline and Integration: Fitting a robust baseline B using “out-of-pulse” samples t (such as low-frequency terms in seasonal trend decomposition, or quantile regression) τ Using a value of 0.2–0.3, subsampling interpolation is performed on the event edges to find the intersection points where the baseline difference is zero; based on this, the event pulse emission is calculated. The peak value and duration are determined by the maximum value within the event segment and the time difference between the start and end of the event, and the contribution rate. .

[0029] (6) Uncertainty assessment: Monte Carlo sampling is performed on the threshold coefficient k, the baseline model residual and the change point boundary. The number of iterations can be 500-1000. The uncertainty is propagated to the event pulse emission and contribution rate, and the corresponding confidence interval and robustness criterion are output.

[0030] In summary, this invention combines statistical robustness with event pattern determination capability, significantly reducing the false negatives and false negatives of the Birch effect, improving boundary certainty and contribution rate estimation accuracy, and is suitable for greenhouse gas observation networks and emission inventory compilation.

[0031] More specifically, this embodiment discloses a supervised method for identifying the Birch effect and pulse emission calculation of soil nitrous oxide, including: Step 1: Collect continuous N2O flux observation data and related covariates (e.g., rainfall events, irrigation / fertilization times, soil moisture content, soil temperature, etc.), and organize all data onto the same time axis to form a time series. Explicitly mark missing data in the observations, and remove physically unreasonable outliers (e.g., negative concentrations or abnormal fluxes caused by sensor malfunctions). Then, establish time indexes and trigger event logs according to observation stations and seasons to record the occurrence times of triggering factors such as rainfall and fertilization. The time resolution of the data is preferably no less than 30 minutes (for high-frequency automatic monitoring), while the method is also compatible with low-frequency data input at the diurnal scale. The above preprocessing ensures the integrity and consistency of the input data.

[0032] Step 2: Detect possible pulse segments using the rolling interquartile range (IQR) adaptive thresholding method. The window length *w* is dynamically set according to the data sampling frequency (e.g., a 24–72 hour window can be used for hourly data, and a 7–21 day window for daily data; for sequences with limited sample size, the window can be appropriately shortened and an odd length is preferred to ensure statistical robustness). Within each sliding window of length *w*, the median *A* and interquartile range (IQR) of the flux values ​​are calculated, and the threshold coefficient *k* is adaptively determined based on the target false alarm rate *α*. The dynamic threshold is defined as *T* = *A* + *k* *IQR* (a negative threshold can also be symmetrically set for possible negative flux). The value of coefficient *k* is calibrated so that the proportion of background normal segments exceeding the threshold is approximately *α* (e.g., *α* is 1%–5%), thereby controlling the false alarm rate while ensuring detection sensitivity. When the continuous frequency is not less than *w*... When the value of a sampling point exceeds a threshold, the continuous interval is marked as a candidate impulse event segment. To reduce the impact of statistical instability at the boundaries, a "reflection fill" process is applied to both ends of the time series during the calculation of rolling statistics. This mirrors and expands the data at the edges of the sequence to ensure that reliable median and IQR estimates are obtained at the beginning and end of the sequence. This adaptive threshold detection strategy can automatically adjust the threshold based on the data noise level and seasonal variations, initially screening out suspicious impulse event segments.

[0033] Step 3: For each candidate segment selected in Step 2, a supervised time series event scoring model (SES) is used for multidimensional feature extraction and event discrimination scoring. First, a series of features characterizing the event's amplitude, morphology, and temporal background are extracted from the candidate segments, including but not limited to: the ratio of peak to baseline level, half-width at half-maximum (the time span corresponding to half the peak height), the maximum slope and ratio of the rising and falling segments, the interquartile range of flux values ​​within the segment (measuring fluctuation amplitude), skewness and kurtosis (reflecting the asymmetry and sharpness of flux distribution), unimodality statistics (determining whether the flux curve has only a single peak), the time lag between the flux segment and the triggering event (the delay of pulse occurrence relative to triggering factors such as rainfall), and the frequency band energy intensity of the flux signal at a typical pulse duration scale. These features constitute a multidimensional feature vector input to the SES model. Next, the trained supervised learning model outputs an event score s (ranging from 0 to 1, representing the credibility of a real pulse event) for each candidate segment. SES models can employ gradient boosting decision tree algorithms (such as LightGBM, XGBoost, etc.) or one-dimensional temporal neural network models (such as temporal convolutional networks TCN, Transformer, etc.). Gradient boosting tree models effectively utilize manually extracted features and have high interpretability, while temporal neural networks can automatically learn temporal patterns from the original sequence. Both types of models contribute to improving the accuracy of impulse event recognition. During model training, it is preferable to first perform weakly supervised pre-training using a large amount of historical data (e.g., using unlabeled samples containing impulses to improve the model's ability to perceive impulse patterns), and then use a small number of manually labeled impulse samples to fine-tune and correct the model, thus balancing the model's generalization ability and specificity. To address the class imbalance problem caused by the relatively scarce impulse event samples, a focal loss function or an adjusted scoring threshold strategy can be introduced during training to reduce the model's bias towards the majority class (non-events). Finally, by setting an appropriate event judgment threshold p... (The score of candidate segments can be screened based on ROC curve analysis, such as maximizing the Yoden index, or by weighing the costs of false positives and false negatives according to the cost matrix. Only segments with scores exceeding the threshold are retained as valid impulse events.)

[0034] Step 4: Further apply Birch pulse morphology constraints to the segments selected as valid events by SES to eliminate interference from segments that do not conform to typical pulse characteristics. Specifically, firstly, limit the duration d of the event to a reasonable range [dmin, dmax] (e.g., the lower limit should not be less than several hours, and the upper limit should not exceed several days) to eliminate noise peaks with excessively short durations or atypical processes with excessively long durations. Secondly, require the event flux curve to exhibit a single-peak shape, that is, there should be only one main peak in an event segment; if there are secondary peaks (commonly known as "shoulder peaks") within the candidate segment, then the segment is split into multiple single-peak sub-events according to the rules of local extrema and valley depth, thereby ensuring that each event has a single morphology and clear physical indicators. Thirdly, constrain the upward and downward slope characteristics of the event, so that the ratio ρ of the maximum rise rate to the maximum fall rate is controlled within the interval [ρmin, ρmax] (e.g., between 0.3 and 3), to avoid judging abnormally steep or overly flat abnormal curves as valid events. Finally, it is stipulated that the time lag between the pulse occurrence and the triggering event must not exceed a pre-set maximum duration Lmax (e.g., 72 hours) to ensure that the identified pulses are indeed associated with reasonable triggering factors. Through these morphological constraints, spurious peaks and composite peaks that do not conform to the typical morphology of the Birch effect can be filtered out, further improving the physical authenticity of the identified pulse events and their robustness to background interference.

[0035] Step 5: For candidate events that have undergone morphological constraints, this step uses dual-channel change point detection to finely calibrate their start and end boundaries. Two methods are employed in parallel to analyze change points in the flux sequence: the Pruned Exact Linear Time (PELT) algorithm and the Bayesian Online Change Point Detection (BOCPD). The PELT method uses a piecewise linear model containing intercept and trend terms as the cost function, searching for the optimal split point in the flux sequence that minimizes the overall cost function, thereby identifying times when the flux mean or trend changes significantly. The BOCPD method, within a Bayesian framework, performs online change point inference for the sequence, assuming that the observation noise follows a heavy-tailed Student's t-distribution (which is more tolerant of outliers than the commonly used Gaussian distribution), and calculates the posterior probability of each new time point as a change point in real time. The two detection methods have different focuses: PELT is sensitive to step or linear trend changes, while BOCPD provides a confidence assessment of each time point as a change point. Combining the results of both methods improves the reliability of boundary determination. When the start or end times of an event given by PELT and BOCPD differ by no more than a preset tolerance Δt, and the posterior probability of the BOCPD change point at the corresponding time is not lower than the threshold p0 (e.g., around 0.7), the time point can be confirmed as the boundary of a pulse event. If the two methods disagree on the boundary (exceeding the time tolerance or insufficient confidence), a local fine-tuning process is triggered in the vicinity of the boundary: by narrowing the analysis window, increasing the temporal resolution, or adjusting the detection parameters, the boundary location with discrepancies is precisely located again. By leveraging the above joint judgment and refinement mechanism, the advantages of both statistical and Bayesian change point detection are fully utilized, ultimately obtaining high-confidence event start time ts and end time te, significantly improving the accuracy and stability of pulse event boundary determination.

[0036] Step 6: After determining the initial event time boundaries, this step estimates the background baseline level of flux and further fine-tunes the event boundary position accordingly. Specifically, a robust baseline B(t) of soil N2O flux is fitted using observation data from all non-pulse periods, reflecting the background variation trend of flux in the absence of pulses. Baseline modeling can employ seasonal trend decomposition methods (such as LOESS smoothing of the original sequence to extract low-frequency trends and seasonal components) or fit the low quantile level of flux based on quantile regression (e.g., selecting quantiles of τ≈0.2–0.3) to represent the background flux. Since pulse events significantly increase flux values, identified pulse segments need to be masked during baseline estimation to remove the influence of these abnormally high values ​​from the fitted data. Baseline fitting can employ an iterative strategy of "mask-fit-update": initially, the baseline is fitted after masking the currently identified pulses. Then, based on the obtained baseline, the influence of any unmasked residual pulses is reassessed, the mask is updated, and the fit is repeated. This process is repeated until the baseline converges to a stable state (convergence criterion can be determined by the baseline change being less than a threshold ε). The resulting B(t) is a robust background flux sequence that varies over time. With the baseline established, to more accurately determine the start and end times of events, this step refines the boundaries of each pulse event's edge region through subsampling interpolation. Specifically, in the vicinity of the event's initiation and termination, a local interpolation curve is constructed based on the original observation data points (e.g., using linear or spline interpolation of adjacent points, or resampling with increased temporal resolution). Then, numerical methods (e.g., solving F(t) = B(t)) are used to precisely locate the intersection point between the pulse segment flux curve and the baseline B(t) curve, thereby accurately determining the start and end times of the pulse event. In other words, the closest time point at which the interpolated flux value F(t) satisfies F(t) = B(t) is used as the time [ts, te] for the pulse to strictly start and fall back to the baseline. This positioning strategy based on the zero-point baseline difference corrects the original coarse boundary truncated by sampling points, making the start and end times of the pulse event more precise and continuous. This facilitates accurate integration of the pulse flux and avoids missing contributions from edge regions.

[0037] Step 7: After baseline correction and boundary clarification, the flux integral for each pulse event needs to be calculated to obtain the pulse emission and its contribution to total emissions, while considering the uncertainty of the results. For a given event, the flux exceeding the baseline is integrated over the time interval [ts, te] to calculate the pulse emission Ai. Simultaneously, the peak value (i.e., the maximum flux value during the pulse event) and duration (i.e., the time interval of the pulse event) are determined by combining the flux curve of the pulse event. Through these calculations, the contribution of the pulse event to total emissions is further evaluated. Since the sampling interval of environmental monitoring data may be uneven, a trapezoidal rule is preferably used for numerical integration to accommodate data with unequal time steps. Then, Ai is compared with the total emissions during the same period to calculate the contribution rate C of the event to the annual or seasonal total emissions (e.g., the proportion of Ai to the total annual flux). If a brief data gap exists within the event time window, this step introduces a local linear state-space model to perform multiple imputation on the missing portion to avoid underestimating or overestimating emissions. This assumes that the flux follows a linear evolution on a small timescale and is superimposed with random perturbations, using adjacent observations to predict and imput the missing segment. Through random sampling of the state-space model, multiple possible imputation sequences can be generated (i.e., Monte Carlo simulation to fill in the missing data), thus obtaining various scenarios for the complete flux curve of the event. Pulse emissions are calculated for each imputation scenario, achieving multiple imputation integration. The results of these different scenarios reflect the impact of the uncertainty of the missing data on the integration calculation. During the uncertainty assessment process, the uncertainty generated by imputation is taken into consideration to ensure that the final output of the event emissions and contribution rate ranges reflects the range of variation caused by the data gap. The event time window is the duration of one N2O emission process, from start to finish.

[0038] Step 8: To fully quantify the uncertainties in pulse identification and emission calculation, this step uses Monte Carlo sampling to conduct multi-source uncertainty propagation analysis on key influencing factors. Specifically, random disturbances are introduced simultaneously from the following three aspects for simulation: (1) Threshold parameter uncertainty: For the threshold coefficient k (and the associated false alarm rate α) in Step 2, sampling values ​​are taken according to its estimated error range to reflect the threshold setting uncertainty in the candidate segment detection link; (2) Baseline estimation uncertainty: For the residual distribution of the robust baseline in Step 6, several possible baseline curves are generated by bootstrap or by directly sampling the fitted residuals. For each simulation, a slightly biased baseline B'(t) is used to recalculate the pulse flux; (3) Boundary determination uncertainty: Using the boundary posterior probability output by BOCPD in Step 5, sampling disturbances are made in the neighborhood of the start and end times of each event to simulate possible deviations in boundary positioning. By combining Monte Carlo sampling with the three factors mentioned above, each simulation yields a different set of thresholds, baselines, and boundaries, based on which the pulse emission Ai' and contribution rate C' are recalculated. Repeating this simulation M times (e.g., M=500~1000), and summing all simulation results, the probability distribution of event emissions and contribution rates can be approximated. Based on this distribution, the required uncertainty indicators can be extracted, such as calculating 95% confidence intervals (commonly using P2.5 and P97.5 quantiles) to quantify the emission uncertainty range at the event level. If necessary, multiple events can be aggregated to obtain an uncertainty estimate of the total annual contribution. Through Monte Carlo uncertainty propagation analysis, uncertainties from various stages, including threshold selection, baseline fitting bias, boundary judgment error, and missing data interpolation, are comprehensively integrated. Finally, the reliable range of emissions and contribution rates for each pulse event is output in the form of confidence intervals, improving the reliability and interpretability of the results.

[0039] Step 9: Following the above steps, the method of this invention summarizes and outputs the identified and calculated results, and provides an external service interface. The system can automatically generate an analysis report, including the total number of detected pulse events, the distribution statistics of the duration and peak amplitude of each event, the contribution rate of each event to total emissions and its confidence interval, and a detailed list of each event (including the start and end times of the event, the corresponding triggering factor type and time lag, the SES model judgment score, the event pulse emission Ai and its confidence interval, etc.). The above results can be used for manual review and can also be used as a reference for subsequent research or management decisions. In addition, this system supports batch data processing and standardized application programming interface (API) calls, which can seamlessly integrate the identification and calculation functions into the data pipeline of the greenhouse gas observation network or the emission inventory accounting platform to achieve large-scale deployment and external service. The entire method has good scalability and universality, providing a reliable tool for consistent Birch pulse identification and emission accounting across sites and regions.

[0040] In one embodiment, this embodiment employs hourly monitoring (sampling interval of 15 minutes), such as... Figure 3 As shown: I. Data and Threshold Parameters: Sampling interval is 15 minutes; rolling window w is 24 hours; target false alarm rate. α The coefficient k is adaptively obtained by setting it to 0.02; the lower limit of the number of consecutive points exceeding the threshold, kmin, is set to 2. The rolling threshold method of "median + k·interquartile range" is used to detect candidate segments, and reflection fill is used at both ends of the window to reduce edge deviation.

[0041] II. Supervised Scorer: Employs a gradient boosting tree model (using LightGBM as an example, 200 trees, maximum depth 6). Input features include peak-to-baseline ratio, full width at half maximum (FWHM), maximum up / down slope and its ratio, intra-segment interquartile range, skewness, kurtosis, unimodality statistics, time lag to the trigger event, and energy intensity in the frequency band corresponding to the typical pulse duration. Class imbalance is addressed through focus loss or threshold shifting. The decision threshold is p. Adjusted based on the "Receiver Operating Characteristic Curve - Youden Index".

[0042] III. Morphological Constraints and Boundary Verification: The lower limit of duration, dmin, is set to 1 hour; the slope ratio range... ρ Set the value to [0.4, 2.5]; the maximum allowable lag Lmax is set to 72 hours. Run PELT and BOCPD in parallel. If the boundary difference between the two does not exceed the set tolerance and the posterior of BOCPD is not lower than p0=0.7, then the boundary is confirmed; otherwise, local fine-tuning is performed in the neighborhood.

[0043] IV. Baseline, Integral, and Uncertainty: Fitting a robust baseline using "out-of-pulse" samples (low-frequency terms or quantile regression of seasonal trend decomposition) τ (≈0.25), perform subsampling interpolation on both ends of the event to solve for the intersection point where the difference between the baseline and the baseline is zero; use the trapezoidal method to complete the time integration to obtain the event pulse emission, peak value and duration; perform Monte Carlo sampling M=1000 times to integrate the threshold parameters, baseline residuals and boundary posterior uncertainties, and output the contribution rate and confidence interval.

[0044] Test results (in this example): ① Number of events detected: 3.

[0045] ②Event details: Event 1: Lasted 2.25 hours, peak concentration 1.76 kg N ha -1 h -1 Pulse emission: 0.102 kg N ha -1 ; Event 2: Lasted 9.25 hours, peak value 2.16 kg N ha -1 h -1 Pulse emission: 0.264 kg N ha -1 ; Event 3: Lasted 3.5 hours, peak concentration 2.04 kg N ha -1 h -1 Pulse emission: 0.139 kg N ha -1 .

[0046] ③ Contribution rate (C): 0.376; Total event pulse emissions: 0.505 kg N ha -1 Total emissions: 1.342 kg Nha -1 .

[0047] In one embodiment, this embodiment uses diurnal scale observations (sampling interval of 1 day), such as... Figure 4 As shown: I. Data and Threshold Parameters: The sampling interval is 1 day; the rolling window w is 14 days; the threshold parameter is based on the target false alarm rate. α Adaptive estimation. Considering water-dominated seasonal differences, the baseline is calculated using quantile regression (…). τ =0.25), and segmented fitting according to soil moisture status (drought / transitional / wet) to reduce the obfuscation of background by season and management practices.

[0048] II. Supervised Scoring and Labeling Strategy: A process of "weakly supervised pre-training + a small amount of fine-tuning with annotations" is adopted; 50 candidate segments are randomly selected monthly for manual review and correction to refresh the model's threshold and feature distribution. To suppress excessive merging, the upper limit of the allowed event duration is d. maxTake 7 days to cover the release of the "long tail".

[0049] III. Boundaries and Integration: PELT and BOCPD are used to jointly verify the boundaries; piecewise linear baselines and subsampling interpolation are used for long events to improve the accuracy of zero-intersection positioning; the event pulse emission and contribution rate are integrated using the trapezoidal method; for periods with fewer missing measurements, linear interpolation is performed using adjacent points, while for continuous missing measurement periods, multiple interpolation is performed using a local linear state-space model and propagated in uncertainty.

[0050] Test results (in this example): ① Number of events detected: 2.

[0051] ②Event details: Event 1: Lasted for 5 days, peak at 1.023 kg N ha -1 d -1 Pulse emission: 2.61 kg N ha -1 ; Event 2: Lasted for 9 days, peak at 1.279 kg N ha -1 d -1 Pulse emission: 4.52 kg N ha -1 .

[0052] ③ Contribution rate (C): 0.561; total event pulse emissions totaled 7.130 kg N ha-1, and total emissions totaled 12.702 kg N ha-1.

[0053] The following examples are for illustrative purposes only and do not limit the scope of protection of the invention in any way. Unless otherwise stated, parameters can be adaptively adjusted according to site data frequency, noise level, and regulatory / accounting requirements.

[0054] This embodiment also provides a supervised identification system for the Birch effect and pulse emission calculation of soil nitrous oxide, comprising: a data management module for acquiring continuous N2O flux observation data and related covariates, and organizing them uniformly onto the same time axis to form a time series; a rolling window judgment and event scoring module for performing rolling window judgment and event discrimination scoring on the time series to generate valid pulse event segments; a morphological constraint module for applying morphological constraints on Birch pulses to exclude pulse event segments that do not conform to typical pulse characteristics; a baseline estimation and integration module for determining the start and end boundaries of events in the morphologically constrained pulse event segments, estimating the background baseline level of flux, fine-tuning the start and end boundaries of events, generating the time interval from the start of the pulse to its fall back to the baseline, and performing unequal time interval integration on the flux exceeding the baseline within the time interval to obtain the pulse emission of the pulse event segment and its contribution rate to the annual or seasonal total emissions; and an uncertainty assessment module for performing uncertainty propagation through Monte Carlo sampling and outputting event-level and annual-level confidence intervals.

[0055] Specifically, this embodiment also provides a supervised identification system for the Birch effect and pulse emission calculation of soil nitrous oxide. The system includes: a data management module, a rolling interquartile range threshold module, a supervised time-series event scoring module, a morphological constraint module, a change point detection and boundary refinement module, a baseline estimation and integration module, an uncertainty assessment module, and a result visualization and interface module. Each module is decoupled and connected using a standardized data structure of "time series + event table," and supports data interface connection with greenhouse gas observation networks / inventory systems. The system includes: a data management module for collecting, verifying, and organizing multi-source data; a rolling threshold module for calculating the median and interquartile range within a rolling window of length w and adaptively determining the coefficient k based on the target false alarm rate to construct a threshold T = median + k·IQR; a supervised time-series event scoring module for extracting amplitude, morphology, and temporal context features from candidate segments and outputting event scores; a morphological constraint module for applying priors such as duration, slope ratio, and time lag to eliminate segments with inappropriate morphology; a change point detection and boundary refinement module for performing joint change point tests based on PELT and BOCPD and refining the start and end boundaries within tolerance; a baseline estimation and integration module for estimating robust baselines on out-of-pulse samples and performing unequal-interval time integration on the out-of-baseline portion to obtain peak value, duration, pulse emission, and contribution rate; an uncertainty assessment module for performing Monte Carlo propagation on the uncertainties of thresholds, baselines, and boundaries and outputting event-level / annual confidence intervals; and a result visualization and interface module for displaying results, generating reports, and connecting to external data / service interfaces.

[0056] The embodiments described above are merely preferred embodiments of the present invention and are not intended to limit the scope of the present invention. Various modifications and improvements made to the technical solutions of the present invention by those skilled in the art without departing from the spirit of the present invention should fall within the protection scope defined by the claims of the present invention.

Claims

1. A supervised method for identifying the Birch effect and pulse emission calculation of soil nitrous oxide, characterized in that, include: Acquire continuous N2O flux observation data and related covariates, and organize them into a time series on the same time axis; The time series is subjected to rolling window judgment and event discrimination scoring to generate valid pulse event segments, and Birch pulse shape constraints are applied to exclude pulse event segments that do not conform to typical pulse characteristics. The start and end boundaries of the events in the morphologically constrained pulse event segment are determined, and the background baseline level of the flux is estimated to fine-tune the start and end boundaries of the events, generating the time interval from the start of the pulse to the fallback to the baseline. The flux exceeding the baseline is integraled over unequal time intervals over the time interval to obtain the pulse emission of the pulse event segment and its contribution to the total annual or seasonal emissions. Uncertainty propagation is then performed using Monte Carlo sampling to output the event-level and annual-level confidence intervals.

2. The supervised identification method for soil nitrous oxide Birch effect and pulse emission calculation method according to claim 1, characterized in that, Generating the valid pulse event segment includes: Set a rolling window of target length, calculate the median and interquartile range of flux values ​​within each rolling window, and adaptively determine the threshold coefficient based on the target false alarm rate to generate a rolling threshold: T = A + k·IQR; where T is the rolling threshold, A is the median, IQR is the interquartile range, and k is the threshold coefficient; when the values ​​of multiple consecutive sampling points in the time series exceed the rolling threshold in the rolling window, they are marked as candidate impulse event segments; Multidimensional feature vectors characterizing the event amplitude, morphology, and temporal background are extracted from the candidate impulse event segments. The multidimensional feature vectors are then input into a supervised temporal event scoring model to obtain event scores. Candidate impulse event segments whose event scores exceed the screening threshold are considered valid impulse events. The supervised temporal event scoring model is obtained by training a gradient boosting decision tree algorithm or a one-dimensional temporal neural network model using a training set. During the training process, historical data from the training set is used for weakly supervised pre-training, and manually labeled pulse samples are used to fine-tune and correct the model. A focus loss function or a scoring threshold adjustment strategy is introduced to reduce the model's bias towards non-events.

3. The supervised identification method for soil nitrous oxide Birch effect and pulse emission calculation method according to claim 1, characterized in that, The morphological constraints applied to the Birch pulse include: The duration of the pulse event segment is limited to a target range, with the event flux curve exhibiting a single-peak shape as the target. If a secondary peak exists in the pulse event segment within the target range, the pulse event segment is divided into multiple single-peak sub-events according to the rules of local extrema and valley depth. The upward and downward slope characteristics and lag duration of the pulse event segment are constrained, i.e., the ratio of the maximum rise rate to the maximum fall rate is controlled within the target range, and the time lag between the pulse occurrence and the triggering event must not exceed the preset maximum duration.

4. The supervised identification method for soil nitrous oxide Birch effect and pulse emission calculation method according to claim 1, characterized in that, Determining the start and end boundaries of the event includes: A parallel pruning precise linear time algorithm and a Bayesian online change point detection method are used to analyze the change points of the flux sequence in the morphologically constrained pulse event segment. If the difference between the start or end time of any pulse event segment given by the pruning precise linear time algorithm and the Bayesian online change point detection method does not exceed a preset tolerance, and the confidence level of the Bayesian online change point detection method at the corresponding time is not lower than the target threshold, then the time point is confirmed as the event start and end boundary of the pulse event segment. If the difference between the start or end time of any pulse event segment given by the pruning precise linear time algorithm and the Bayesian online change point detection method exceeds the time tolerance or the confidence level is insufficient, then local fine-tuning is triggered in the vicinity of the boundary, and secondary positioning is performed on the boundary positions with divergence to obtain the event start and end boundary of the pulse event segment.

5. The supervised identification method for soil nitrous oxide Birch effect and pulse emission calculation method according to claim 1, characterized in that, Generating the time interval includes: The baseline of the soil N2O flux observation data is fitted using observation data during non-pulse periods. That is, the baseline is fitted after masking the currently identified pulses. Based on the baseline, the influence of unmasked residual pulses is reassessed. The mask is updated and fitted again until the baseline converges to a stable state. In the vicinity of the event's initiation and termination within the event's start and end boundaries, a local interpolation curve is constructed based on the original observation data points to determine the intersection point of the pulse segment flux curve and the baseline curve. This intersection point is then used as the time interval between the pulse's start and its return to the baseline.

6. The supervised identification method for soil nitrous oxide Birch effect and pulse emission calculation method according to claim 1, characterized in that, During the process of integrating fluxes exceeding the baseline at unequal intervals over the time interval, if there is a brief data gap within the event time window, a local linear state-space model is introduced to perform multiple interpolation on the data gap. That is, the data gap is predicted and interpolated using adjacent observations, and multiple interpolation sequences are generated by random sampling through the local linear state-space model to obtain multiple interpolation scenarios for the complete flux curve of the event. For each interpolation scenario, the pulse emission is calculated and multiple interpolation integration is performed.

7. The supervised identification method for soil nitrous oxide Birch effect and pulse emission calculation method according to claim 1, characterized in that, The output of the event-level and year-level confidence intervals includes: Joint Monte Carlo sampling of threshold coefficients, baselines, and event start and end boundaries is performed. Each simulation sampling yields a different set of threshold coefficients, baselines, and event start and end boundaries. The pulse emission and contribution rate are recalculated. The simulation is repeated a target number of times. All simulation results are summarized to obtain the probability distribution of event emission and contribution rate. Finally, the confidence range of emission and contribution rate for each pulse event is output in the form of confidence intervals.

8. The supervised identification method for soil nitrous oxide Birch effect and pulse emission calculation method according to claim 7, characterized in that, The joint Monte Carlo sampling of the threshold coefficient, the baseline, and the event start and end boundaries includes: Based on the estimated error range of the threshold coefficient, sampling values ​​are taken, and several baseline curves are generated through the residual distribution of the baseline. These curves are used to recalculate the pulse flux for each simulation using a slightly deviated baseline. Based on the boundary posterior probability, the start and end times of each pulse event segment are sampled and perturbed in its neighborhood.

9. A supervised identification system for soil nitrous oxide Birch effect and pulse emission calculation system implemented according to any one of claims 1-8, characterized in that, include: The data management module is used to acquire continuous N2O flux observation data and related covariates, and organize them into a time series on the same time axis. The rolling window judgment and event scoring module is used to perform rolling window judgment and event discrimination scoring on the time series to generate valid pulse event segments; The shape constraint module is used to apply shape constraints to Birch pulses and to exclude pulse event segments that do not conform to typical pulse characteristics. The baseline estimation and integration module is used to determine the start and end boundaries of the event in the morphologically constrained pulse event segment and estimate the background baseline level of the flux. It is used to fine-tune the start and end boundaries of the event, generate the time interval from the start of the pulse to the fallback to the baseline, and perform unequal time integration on the flux exceeding the baseline in the time interval to obtain the pulse emission of the pulse event segment and its contribution rate to the annual or seasonal total emissions. The uncertainty assessment module is used to propagate uncertainty through Monte Carlo sampling and output confidence intervals at the event and year levels.

Images