A multi-factor coupling evaluation method for ecosystem service function of wine grape planting area

Abstract

The application discloses a kind of multi-factor coupling evaluation methods of wine grape planting area ecosystem service function, belong to the technical field of agricultural ecological environment evaluation.Ground monitoring, remote sensing, terrain, soil, weather and cultivation management data are obtained, and grapeyard management unit is divided according to phenophase spatial-temporal alignment;For each management unit, counterfactual reference matched with environmental background is established, and the net effect of ecological system services such as supply, regulation and habitat maintenance is calculated;First-order contribution and second-order non-additive coupling contribution are determined by single-factor replacement and double-factor joint replacement;According to the model error, data missing rate and scale consistency correction weight, and introducing ecological pressure constraint, the multi-factor coupling evaluation value is obtained, the limiting service and main factor combination are identified, and the ecological conservation, maintenance promotion, regulation improvement and key repair partition are output.The application can reduce the influence of environmental background difference on the evaluation result, and improve the interpretability of the evaluation result.

Description

Technical Field

[0001] This invention relates to the field of agricultural ecological environment monitoring and ecosystem service evaluation technology, specifically to a multi-factor coupled evaluation method for the ecosystem service function of wine grape growing areas. Background Technology

[0002] In addition to supplying grape products, wine grape growing areas also influence ecological processes such as soil conservation, windbreak and sand fixation, net carbon sequestration, and habitat maintenance through vegetation cover, root stabilization, windbreaks, grass cover, and soil organic matter accumulation. Especially in arid, semi-arid, and desertified regions, wine vineyards are typically formed from wasteland, Gobi desert, or low-coverage land, and their ecological effects are simultaneously influenced by factors such as precipitation, irrigation, groundwater, soil salinity, topography, planting age, mulching methods, and water and fertilizer management. Therefore, when assessing the ecosystem services of grape-growing areas, it is necessary to distinguish between changes caused by grape cultivation activities and differences caused by natural background conditions within a unified spatial unit and phenological timescale.

[0003] In existing ecological assessments of wine grape growing areas, common practices include standardizing and weighting remote sensing indices, ecosystem service model outputs, or landscape pattern indicators. Some methods also directly compare service function values ​​across different plots or years. When soil, topography, and climate conditions vary significantly across plots, these direct comparisons can easily introduce environmental background differences into the ecological effects of grape cultivation activities. Furthermore, factors such as irrigation, cover, fertilization, and windbreak configuration often collectively influence the same ecosystem service function. Judging driving forces solely based on correlation or individual weights makes it difficult to distinguish between the effects of individual factors and the non-additive effects of their combined influence. When there are differences in the measured support levels and model calculation accuracy for different service functions, directly using fixed weights for aggregation may amplify the impact of low-reliability indicators on the overall results. Therefore, existing methods still have room for improvement in environmental background correction, identification of factor combination effects, and the conversion of assessment results into specific plot management measures. Summary of the Invention

[0004] The purpose of this invention is to provide a multi-factor coupled assessment method for ecosystem service functions in wine grape growing areas, comprising the following steps:

[0005] S1. Obtain ground monitoring data, remote sensing data, topographic data, soil data, meteorological data, grape growth data, and cultivation management data of the wine grape planting area to be evaluated. Divide the data into multiple phenological time windows according to the phenological period of wine grapes, and perform spatial registration and temporal alignment.

[0006] S2. Based on the boundaries of vineyard plots, slope aspect change lines, irrigation and drainage zone boundaries, windbreaks, and roads and ditches, the wine grape planting area to be evaluated is divided into multiple vineyard management units. A candidate reference record set is established for each vineyard management unit. The candidate reference records include at least one of the following: data records from the same spatial location during historical periods before vineyard establishment and data records from unplanted areas within the corresponding phenological time window. For the candidate reference records corresponding to the unplanted areas, a matching distance is determined based on the background factors of topography, soil, and climate, as well as the spatial distance between the unplanted area and the vineyard management unit, wherein the spatial distance is a component of the matching distance. A single counterfactual reference record is determined based on the matching distance, soil type, landform type, climate zone category constraints, and the results of the background factor balance test. When there is no single candidate reference record that simultaneously satisfies the category constraints, the matching distance threshold, and the preset balance threshold, a synthetic counterfactual reference value is generated based on multiple candidate reference records.

[0007] The matching and balance checks of the candidate reference records are performed in the following manner:

[0008] ;

[0009] In the formula, For the first The first vineyard management unit and the first The matching distance between the candidate reference records; The number of continuous background factors participating in the matching; The background factor number; and The first The first vineyard management unit and the first The first candidate reference record corresponding to the Background factor; For the first Interquartile range of the background factor; In order to be with the first The background factors have the same dimensions and are used to prevent positive numbers with a denominator of zero; Spatial distance; As a spatial distance reference value and ; and These are the background factor matching weights and spatial distance matching weights, respectively.

[0010] The and All are non-negative numbers and satisfy:

[0011] ;

[0012] For the first Matching weights of the background factors, The matching weights are spatial distances, and both are non-negative numbers;

[0013] The relative difference between a vineyard management unit and its corresponding counterfactual reference on the same background factor is defined as standardized difference and determined according to the following formula:

[0014] ;

[0015] In the formula, For the first The vineyard management unit in the first Standardized differences on the background factors; For the first The reference for the balance test of the vineyard management unit is the first... The background factor value; when the reference to be tested for balance is a single candidate reference record. For the single candidate reference record corresponding to the first The background factor value; when the reference to be tested for equilibrium is a synthetic reference. The first candidate reference record is obtained by weighting multiple candidate reference records according to their corresponding synthesis weights. Background factor value; and The meaning is the same as the definition in the aforementioned matching distance calculation;

[0016] When the standardized differences of each continuous background factor are not greater than the corresponding preset balance threshold, the reference to be tested for balance is determined to have passed the continuous background factor balance test.

[0017] Among the candidate reference records that simultaneously meet the constraints of soil type, landform type and climate zone category, matching distance threshold and preset balance threshold, the candidate reference record with the smallest matching distance is selected as the single counterfactual reference record.

[0018] S3. Determine the ecosystem service function mass of the vineyard management unit and its corresponding counterfactual reference service function, wherein the counterfactual reference service function mass is the ecosystem service function mass corresponding to a single counterfactual reference record, or the ecosystem service function mass calculated by multiple candidate reference records according to the composite weight; form the service function vector of the vineyard management unit and the counterfactual reference service function vector, and obtain the net ecosystem service effect vector after reducing the impact of environmental background differences based on the mass difference between the two.

[0019] S4. Establish a set of influencing factors consisting of natural background factors, ecological state factors, and cultivation management factors. Replace each individual influencing factor in the vineyard management unit with its corresponding counterfactual reference influencing factor value, and calculate the first-order contribution of each influencing factor to different ecosystem service functions. Then, replace both influencing factors with their corresponding counterfactual reference influencing factor values ​​simultaneously, and calculate the second-order non-additive coupling contribution of the two influencing factors. The counterfactual reference influencing factor value is derived from a single counterfactual reference record or generated from multiple candidate reference records according to a composite weight.

[0020] S5. Convert the net effect of the ecosystem services into a non-negative evaluation value, and determine the basic weight of each ecosystem service function based on the non-negative evaluation values ​​of each vineyard management unit within the same phenological time window.

[0021] The reliability coefficients of each ecosystem service function are determined according to the following formula:

[0022] ;

[0023] In the formula, For the first Within the first phenological time window Credibility coefficient of ecosystem service function; This corresponds to the missing data rate; The root mean square error between the model-calculated values ​​and the ground-measured values; The standard deviation of the measured values ​​on the ground; The scale consistency coefficient is the result of summarizing the calculation results at different spatial resolutions into the same vineyard management unit.

[0024] In order to be with the first Within the first phenological time window The quantities of ecosystem service functions have the same dimensions and are used to prevent positive numbers with a denominator of zero; the data missing rate satisfy The scale consistency coefficient satisfy When the first Ecosystem service functions in the first When no valid measured samples are available within a phenological time window, the reliability coefficient is... Set to 0;

[0025] The basic weights are then corrected using the credibility coefficient.

[0026] S6. The first-order contribution and the second-order non-additive coupling contribution of the model are dimensionless. The factor coupling correction of each ecosystem service function is determined based on the dimensionless first-order contribution and the second-order non-additive coupling contribution. The ecological pressure constraint is determined based on the measured values ​​of the ecological pressure factors and the corresponding pressure thresholds and upper limits. The net effect vector of ecosystem services, the factor coupling correction, the corrected weights, and the ecological pressure constraint are then fused to obtain the multi-factor coupling evaluation value of each vineyard management unit in different phenological time windows.

[0027] S7. Based on the corrected weights, determine the comprehensive first-order contribution intensity of each influencing factor and the comprehensive coupling contribution intensity of each combination of influencing factors. Then, determine the controlling influencing factor and the combination of controlling factors based on the comprehensive first-order contribution intensity and the comprehensive coupling contribution intensity. Determine the allowable adjustment set of cultivation management factors based on grape cultivation standards, the control range of management equipment, and historical effective management records. Perform influencing factor replacement simulations within the allowable adjustment set. Determine the ecological recoverability index based on the increment of the multi-factor coupling evaluation value before and after adjustment. Determine the evaluation zoning threshold and the recovery zoning threshold based on the multi-factor coupling evaluation value and the ecological recoverability index of each vineyard management unit within the same phenological time window. Output the spatial zoning results based on the high-low relationship between the multi-factor coupling evaluation value and the ecological recoverability index relative to the evaluation zoning threshold and the recovery zoning threshold.

[0028] Preferably, the phenological time window includes at least the dormancy period, the budding to flowering period, the young fruit development period, the color-changing period, and the ripening and harvesting period;

[0029] The remote sensing data includes one or more of the following: vegetation cover, leaf area index, net primary productivity, surface temperature, evapotranspiration, and remote sensing ecological indices.

[0030] The cultivation management data includes one or more of the following: irrigation amount, fertilizer amount, pesticide application amount, inter-row coverage ratio, soil tillage method, pruning intensity, windbreak belt parameters, and orchard establishment year.

[0031] Preferably, when there is no single candidate reference record that simultaneously satisfies the soil type, landform type, and climate zone category constraints, the matching distance threshold, and the preset balance threshold, multiple candidate reference records are selected from the candidate reference records that satisfy the soil type, landform type, and climate zone category constraints in ascending order of matching distance, and a synthetic counterfactual reference value is generated according to the following formula:

[0032]

[0033]

[0034] In the formula, To avoid matching dimensionless positive numbers with a distance of zero; For candidate reference record index, For the first A set of candidate reference records corresponding to each vineyard management unit; For the first The composite weight of each candidate reference record;

[0035] For the first The vineyard management unit in the first Within the first phenological time window The synthetic counterfactual reference quality of ecosystem service functions; For the first The candidate reference record in the first Within the first phenological time window The quality of ecosystem service functions; Ecosystem service function sequence number. This refers to the sequence number of the phenological time window;

[0036] For continuous impact factors, the composite counterfactual reference impact factor value is generated according to the following formula:

[0037]

[0038] In the formula, For the first The vineyard management unit in the first Within the first phenological time window A composite counterfactual reference value for each influencing factor; For the first The candidate reference record in the first Within the first phenological time window The values ​​of the impact factors, The impact factor number;

[0039] For categorical impact factors, the category state corresponding to the candidate reference record with the largest composite weight is used as the composite counterfactual reference state.

[0040] When there is no management activity corresponding to the cultivation management factor to be replaced in the candidate reference record, the counterfactual reference state corresponding to the cultivation management factor is determined as the preset no management baseline state.

[0041] For continuous background factors of synthetic counterfactual references, the corresponding synthetic background factor values ​​are obtained by weighted summation according to the synthetic weights of each candidate reference record; for categorical background factors of soil type, landform type and climate zone, the category state corresponding to the candidate reference record that satisfies the category constraints is adopted.

[0042] The generated synthetic counterfactual references are subjected to a balance test. When the standardized difference of any continuous background factor exceeds the corresponding preset balance threshold, a new candidate reference record is selected and the synthetic counterfactual reference value is recalculated. When the standardized difference of each continuous background factor is not greater than the corresponding preset balance threshold, the corresponding synthetic counterfactual reference value is determined as the counterfactual reference of the vineyard management unit.

[0043] Preferably, the ecosystem service function vector includes at least product supply function, water regulation function, soil conservation function, windbreak and sand fixation function, net carbon sequestration function and habitat maintenance function;

[0044] Among them, the soil conservation function is determined based on the difference between potential soil erosion and actual soil erosion; the windbreak and sand fixation function is determined based on the difference between potential wind erosion and current wind erosion under bare land conditions; and the net carbon sequestration function is determined jointly based on changes in vegetation carbon storage, changes in soil organic carbon, and carbon emissions from cultivation and management.

[0045] Preferably, the difference in ecosystem service functions is first determined according to the following formula:

[0046]

[0047] The difference in ecosystem service functions is then robustly standardized using the following formula:

[0048]

[0049] In the formula, For the first The vineyard management unit in the first The difference in quantity of the i-th ecosystem service function relative to the counterfactual reference within a phenological time window; Quality of service functions for vineyard management units; The corresponding counterfactual reference service function quality is either the service function quality corresponding to a single counterfactual reference record or the service function quality synthesized from multiple candidate reference records. and Each of the vineyard management units in the 1st The 95th and 5th quantiles of the difference in the i-th ecosystem service function within a phenological time window; In order to be with the first Within the first phenological time window The quantities of each ecosystem service function have the same dimensions and are used to prevent positive numbers with a denominator of zero;

[0050] in, The truncation function has the following expression:

[0051] .

[0052] Preferred, the first The impact factor on the first The vineyard management unit in the first Within the first phenological time window The first-order contribution of an ecosystem service function model is determined according to the following formula:

[0053]

[0054] No. The first influence factor and the first The second-order non-additive coupling contribution of each influencing factor is determined according to the following formula:

[0055]

[0056] In the formula, For the first Within the first phenological time window Computational models for ecosystem service functions; For the first The vineyard management unit in the first Vector of influencing factors within a phenological time window; , The index of each impact factor in the set of impact factors; Indicates the first The impact factor vector is formed by replacing each impact factor with its corresponding counterfactual reference value; Indicates the first The impact factor vector is formed by replacing each impact factor with its corresponding counterfactual reference value; Indicates that at the same time the first The first impact factor and the first The impact factor vector is formed by replacing each impact factor with its corresponding counterfactual reference value;

[0057] when When, the corresponding combination of influencing factors is determined as a synergistic enhancement combination; when When this happens, the corresponding combination of influencing factors is determined to be a mutually restrictive combination;

[0058] To eliminate the influence of dimensional differences between different ecosystem service functions on the comparison results of contribution intensity, the first-order contribution and the second-order non-additive coupling contribution of the model were dimensionless according to the quantile of the difference in the corresponding ecosystem service functions:

[0059]

[0060]

[0061] In the formula, This represents the first-order contribution of the dimensionless model. This represents the contribution of the second-order non-additive coupling to the dimensionless model. and The first Within the first phenological time window The 95th and 5th percentiles of the differences in ecosystem service functions; In order to be with the first Within the first phenological time window The quantities of each ecosystem service function have the same dimensions and are used to prevent positive numbers with a denominator of zero.

[0062] Preferably, the net ecosystem service effect corresponding to each vineyard management unit is considered. The values ​​were converted to non-negative evaluation values, and the basic weights of each ecosystem service function were determined using the entropy weight method based on the non-negative evaluation values ​​of each vineyard management unit within the same phenological time window. ;

[0063] The corrected weights are determined according to the following formula:

[0064]

[0065] In the formula, For the first Within the first phenological time window Adjusted weights for ecosystem service functions; These serve as the basic weights for corresponding ecosystem service functions; Summation index for ecosystem service functions ; The number of ecosystem service functions;

[0066] when

[0067]

[0068] At that time, the first If no ecosystem service function with valid credibility exists within a phenological time window, the multi-factor coupling assessment value will not be calculated for that phenological time window.

[0069] Preferably, the first step is to calculate the following formula: The vineyard management unit in the first Factor coupling correction for the i-th ecosystem service function within a phenological time window:

[0070]

[0071] in, , Let be the index of each different impact factor in the set of impact factors, and , A dimensionless positive number used to prevent the denominator from being zero during the calculation of factor coupling corrections;

[0072] Then calculate the multi-factor coupling evaluation value according to the following formula:

[0073]

[0074] In the formula, For the first The vineyard management unit in the first Multi-factor coupled evaluation values ​​within a phenological time window; The coefficient contributes to the coupling adjustment. This represents the ecological pressure constraint coefficient. This is a constraint on ecological pressure.

[0075] Preferably, in adjusting the weights In the ecosystem service functions, the first The vineyard management unit in the first Weight adjustment within each phenological time window Net effect of ecosystem services The ecosystem service function corresponding to the minimum value of the product is determined as the limited ecosystem service function of the vineyard management unit within the phenological time window;

[0076] Calculate the first according to the following formula. The combined first-order contribution strength of the influencing factors:

[0077]

[0078] In the formula, For the first The vineyard management unit in the first Within the first phenological time window The combined first-order contribution intensity of each influencing factor; The number of ecosystem service functions; For the first Within the first phenological time window Adjusted weights for ecosystem service functions; This represents the first-order contribution of the dimensionless model.

[0079] Calculate the first according to the following formula. The first impact factor and the first The combined coupling contribution strength of the influencing factors:

[0080]

[0081] In the formula, For the first The vineyard management unit in the first Within the first phenological time window The first impact factor and the first The combined coupling contribution strength of each influencing factor; This represents the contribution of the second-order non-additive coupling to the dimensionless model.

[0082] The overall first-order contribution intensity will be considered. The first corresponding to the maximum value One influencing factor was identified as the controlling influencing factor, and the combined contribution intensity was determined. The first corresponding to the maximum value The first impact factor and the first One influencing factor is determined to be the combination of the main controlling factors;

[0083] Within the allowable adjustment set of irrigation amount, fertilizer amount, inter-row cover ratio, soil tillage method, and windbreak belt parameters, an influencing factor replacement simulation was conducted, and the ecological restoreability index of the vineyard management unit was calculated according to the following formula:

[0084]

[0085] In the formula, For the first The vineyard management unit in the first Ecological recoverability indicators within a phenological time window; The set of permissible adjustments is defined by grape cultivation standards, the scope of control of management equipment, and historically valid management records. To allow adjustment of the vector of influence factors after replacement within the set; The multi-factor coupling evaluation value is calculated using the replaced influence factor vector; The values ​​are the multi-factor coupling assessment values ​​before adjustment;

[0086] In calculating the At that time, only the selected cultivation management factors within the allowed adjustment set are replaced with the adjusted values, while the unselected influencing factors and the corresponding ecological pressure constraints of the vineyard management unit within that phenological time window are not considered. The ecological pressure constraint is not recalculated with each single influencing factor replacement simulation, keeping the original values ​​unchanged.

[0087] Based on the multi-factor coupled assessment values ​​and ecological restoreability indicators of each vineyard management unit within the same phenological time window, the [number]th [period] was determined. The evaluation zone threshold corresponding to each phenological time window and restorative partition threshold ;

[0088] when and At that time, the corresponding vineyard management unit will be divided into an ecological conservation area;

[0089] when and At that time, the corresponding vineyard management unit will be divided into a maintenance and improvement zone;

[0090] when and At that time, the corresponding vineyard management unit will be divided into a regulation and improvement zone;

[0091] when and At that time, the corresponding vineyard management unit will be divided into key restoration areas.

[0092] Preferably, the vineyard management unit is divided into independent calibration units and verification units according to a preset ratio. The parameters of the corresponding ecosystem service function calculation model are corrected by using the grape yield and quality test data, soil erosion monitoring data, wind and sand transport monitoring data, soil moisture content and irrigation and drainage records, vegetation and soil carbon storage, cultivation management energy consumption and input usage data and biological survey data of the calibration unit.

[0093] A matching equilibrium test is performed on the topographic, soil, and climate background factors between the vineyard management unit and a single counterfactual reference record or a synthetic counterfactual reference value. For continuous background factors, if the standardized variance of any continuous background factor exceeds a preset equilibrium threshold, the single counterfactual reference record is re-screened or the synthetic counterfactual reference value is recalculated. For categorical background factors, the test is performed according to the corresponding category constraints. When the standardized variance of each continuous background factor is not greater than the equilibrium threshold and each categorical background factor meets the corresponding category constraints, the model error of each ecosystem service function is calculated using the validation unit, and the corresponding confidence coefficient is determined accordingly. The validation unit is only used for model error and evaluation result stability testing, does not participate in the correction of ecosystem service function calculation model parameters, and keeps the counterfactual reference relationship that has passed the equilibrium test and the ecosystem service function calculation model parameters determined by the calibration unit unchanged.

[0094] Compared with existing technologies, this invention provides a multi-factor coupled assessment method for ecosystem service functions in wine grape growing areas, which has the following beneficial effects:

[0095] 1. This invention utilizes counterfactual references with similar environmental backgrounds to calculate the net effect of vineyard ecosystem services. This can reduce the interference of differences in topography, soil, climate, and regional location on the judgment results of the ecological effect of grape cultivation, and avoid directly classifying plots with better natural conditions as plots with higher levels of ecological management.

[0096] 2. This invention calculates the first-order contribution and the second-order non-additive coupling contribution of the model by single-factor substitution and dual-factor joint substitution. It can compare the degree of influence of different influencing factors on the calculation results of ecosystem service functions, and identify the synergistic enhancement or mutual constraint relationship of different combinations of influencing factors in the model output, so that the analysis of influencing factors no longer depends solely on the degree of correlation or single weight.

[0097] 3. This invention organizes monitoring data according to phenological time windows such as dormancy period, budding to flowering period, young fruit development period, color change period and ripening and harvesting period, which can identify the differences in the effects of the same ecological factor at different growth stages and reduce the problem of the annual average value masking the stage-specific ecological pressure.

[0098] 4. This invention modifies the basic weights of ecosystem service functions based on model error, data missing rate, and scale consistency, so that the weights of ecosystem service functions with insufficient measured support or low scale consistency in the comprehensive evaluation are constrained by the credibility coefficient, thereby improving the stability of multi-model fusion results.

[0099] 5. Based on the results of ecosystem service evaluation, this invention further identifies limiting ecosystem service functions and combinations of controlling factors, and combines ecological restoreability indicators to perform spatial zoning, which can provide a basis for irrigation control, inter-row cover, soil improvement, windbreak belt optimization, and old garden renovation corresponding to specific plots. Attached Figure Description

[0100] Figure 1 A diagram illustrating the division of vineyard management units and counterfactual reference matching;

[0101] Figure 2 A schematic diagram of the comprehensive evaluation process for single-factor replacement, two-factor joint replacement, and multi-factor coupling.

[0102] Figure 3 A comparison chart of the average standardized difference of background factors and equilibrium threshold before and after counterfactual reference matching;

[0103] Figure 4 Summary chart of independent verification results for each ecosystem service function;

[0104] Figure 5 A graph showing the calculation results of single-factor and two-factor coupling contributions for management unit M028;

[0105] Figure 6A diagram showing the composition of the calibration set and independent validation set for the vineyard management unit. Detailed Implementation

[0106] refer to Figures 1-3 Example 1: Multi-source data acquisition and phenological time alignment

[0107] The experimental area for this embodiment is located in the grape-producing region on the eastern foothills of the Helan Mountains in Ningxia. The evaluation area consisted of continuously distributed vineyards and their surrounding wasteland, Gobi Desert, and grasslands. The collected data included:

[0108] Ground monitoring data include soil moisture content, soil temperature, electrical conductivity, pH value, organic matter, total nitrogen, total phosphorus, wind speed, air temperature and humidity, precipitation, groundwater depth and irrigation amount;

[0109] Remote sensing data, including multispectral satellite imagery, UAV multispectral imagery, vegetation cover, leaf area index, land surface temperature, evapotranspiration, net primary productivity, and remote sensing ecological indices;

[0110] Grape growth data, including bud break date, flowering date, veraison date, ripening date, yield, soluble solids content, total acid and phenolic content;

[0111] Management data includes irrigation volume, fertilizer application volume, pesticide application volume, pruning intensity, inter-row grass or mulch ratio, amount of branches returned to the field, burial method for winter protection, and windbreak structure.

[0112] In a specific set of experiments, continuous monitoring data from the Helan Mountain East Foothills grape-producing region in Ningxia from 2016 to 2023 were used, and a total of 64 vineyard management units were established. These units were stratified according to establishment year, soil type, slope grade, and irrigation method. Forty-four vineyard management units (M001–M044) were used as calibration units, and twenty vineyard management units (M045–M064) were used as validation units. These units remained fixed during the monitoring period. The calibration units were used to determine the parameters of the ecosystem service function calculation model, as well as the weights and adjustment parameters requiring calibration. The validation units were spatially independent of the calibration units and did not participate in the determination of these parameters; they were only used to verify model calculation errors and the stability of evaluation results.

[0113] The configuration of the calibration unit and the independent verification unit is shown in the following figure. Figure 6 Of these, there are 44 calibration units, accounting for 68.75% of the total number of vineyard management units; and 20 independent verification units, accounting for 31.25% of the total number of vineyard management units.

[0114] The dormancy period, budding to flowering period, young fruit development period, color-changing period, and ripening and harvesting period are divided according to the actual measured phenological dates of each year. For hourly sensor data, the mean, cumulative value, extreme value, or coefficient of variation are used to summarize the data within the phenological time window, based on the nature of the corresponding monitoring indicators. For satellite or UAV data, data that is close to the phenological nodes and meets the image quality requirements are used. For detection data that can reflect changes in phenological stages, the data are assigned to the corresponding phenological time window according to the sampling date. For data with slow annual changes, such as the number of years of orchard establishment, soil organic carbon, and habitat pattern, the most recent valid observation results of the current year are used, and no artificial interpolation is made based on the division of phenological time windows.

[0115] Synchronous ground calibration was performed on remotely sensed surface temperature data. During the 2022 growing season, 12 sets of remote sensing observation data were acquired, with 20 ground sampling points set up in each set, resulting in 240 sets of synchronous observation data. Ground temperature was measured using an infrared thermometer on the day of satellite transit, with the time difference between ground sampling and satellite transit controlled within 30 minutes. Based on all 240 sets of synchronous observation data, the root mean square error between remotely sensed surface temperature and measured ground temperature was calculated to be 0.65℃, and the coefficient of determination was... It is 0.90.

[0116] When different data have inconsistent spatial resolutions, separate data layers are retained at the vineyard row scale, plot scale, and landscape scale. Different scales are primarily used to organize corresponding types of input data and do not imply that the same ecosystem service function is calculated only at a single scale. For the same ecosystem service function participating in the scale consistency test, calculations are performed at at least two available spatial resolutions, and the results are then aggregated into the same vineyard management unit for comparison and evaluation. Specifically, the vineyard row scale is mainly used for organizing vegetation cover, soil moisture, and management input data; the plot scale is mainly used for organizing soil conservation, product supply, and carbon storage data; and the landscape scale is mainly used for organizing habitat connectivity and landscape fragmentation data.

[0117] Example 2: Vineyard Management Unit Division and Counterfactual Reference Establishment

[0118] First, extract the boundaries of the vineyard plots. For larger plots, subdivide them according to the following boundaries:

[0119] 1. Abrupt change in slope aspect exceeding the preset value;

[0120] 2. Irrigation and drainage boundaries belonging to different irrigation branch pipes or valve control ranges;

[0121] 3. The boundaries of the influence area on both sides of the windbreak forest belt;

[0122] 4. Main roads, drainage ditches, and bare ground isolation strips exceeding the preset width;

[0123] 5. A dividing line where soil salinity or groundwater depth changes abruptly.

[0124] For plots of land with significant differences in grape varieties, years of establishment, or main management methods, further subdivisions should be made according to the corresponding difference boundaries to ensure that grape varieties, years of establishment, and main management methods remain consistent within the same vineyard management unit.

[0125] Counterfactual references prioritize valid historical data records from the same spatial location prior to park establishment. When using historical data records, meteorological factors are re-statistically analyzed according to corresponding phenological time windows, and historical records with precipitation, accumulated temperature, and extreme weather conditions exceeding the matching threshold for the period to be evaluated are removed. When historical records that meet the continuity and matching requirements are lacking, spatial candidate reference records are selected from surrounding uncultivated wasteland, Gobi desert, or natural grassland areas according to the aforementioned matching distance, category constraints, and background factor balance test rules. The background factors, ecosystem service function data, and influencing factor data of the corresponding area constitute the candidate reference records.

[0126] The baseline factors used for matching include elevation, slope, aspect, soil texture, initial soil organic matter, soil salinity, groundwater depth, annual precipitation, accumulated temperature ≥10℃, and distance from mountains or rivers.

[0127] When the minimum matching distance still exceeds the preset threshold, instead of using a single counterfactual reference record, multiple candidate reference records are selected, and the synthesis weight is determined according to the reciprocal of the matching distance to generate a synthetic counterfactual reference value.

[0128] Specifically, no. The first vineyard management unit and the first The matching distance between candidate reference records is determined by the following formula:

[0129]

[0130] In the formula, For the first The first vineyard management unit and the first The matching distance between the candidate reference records; The number of continuous background factors participating in the matching; The background factor number; and The first The first vineyard management unit and the first The first candidate reference record corresponding to the Background factor; For the first Interquartile range of the background factor; In order to be with the first The background factors have the same dimensions and are used to prevent positive numbers with a denominator of zero; For the first The first vineyard management unit and the first The spatial distance between the corresponding spatial locations of the candidate reference records; when the candidate reference records originate from the same spatial location during a historical period before the park's construction. ; As a spatial distance reference value and ; and These are the background factor matching weight and the spatial distance matching weight, respectively.

[0131] The matching weights of the background factor and the spatial distance matching weights satisfy the following:

[0132]

[0133] in, For the first Matching weights of the background factors, represents the matching weight for spatial distance, and both are non-negative numbers.

[0134] When performing balance tests on the matching results of continuous background factors, the standardized differences between the vineyard management unit and the reference to be tested on each background factor are calculated:

[0135]

[0136] In the formula, For the first The vineyard management unit in the first Standardized differences on the background factors; For the first The reference for the balance test of the vineyard management unit is the first... The baseline factor value; when the reference for the balance test is a single candidate reference record. For the single candidate reference record corresponding to the first Background factor values; when the reference for the equilibrium test is a synthetic reference. The first candidate reference record is obtained by weighting multiple candidate reference records according to their corresponding synthesis weights. Background factor value; and The meaning is the same as the definition in the aforementioned matching distance calculation.

[0137] Among candidate reference records that simultaneously meet the constraints of soil type, landform type, and climate zone category, matching distance threshold, and preset balance threshold, the candidate reference record with the smallest matching distance is determined as the single counterfactual reference record.

[0138] When no single candidate reference record simultaneously satisfies the category constraint, the matching distance threshold, and the preset balance threshold, multiple candidate reference records are selected from those satisfying the category constraint in ascending order of matching distance, and the corresponding synthetic counterfactual reference service function material quality is determined according to the following formula:

[0139]

[0140] in:

[0141]

[0142] In the formula, To avoid matching dimensionless positive numbers with a distance of zero; Index for candidate reference records; For the first A set of candidate reference records corresponding to each vineyard management unit; For the first The composite weight of each candidate reference record; For the first The vineyard management unit in the first Within the first phenological time window The synthetic counterfactual reference quality of ecosystem service functions; For the first The candidate reference record in the first Within the first phenological time window The quality of ecosystem service functions.

[0143] For continuous impact factors, their composite counterfactual reference impact factor values ​​are determined according to the following formula:

[0144]

[0145] In the formula, For the first The vineyard management unit in the first Within the first phenological time window A composite counterfactual reference value for each influencing factor; For the first The candidate reference record in the first Within the first phenological time window The values ​​of the impact factors, This is the impact factor number.

[0146] For categorical influencing factors, the category state corresponding to the candidate reference record with the largest composite weight is used as the composite counterfactual reference state; when there is no management activity corresponding to the cultivation management factor to be replaced in the candidate reference record, the counterfactual reference state of the corresponding cultivation management factor is determined as the preset no-management baseline state.

[0147] For continuous background factors of synthetic counterfactual references, the corresponding synthetic background factor values ​​are obtained by weighted summation according to the synthetic weights of each candidate reference record; for categorical background factors of soil type, landform type and climate zone, the category status corresponding to the candidate reference record that satisfies the category constraint is adopted.

[0148] To test the counterfactual reference matching results, 20 validation units were used as subjects to compare the standardized differences before and after background factor matching. Each validation unit was tested on the [number missing]th [number missing]. The average standardized difference over the background factors is determined by the following formula:

[0149]

[0150] In the formula, For the set of verification units, To verify the number of units, this embodiment... As a control, the nearest unplanted area was used as a reference before matching, without screening for background factors such as topography, soil, and climate; after matching, a counterfactual reference was determined according to the aforementioned matching distance, category constraints, and balance test. In this embodiment, the preset balance threshold was set to 0.2, and the results are as follows: Figure 3 As shown.

[0151] Figure 3 To reflect the overall matching effect of the counterfactual references before and after matching for the 20 validation units, the balance test for each vineyard management unit was still determined item by item according to the aforementioned standardized differences and preset balance thresholds. After matching, the average standardized differences of each background factor were no greater than 0.08. For counterfactual references that passed the balance test, the same counterfactual reference relationship continued to be used in the subsequent calculations of net ecosystem service effects and impact factor replacements.

[0152] Example 3: Calculation of Ecosystem Service Functional Masses

[0153] The vineyard management unit is divided into two independent calibration units and validation units according to a preset ratio. The calibration unit is used to determine the parameters of the ecosystem service function calculation model and the weights and adjustment parameters that need to be calibrated within the ecosystem service function calculation model. The validation unit is only used to verify the stability of model errors and evaluation results and does not participate in the correction of model parameters.

[0154] 1. Product supply function

[0155] The product supply function can be determined based on the yield per unit area and grape quality correction coefficient:

[0156]

[0157] In the formula, Grape yield per unit area It is a quality correction factor determined based on soluble solids, total acid, phenolic substances, or brewing quality grade.

[0158] For counterfactual reference records originating from unplanted areas determined according to the aforementioned matching rules, when wine grapes are not planted in the corresponding area, the quality of the wine grape product supply function is determined to be 0, and no grape quality correction coefficient is calculated for the counterfactual reference record; for counterfactual reference records originating from the same spatial location during the historical period before the establishment of the vineyard, when wine grapes were not planted during the corresponding historical period, the quality of the wine grape product supply function is also determined to be 0.

[0159] When using unplanted areas as a counterfactual reference, the difference in product supply function reflects the change in supply services resulting from the planting status of wine grapes relative to the unplanted status. The impact of different cultivation management practices on product supply function is determined separately in subsequent cultivation management factor replacement calculations.

[0160] The product supply function index is used to evaluate the supply services of vineyards. It does not directly participate in the comprehensive ecological function score based on market prices, in order to reduce the impact of price fluctuations in different years on the ecological function evaluation results.

[0161] Soluble solids, total acid, and phenolic compounds were dimensionless according to the suitable ranges used for evaluating the quality of wine grapes; values ​​within the suitable range were recorded as 1, and values ​​deviating from the suitable range decreased according to their distance from the nearest range boundary. Quality correction coefficient. This is the weighted average of all dimensionless quality indicators, ranging from 0 to 1. The weight of each quality indicator is determined using calibration unit data.

[0162] 2. Moisture regulation function

[0163] Calculate the root zone water regulation amount based on the water balance within the phenological time window:

[0164]

[0165] In the formula, For effective precipitation, For effective irrigation volume, This is the actual evaporation rate. Surface runoff, This refers to the amount of deep leakage.

[0166] Effective irrigation volume As an input to the root zone water balance, it is used in the calculation of water regulation and is used to characterize the actual amount of irrigation water obtained by the grape root zone within the corresponding phenological time window. Considering that excessive irrigation can increase root zone water input, it may also exacerbate regional water resource consumption. Therefore, irrigation water intake pressure and groundwater level decline are not considered as additional ecosystem service gains, but are deducted as ecological pressure constraints in the subsequent comprehensive assessment.

[0167] The water regulation value represents the net water regulation status of the grape root zone within the corresponding phenological time window; an increase in its value indicates improved water conditions in the root zone. Water consumption caused by irrigation and groundwater pressure are not included in the water regulation function gain and will be calculated separately in the subsequent ecological pressure constraint.

[0168] 3. Soil conservation function

[0169] Soil retention capacity is calculated using the difference between potential erosion and actual erosion.

[0170]

[0171] in:

[0172]

[0173]

[0174] In the formula, As the erosivity factor of rainfall, As a soil erodibility factor, For slope length and slope factor, For vegetation cover and management factors, Factors related to soil and water conservation measures.

[0175] In this embodiment, the soil conservation function is used to characterize the amount of soil retained under hydraulic erosion conditions; the wind erosion reduction under wind erosion conditions is separately included in the windbreak and sand fixation function and is not included in the soil conservation function.

[0176] 4. Windbreak and sand-fixing function

[0177] Determined based on the difference between the potential wind erosion under bare land conditions and the existing wind erosion of the vineyard:

[0178]

[0179] In the formula, Quality of materials for windbreak and sand fixation; This represents the potential wind erosion under bare land conditions. This represents the current wind erosion level of the vineyard.

[0180] Potential wind erosion is calculated based on wind speed, soil particle composition, soil moisture content, and bare land length; current wind erosion also takes into account grape canopy coverage, inter-row coverage, windbreak height, windbreak permeability, and belt spacing.

[0181] Potential wind erosion and current wind erosion are calculated using the same caliber. When calculating potential wind erosion, vegetation and protection conditions are set according to the bare land condition. When calculating current wind erosion, the actual canopy cover, inter-row cover, and windbreak parameters of the vineyard are used. All other corresponding inputs maintain the same data caliber.

[0182] 5. Net carbon fixation function

[0183]

[0184] In the formula, Changes in carbon storage of grapevines, inter-row vegetation, and windbreak vegetation; Changes in soil organic carbon storage; The equivalent carbon emissions from irrigation energy use, mechanized operations, fertilizer input, and pesticide use.

[0185] The portion of branches that are shredded and returned to the field is included in the soil carbon input and the current decomposition release based on the carbon content and decomposition coefficient of the branches.

[0186] 6. Habitat maintenance function

[0187] Based on the spatial distribution of vineyards, inter-row vegetation, windbreaks, natural patches, and water bodies, habitat quality, patch connectivity, and disturbance intensity are calculated. For vineyards connected by windbreaks or natural vegetation, the contribution to habitat connectivity is increased; for areas with roads, continuous bare land, and high-intensity pesticide use, the habitat degradation coefficient is increased.

[0188] Habitat quality, patch connectivity, and disturbance intensity were normalized to the range of 0 to 1, and habitat maintenance function was determined according to the following formula:

[0189]

[0190] In the formula, To maintain habitat function; The normalized habitat quality; The normalized patch connectivity; The normalized interference intensity; , and The weights for the corresponding indicators are all non-negative numbers and satisfy the following conditions:

[0191]

[0192] In this embodiment, , and The same weights are applied; when calibration data is available, the weights can also be determined based on the calibration data.

[0193] All service functions subsequently involved in the calculation of the net effect of ecosystem services, the first-order contribution of the model, and the second-order non-additive coupled contribution of the model are processed in the direction that the increase of the numerical value indicates the enhancement of the corresponding service function.

[0194] When value accounting is required, a separate value layer can be established after the above-mentioned quality assessment, based on the market value method, substitution cost method, or recovery cost method. The value does not participate in the calculation of service function quality and coupling contribution.

[0195] After obtaining the biomass of each ecosystem service function corresponding to the vineyard management unit and its counterfactual reference, the two are differencing to characterize the net effect of ecosystem services relative to the environmental background conditions. The vineyard management unit in the first Within the first phenological time window The difference in the quantity of each ecosystem service function is determined by the following formula:

[0196]

[0197] In the formula, For the first The vineyard management unit in the first Within the first phenological time window The difference in the quantity of an ecosystem service function relative to a counterfactual reference; The quality of ecosystem service functions corresponding to the vineyard management unit; The quality of the corresponding counterfactual reference service function.

[0198] To reduce the impact of outlier observations on comparisons between different vineyard management units, the differences in ecosystem service functions were robustly standardized using the 5th and 95th quantiles.

[0199]

[0200] In the formula, Net effect of ecosystem services after robust standardization; and Each of the vineyard management units in the 1st Within the first phenological time window The 95th and 5th percentiles of the differences in ecosystem service functions; In order to be with the first Within the first phenological time window The quantities of each ecosystem service function have the same dimensions and are used to prevent positive numbers with a denominator of zero.

[0201] The truncation function is:

[0202]

[0203] After the above treatment, The value ranges from -1 to 1. A positive value indicates that the corresponding ecosystem service function has a net gain relative to the counterfactual reference, while a negative value indicates a net loss.

[0204] Example 4: Multi-factor Coupling Contribution Calculation

[0205] Establish an impact factor vector:

[0206]

[0207] In the formula, The total number of impact factors. For the first One influencing factor, .

[0208] Impact factors are divided into three categories:

[0209] Natural background factors include precipitation, temperature, wind speed, slope, soil texture, and groundwater depth;

[0210] Ecological state factors include vegetation cover, leaf area index, soil organic matter, soil salinity, and habitat connectivity;

[0211] Cultivation management factors include irrigation amount, fertilizer amount, inter-row cover ratio, soil tillage method, pruning intensity and windbreak parameters.

[0212] For continuous cultivation management factors, when using a single counterfactual reference record, the impact factor value corresponding to that record is used as the counterfactual reference value; when using a composite counterfactual reference, the composite counterfactual reference value of the continuous impact factor is calculated according to the composite weight of the candidate reference records; when the reference record does not have a corresponding management activity, the baseline value corresponding to the state of no management input is used. For categorical cultivation management factors such as soil tillage methods and pruning methods, categorical coding is used. The coding results are only used to distinguish different categorical states and do not represent the continuous change relationship between different categories based on the numerical value of the coded value; when using a single counterfactual reference record, the categorical state corresponding to that record is used as the counterfactual reference state; when using a composite counterfactual reference, the categorical state corresponding to the candidate reference record with the largest composite weight is used as the counterfactual reference state; when the corresponding candidate reference record does not have that type of management activity, the pre-coded baseline state of no management is used. All replacement values ​​are within the valid input range of the corresponding ecosystem service function calculation model.

[0213] When a certain influencing factor is not used as input to the corresponding ecosystem service function calculation model, the replacement of the influencing factor does not change the corresponding model output, and its first-order contribution is 0; when the replacement of a two-factor containing the influencing factor does not produce additional output changes, the corresponding second-order non-additive coupling contribution is 0.

[0214] For the The vineyard management unit in the first Vector of influencing factors within a phenological time window , No. The impact factor on the first The first-order contribution of an ecosystem service function model is determined according to the following formula:

[0215]

[0216] No. The first influence factor and the first The second-order non-additive coupling contribution of each influencing factor is determined according to the following formula:

[0217]

[0218] In the formula, For the first Within the first phenological time window Computational models for ecosystem service functions; Indicates the first The impact factor vector is formed by replacing each impact factor with its corresponding counterfactual reference value; Indicates the first The impact factor vector is formed by replacing each impact factor with its corresponding counterfactual reference value; Indicates that at the same time the first The and the first The resulting vector of influence factors is formed by replacing each influence factor with its corresponding counterfactual reference value.

[0219] when When the corresponding combination of influencing factors is determined to be a synergistic enhancement combination; when When this occurs, the corresponding combination of influencing factors is determined to be a mutually restrictive combination.

[0220] When the second-order non-additive coupling contribution of the model is equal to 0, it means that the corresponding combination of influencing factors has not produced a second-order non-additive effect.

[0221] To eliminate the problem of incomparable contribution intensity caused by different units of measurement for different ecosystem service functions, the first-order contribution and the second-order non-additive coupled contribution of the model are dimensionless according to the quantile of the difference in the corresponding ecosystem service functions:

[0222]

[0223]

[0224] In the formula, The meaning is the same as the definition in the aforementioned robust standardization process of net ecosystem service effects. This represents the first-order contribution of the dimensionless model. This represents the contribution of the second-order non-additive coupling of the dimensionless model.

[0225] The first-order contribution and the second-order non-additive coupling contribution of the model are determined based on the difference in model output before and after the replacement of the influencing factors. They are used to compare the degree of influence of different influencing factors and combinations of influencing factors on the calculation results of the corresponding ecosystem service functions.

[0226] Taking the water regulation function as an example, keeping other factors unchanged, the irrigation amount of the vineyard is replaced with the corresponding counterfactual reference value to obtain the first-order contribution of the irrigation factor in the model; then keeping other influencing factors unchanged, the irrigation amount and groundwater depth are replaced with the corresponding counterfactual reference values ​​to obtain the joint replacement result of the two factors.

[0227] Management unit M028 is located in the Helan Mountain East Foothills grape-growing region of Ningxia, corresponding to operation area 28, covering an area of ​​approximately 1.1 hm². Its central location is at 106.02307°E, 38.49719°N. The planted variety is Chardonnay, and the vineyard was established in 2010. The soil type is aeolian sandy soil, and the plot has a slope of approximately 2.5°. Drip irrigation is used. This management unit is included in the continuous monitoring scope from 2016 to 2023.

[0228] During the young fruit development period in 2023, three rounds of field surveys and soil sampling were conducted in management unit M028 on June 18, June 30, and July 10, 2023. The inter-row cover ratio was investigated using the quadrat method. Five 1 m × 1 m quadrats were set up within each management unit, and the inter-row cover was investigated within each quadrat using a 10 cm × 10 cm grid. Each quadrat was read twice. The results from each quadrat survey were then summarized to form a representative value for the inter-row cover ratio during the young fruit development period of that management unit.

[0229] Soil samples were collected at a depth of 0–20 cm. Five sampling points were set up in a quincunx pattern within the M028 management unit for each period. Equal volumes of soil samples from each point were mixed to form one composite sample. Soil organic matter was analyzed using the potassium dichromate oxidation-external heating method, according to NY / T 1121.6—2006. During the analysis, a BSA224S electronic analytical balance with a range of 220 g and a graduation of 0.1 mg was used for sample weighing. A digital display constant-temperature oil bath was used for external heating, and titration was performed using a 25 mL acid burette. Each composite soil sample was analyzed in triplicate, and the average of the three results was used as the measured soil organic matter value within the corresponding phenological time window for that management unit. The above field investigation and testing results served as the model input data for the corresponding phenological time window of the M028 management unit.

[0230] Taking the data from the fruit development stage of management unit M028 in 2023 as an example, the impact of inter-row cover and soil organic matter on soil conservation function was calculated. Based on the aforementioned field survey and soil testing results, the current value of inter-row cover ratio during the fruit development stage in this management unit was 20%, and the current value of soil organic matter was 7.2 g / kg. According to the matching results of the corresponding counterfactual reference unit, the counterfactual reference value of inter-row cover ratio was 0, and the counterfactual reference value of soil organic matter was 5.8 g / kg. The soil conservation function calculation model output was 14.2 t / hm² under the original state; after replacing only the inter-row cover ratio with the counterfactual reference value, the model output was 9.8 t / hm², corresponding to a first-order contribution of 4.4 t / hm²; after replacing only the soil organic matter with the counterfactual reference value, the model output was 10.5 t / hm², corresponding to a first-order contribution of 3.7 t / hm²; and after replacing both influencing factors simultaneously, the model output was 6.7 t / hm². The calculation results of soil conservation function under different replacement states are shown below. Figure 5 According to the aforementioned formula for calculating the second-order non-additive coupling contribution of the model, the second-order non-additive coupling contribution of this combination of influencing factors is 0.6t / hm². Therefore, this combination is determined to be a synergistic enhancement combination.

[0231] The aforementioned 20% inter-row cover ratio and 7.2 g / kg soil organic matter content were derived from field survey records and soil testing results within the corresponding phenological time window of management unit M028, respectively. 0 and 5.8 g / kg were the corresponding reference values ​​determined through counterfactual reference matching. 14.2 t / hm², 9.8 t / hm², 10.5 t / hm², and 6.7 t / hm² represent the soil conservation function model outputs obtained under the following conditions: keeping other input conditions unchanged, under the original state, replacing only the inter-row cover ratio, replacing only the soil organic matter, and replacing both influencing factors simultaneously.

[0232] If the second-order non-additive coupling contribution of irrigation volume and groundwater depth is negative, then according to the sign determination rule of the second-order non-additive coupling contribution of the model, the corresponding combination of influencing factors is determined as a mutually restrictive combination.

[0233] Taking soil conservation function as an example, the first-order contributions of inter-row cover, soil organic matter, and tillage methods can be calculated separately. When the second-order non-additive coupling contribution of the "inter-row cover-soil organic matter" model is positive, the corresponding combination of influencing factors is determined as a synergistic enhancement combination according to the sign determination rule of the second-order non-additive coupling contribution of the model.

[0234] When calculating the net carbon sequestration (NCS) function of the same management unit, the current irrigation amount is 45 m³ / mu, corresponding to a counterfactual reference value of 0; the current fertilizer application, calculated as pure nitrogen, is 28 kg / mu, also corresponding to a counterfactual reference value of 0. Under the original conditions, the NCS calculation model output is 5.72 tC / hm²; after replacing only the irrigation amount, the model output is 4.81 tC / hm², corresponding to a first-order contribution of 0.91 tC / hm²; after replacing only the fertilizer amount, the model output is 4.65 tC / hm², corresponding to a first-order contribution of 1.07 tC / hm²; and after replacing both irrigation and fertilizer amounts, the model output is 3.98 tC / hm². Therefore, the second-order non-additive coupling contribution of the model is 0.24 tC / hm², and the corresponding influencing factor combination is determined to be a synergistic enhancement combination.

[0235] The first-order contribution and the second-order non-additive coupled contribution of the model were dimensionless according to the quantile of the difference in the corresponding ecosystem service functions, so as to eliminate the dimensional differences between product supply, soil retention, net carbon sequestration and other ecosystem service functions.

[0236] For the The vineyard management unit in the first The first-order contributions of each dimensionless model within a phenological time window are summed by weighted absolute values ​​according to the corrected weights of the corresponding ecosystem service functions to obtain the comprehensive first-order contribution intensity of each influencing factor; the second-order non-additive coupling contributions of each dimensionless model are summed by weighted absolute values ​​according to the corresponding corrected weights to obtain the comprehensive coupling contribution intensity of each combination of influencing factors.

[0237] The influencing factor with the largest overall first-order contribution intensity was determined as the main controlling influencing factor, and the combination of influencing factors with the largest overall coupled contribution intensity was determined as the main controlling factor combination. This avoids directly comparing product supply, soil retention, windbreak and sand fixation, and net carbon sequestration under different units of measurement.

[0238] The combined first-order contribution intensity and the combined coupling contribution intensity are used to compare the degree of influence of different influencing factors and combinations of influencing factors; the specific direction of influence is still determined based on the positive or negative sign of the corresponding model's first-order contribution and the model's second-order non-additive coupling contribution.

[0239] Example 5: Credibility Correction and Comprehensive Evaluation

[0240] This embodiment uses field-measured data from 20 independent verification units (M045–M064) to verify the calculated results of corresponding ecosystem service functions. Units M045–M064 are all located in the grape-growing region of the eastern foothills of the Helan Mountains in Ningxia and remained fixed during the monitoring period; the location and boundaries of each verification unit remained unchanged throughout the corresponding monitoring year. The verification units are spatially independent of the calibration units (M001–M044) used for model parameter determination. The numbering and division relationship between the calibration units (M001–M044) and the independent verification units (M045–M064) is detailed below. Figure 6 .

[0241] For the water regulation function, soil moisture content was observed using a TRIME-PICO IPH time-domain reflectometry soil moisture meter from 2018 to 2023, and corresponding irrigation water volume was recorded using an LXS-50E rotary wet scrubber. Each validation unit had one representative soil moisture profile observation point, with stratified measurements taken at soil layers of 0–20 cm, 20–40 cm, and 40–60 cm. Each layer was read twice, and the average value was taken. Observations were conducted annually according to five phenological time windows: dormancy, budding to flowering, young fruit development, color change, and ripening / harvesting. Corresponding irrigation and drainage data were also recorded.

[0242] To verify the soil conservation function, ground verification was conducted using fixed erosion observation plots from 2020 to 2023. Each verification unit had one fixed erosion observation plot, with three slope positions (upper, middle, and lower) within each plot. Three φ5 mm × 30 cm steel erosion needles were installed at each slope position, for a total of nine needles. A 602-01S type vernier caliper with a range of 0–150 mm and a graduation of 0.02 mm was used to record water erosion and deposition using 20 m × 5 m runoff plots and collection buckets. One observation was conducted after each effective precipitation event, with each needle reading repeated twice and the average value taken. The observation period was the erosion season from May to October each year, and observations were completed within 24 hours of an effective precipitation event with a minimum of 10 mm.

[0243] To assess habitat maintenance function, field surveys were conducted in verification units M045–M064 from June to August 2023. One 50-meter transect was established in each verification unit, between the grapevines and in the windbreak belt, for a total of two transects per unit, resulting in 40 transects across the 20 verification units (M045–M064). A 50-meter measuring tape was used for transect placement, and a Hunter II 8×42 binoculars were used for bird observation. Insect and vegetation surveys were conducted using Magnet nets and 1-meter × 1-meter transect frames. Field surveys were conducted once each in June, July, and August for each transect, for a total of three surveys.

[0244] Multiple observations within the same year for the same validation unit are aggregated according to the data processing method of the corresponding service function. The aggregated result is used as the measured value for that unit in that year, and is paired with the model calculation value for the same year to form a set of paired data. Therefore, the water regulation function consists of 120 sets of paired data from 20 validation units over 6 years from 2018 to 2023; the soil conservation function consists of 80 sets of paired data from 20 validation units over 4 years from 2020 to 2023; and the habitat maintenance function consists of 20 sets of paired data from 20 validation units in 2023.

[0245] The product supply function was verified using actual grape yield per unit area, soluble solids, total acid and phenolic substances test data.

[0246] For soil conservation, data from ground erosion needles or water erosion sedimentation were used for verification; for windbreak and sand fixation, data from sand collection troughs, near-surface sand transport, and wind speed profiles were used for verification; for water regulation, soil moisture profiles and irrigation and drainage records were used for verification; for net carbon sequestration, plant biomass, soil organic carbon, cultivation management energy consumption, and input usage records were used for verification; and for habitat maintenance, data from transect surveys, insect surveys, or bird surveys were used for verification.

[0247] The measured data from the aforementioned independent validation units were used to verify the calculation results of some ecosystem service functions. Each set of validation data was paired one-to-one with the corresponding model calculation values ​​according to the validation unit number and year; the validation data were not used in determining the corresponding model parameters.

[0248] For the water regulation function, 120 sets of annual paired data from 20 validation units (M045–M064) from 2018 to 2023 were used. The RMSE between calculated and measured values ​​was 3.57 mm, the MAE was 2.81 mm, and the coefficient of determination (R²) was 0.86. For the soil conservation function, 80 sets of annual paired data from 20 validation units from 2020 to 2023 were used. The RMSE was 1.21 t / hm², the MAE was 0.95 t / hm², and the coefficient of determination (R²) was 0.81. For the habitat maintenance function, 20 sets of annual paired data from 20 validation units in 2023 were used. The RMSE was 0.059, the MAE was 0.047, and the coefficient of determination (R²) was 0.72. A summary of the independent validation results for water regulation, soil conservation, and habitat maintenance functions can be found in [link to relevant documentation]. Figure 4 .

[0249] The above cross-year validation results are used to examine the overall calculation error of the corresponding ecosystem service function calculation model on independent validation data, and are not directly used as parameters for calculating the confidence coefficient of a single phenological time window. (Calculation of the...) Within the first phenological time window When determining the confidence coefficient of an ecosystem service function, the root mean square error between the model-calculated value and the ground-measured value, as well as the standard deviation of the ground-measured value, are determined based on the paired data of the validation unit within the phenological time window.

[0250] For the Within the first phenological time window For each ecosystem service function, a reliability coefficient is determined based on data missing rate, model calculation error, and scale consistency:

[0251]

[0252] In the formula, This is the credibility coefficient; This corresponds to the missing data rate; The root mean square error between the model-calculated values ​​and the ground-measured values; The standard deviation of the measured values ​​on the ground; The scaling consistency coefficient; In order to be with the first Within the first phenological time window The quantities of each ecosystem service function have the same dimensions and are used to prevent positive numbers with a denominator of zero.

[0253] The data missing rate is determined by the proportion of missing records and records not included in the calculation due to data quality issues to the total number of records that should be included in the calculation within the corresponding ecosystem service function phenological time window.

[0254] In the data processing of this embodiment, the data missing rate is calculated based on the original records before interpolation. For data with a missing rate of less than 10%, in order to complete the calculation within the corresponding phenological time window, within the data range of the calibration unit and the verification unit, the three management units with valid records that are spatially closest within the same phenological time window are selected, and the missing values ​​are supplemented using an inverse distance weighting method; the supplemented data does not change the original statistical result of the missing rate. When the missing rate of the corresponding data is not less than 10%, the corresponding indicator of the management unit within the corresponding phenological time window is marked as invalid data and is not used for this model calculation.

[0255] Anomalies were initially screened using a 3-interquartile range. When the data met the criteria... or At that time, the corresponding original records are reviewed, among which Records confirmed upon review to be due to instrument malfunction, communication abnormality, data entry error, or detection anomaly will be removed and included in the data missing rate; records verified to be due to actual environmental changes will be retained.

[0256] The basic weights of ecosystem service functions were determined for each phenological time window. The net effect of ecosystem services after robust standardization was then considered. The value range is from -1 to 1. First, convert it to the range of 0 to 1 using the following formula:

[0257]

[0258] Then, using all vineyard management units within the same phenological time window... As input data for the entropy weight method, the basic weights of each ecosystem service function are calculated. .when When the corresponding entropy value is calculated, it is treated as 0.

[0259] When the Ecosystem service functions in the first All vineyard management units within a phenological time window When the sum of the entropy weight differences of all ecosystem service functions is 0, the entropy weight difference coefficient of that ecosystem service function is set to 0; when the entropy weight difference coefficients of all ecosystem service functions within the phenological time window are all 0, the basic weights of each ecosystem service function are all set to 0. .

[0260] After obtaining the basic weights Then, the basic weights are adjusted using the corresponding credibility coefficients:

[0261]

[0262] In the formula, For the first Within the first phenological time window Adjusted weights for ecosystem service functions; These serve as the basic weights for corresponding ecosystem service functions; Summation index for ecosystem service functions ; The number of functions that serve the ecosystem.

[0263] when

[0264]

[0265] When, it indicates the first If no ecosystem service function with valid credibility exists within a phenological time window, the multi-factor coupling assessment value will not be calculated for that phenological time window.

[0266] When there are few valid measured records of an ecosystem service function, and it exhibits a high data missing rate or large model calculation errors, its reliability coefficient decreases accordingly. Even if the indicator obtains a large base weight through the entropy weight method, its adjusted weight is still constrained by the reliability coefficient.

[0267] In this embodiment, the base weights are used to reflect the degree of differentiation between different vineyard management units for each ecosystem service function, and are not used as a direct ranking result of the ecological value of each ecosystem service function. When there are no effective differences among all ecosystem service functions, equal weighting is applied.

[0268] For ecosystem service functions participating in the scale consistency test, their calculation results at at least two available spatial resolutions are aggregated into the same vineyard management unit. The rank correlation coefficients of the evaluation results at any two available spatial resolutions are calculated, and the average of each rank correlation coefficient is truncated to the range of 0 to 1 as the scale consistency coefficient. When an ecosystem service function has no valid measured samples within the corresponding phenological time window, its confidence coefficient is set to 0, and this service function is not included in the comprehensive score for that phenological time window.

[0269] The reliability coefficients for each ecosystem service function were determined separately according to phenological time windows. For the first... Within the first phenological time window For each ecosystem service function, the missing data rate is determined based on the original missing records within the phenological time window and the records removed after data quality verification. The root mean square error between the model calculated values ​​and the ground-measured values, as well as the standard deviation of the ground-measured values, are determined based on paired data from the validation units within the phenological time window. The scale consistency coefficient is determined based on the rank correlation results of the calculation results at different spatial resolutions after being aggregated into the same vineyard management unit. Thus, the data completeness, model calculation error, and scale stability of each phenological stage are reflected in their respective reliability coefficients.

[0270] The dimensionless conversion of each factor's contribution is further processed. For the same ecosystem service function, the proportion of each combination in the total second-order non-additive effect is first determined based on the absolute value of the second-order non-additive coupling contribution of each factor combination. Then, combining the second-order non-additive coupling contribution of the combination with the first-order contributions of the two corresponding factors, the coupling correction effect of the combination on the service function is calculated. Finally, the calculation results of each factor combination are summed to obtain the factor coupling correction amount. .

[0271]

[0272] In the formula, , and , All are impact factor indices in the impact factor set, where, , Used for summing factor combinations. , Used for normalized summation in the denominator, and , The total number of impact factors; To prevent dimensionless positive numbers with zero denominators during the calculation of factor coupling corrections, the absolute value of the second-order non-additive coupling contribution is used to determine the proportion of each influencing factor combination. This ensures that combinations with significant coupling effects have a corresponding degree of influence when summarizing, while retaining the positive and negative directions of the second-order non-additive coupling contribution itself.

[0273] For a given combination of influencing factors, when its second-order non-additive coupling contribution is larger than the first-order contributions of each of the two influencing factors, the influence of the combination on the factor coupling correction amount increases accordingly; when the second-order non-additive coupling contribution is close to zero, the influence of the combination on the correction result also decreases.

[0274] After summarizing the combinations of various impact factors, A positive value indicates that the overall second-order non-additive effects under this ecosystem service function are positively modified. A negative value indicates an overall negative correction. This correction is not considered a standalone result in the ecosystem service assessment but is included in the multi-factor coupled assessment along with the net ecosystem service effect; the coupling contribution adjustment coefficient... Used to control the extent to which the factor coupling correction plays a role in the overall evaluation.

[0275] In this embodiment, ecological stress factors include one or more of the following: groundwater level decline, soil salinity accumulation, fertilizer and pesticide application intensity, and habitat fragmentation.

[0276] Ecological pressure constraints are calculated using the following formula:

[0277]

[0278] In the formula, This refers to the sequence number of the ecological stress factor. The number of ecological pressure factors involved in the calculation of ecological pressure constraints; For the first Ecological stress factors in the first The first vineyard management unit, the first Measured values ​​within a phenological time window; For the first The threshold at which an ecological stress factor begins to generate significant ecological stress; For the first The upper limit of the ecological stress factor, and ; For the first The weights of each ecological pressure factor are non-negative and satisfy the following conditions:

[0279]

[0280] when When, the normalized stress value of this ecological stress factor is taken as 0; when At that time, the normalized stress value of the ecological stress factor is taken as 1.

[0281] Ecological pressure constraints are used to reflect resource consumption and environmental pressure. They are deducted during the multi-factor coupling assessment stage and are not counted repeatedly as ecosystem service function mass.

[0282] Coupling contribution adjustment coefficient and ecological pressure constraint coefficient All parameters are taken within the range of 0 to 1, and the parameter combination is calibrated using a calibration unit. Once the parameter combination is determined, it remains unchanged, and then the calculation error of the corresponding ecosystem service function and the stability of the final evaluation results are verified using a validation unit.

[0283] After obtaining the net effect of ecosystem services, factor coupling correction, correction weights, and ecological pressure constraints, the first... The vineyard management unit in the first The multi-factor coupling evaluation value within a phenological time window is determined according to the following formula:

[0284]

[0285] In the formula, For the first The vineyard management unit in the first Multi-factor coupled evaluation values ​​within a phenological time window; Net effect on ecosystem services; This is the factor coupling correction amount; The revised weights for ecosystem service functions; The coefficient contributes to the coupling adjustment. This represents the ecological pressure constraint coefficient. This is a constraint on ecological pressure.

[0286] When a management unit has high net carbon sequestration, soil conservation, and product supply functions, but the groundwater level continues to decline or the soil salinity exceeds the threshold, its comprehensive assessment value will be reduced accordingly to avoid classifying high water-consuming vineyards as high ecological benefit areas solely based on vegetation cover and yield.

[0287] Example 6: Restricted Service Identification and Spatial Partitioning

[0288] For the The vineyard management unit in the first The various ecosystem service functions within a phenological time window are first screened and their weights adjusted. The ecosystem service functions were evaluated, and corrected weights were calculated for each. Net effect of ecosystem services The product of the products; the ecosystem service function corresponding to the minimum value in the product is determined as the limited ecosystem service function of the vineyard management unit within the corresponding phenological time window.

[0289] The smaller the product of the adjusted weight and the net effect of ecosystem services, the weaker the net effect of the service function after weight adjustment. When the product is negative, it indicates a net loss relative to the counterfactual reference. Therefore, the ecosystem service function corresponding to the minimum product value is identified as a restricted ecosystem service function.

[0290] No. The vineyard management unit in the first The combined first-order contribution intensity of the h-th influencing factor within a phenological time window is determined according to the following formula:

[0291]

[0292] No. The first influence factor and the first The combined coupling contribution intensity of the influencing factors is determined according to the following formula:

[0293]

[0294] In the formula, To synthesize the first-order contribution intensity; This contributes to the overall coupling strength. The influence factor corresponding to the maximum value is determined as the controlling influence factor. The combination of influencing factors corresponding to the maximum value is determined as the combination of controlling factors.

[0295] Replacement simulations were performed on cultivation management factors within the allowable adjustment set, and the ecological recoverability index was determined according to the following formula:

[0296]

[0297] In the formula, For the first The vineyard management unit in the first Ecological recoverability indicators within a phenological time window; The set of permissible adjustments is defined by grape cultivation standards, the scope of control of management equipment, and historically valid management records. To allow adjustment of the vector of influence factors after replacement within the set; The multi-factor coupling evaluation value is calculated using the replaced influence factor vector; The values ​​are the multi-factor coupling evaluation values ​​before adjustment.

[0298] In the ecological recoverability simulation, each simulation only changes the selected cultivation management factors within the allowable adjustment set, while the remaining influencing factors remain in their pre-adjustment state; the multi-factor coupling evaluation value after the replacement is calculated. At that time, the ecological pressure constraint within the corresponding phenological time window The calculation results before adjustment are used, and the simulation is not recalculated with each replacement of the impact factor.

[0299] The ecological recoverability index in this embodiment is used to evaluate the potential for improvement in the assessed value obtained by adjusting cultivation management factors within the allowable adjustment set under current ecological stress conditions. Single-replacement calculations do not predict long-term changes in relatively slow-changing conditions such as groundwater and soil salinity.

[0300] For each vineyard management unit, the multi-factor coupling assessment value, comprehensive first-order contribution intensity, comprehensive coupling contribution intensity, and ecological recoverability index are calculated within each phenological time window. Restrictive services are identified, controlling factors are identified, and spatial partitioning is performed within the same phenological time window.

[0301] In this embodiment, the 50th percentile of the multi-factor coupled evaluation values ​​of all vineyard management units within the same phenological time window is used as the evaluation partition threshold. The 50th percentile of the ecological recoverability index of all vineyard management units within the same phenological time window was used as the threshold for recovery zoning. In other implementations, assessment zoning thresholds and restoration zoning thresholds may be predetermined based on historical evaluation samples or regional ecological management objectives.

[0302] When using the above quantiles to determine the zoning threshold, the resulting zoning is used to compare the relative status and management priority of different vineyard management units within the same evaluation area and the same phenological time window, and is not used as the result of determining whether the ecological environment quality meets the corresponding standards.

[0303] when and At that time, the corresponding vineyard management unit will be designated as an ecological conservation zone. These management units currently have a high level of integrated ecosystem services, and there is relatively little room for further improvement of the multi-factor coupling assessment value within the allowable adjustment set. Therefore, existing inter-row cover, windbreaks, and water and fertilizer management practices will be maintained.

[0304] when and At that time, the corresponding vineyard management unit was divided into a maintenance and enhancement zone. This type of management unit currently has a high level of integrated ecosystem services, while still having some room for management and adjustment. The corresponding cultivation management parameters were adjusted according to the main influencing factors and combinations of main factors.

[0305] when and At that time, the corresponding vineyard management unit is divided into a regulation and improvement zone. The current integrated ecosystem service level of this type of management unit is relatively low, but a large increase in the assessment value can be obtained by adjusting cultivation management factors. Based on the main influencing factors and the combination of main control factors, priority is given to adjusting irrigation amount, fertilization amount, inter-row cover ratio, soil tillage method or windbreak parameters.

[0306] when and At that time, the corresponding vineyard management units were designated as key restoration areas. These management units currently have a low integrated ecosystem service level, and the incremental assessment values ​​that can be obtained by adjusting cultivation management factors within the allowable adjustment set are relatively small. Based on the limiting ecosystem service functions and sources of ecological stress, adjustments were made to windbreaks, irrigation and drainage systems, and soil management measures.

[0307] The zoning results are based on the vineyard management unit as the smallest output object. Each management unit must output at least the multi-factor coupled assessment value, spatial zoning category, limiting ecosystem service function, main controlling influencing factor, combination of main controlling factors, ecological restoreability index, sources of ecological stress and recommended control measures within each phenological time window.

Claims

1. A multi-factor coupled assessment method for ecosystem service functions in wine grape growing areas, characterized in that, Includes the following steps: S1. Obtain ground monitoring data, remote sensing data, topographic data, soil data, meteorological data, grape growth data, and cultivation management data of the wine grape planting area to be evaluated. Divide the data into multiple phenological time windows according to the phenological period of wine grapes, and perform spatial registration and temporal alignment. S2. Based on the boundaries of vineyard plots, slope aspect change lines, irrigation and drainage zone boundaries, windbreaks, and roads and ditches, the wine grape planting area to be evaluated is divided into multiple vineyard management units. A candidate reference record set is established for each vineyard management unit. The candidate reference records include at least one of the following: data records from the same spatial location during historical periods before vineyard establishment and data records from unplanted areas within the corresponding phenological time window. For the candidate reference records corresponding to the unplanted areas, a matching distance is determined based on the background factors of topography, soil, and climate, as well as the spatial distance between the unplanted area and the vineyard management unit, wherein the spatial distance is a component of the matching distance. A single counterfactual reference record is determined based on the matching distance, soil type, landform type, climate zone category constraints, and the results of the background factor balance test. When there is no single candidate reference record that simultaneously satisfies the category constraints, the matching distance threshold, and the preset balance threshold, a synthetic counterfactual reference value is generated based on multiple candidate reference records. The matching and balance checks of the candidate reference records are performed in the following manner: ； In the formula, For the first The first vineyard management unit and the first The matching distance between the candidate reference records; The number of continuous background factors participating in the matching; The background factor number; and The first The first vineyard management unit and the first The first candidate reference record corresponding to the Background factor; For the first Interquartile range of the background factor; In order to be with the first The background factors have the same dimensions and are used to prevent positive numbers with a denominator of zero; Spatial distance; As a spatial distance reference value and ; and These are the background factor matching weights and spatial distance matching weights, respectively. The and All are non-negative numbers and satisfy: ； For the first Matching weights of the background factors, The matching weights are spatial distances, and both are non-negative numbers; The relative difference between a vineyard management unit and its corresponding counterfactual reference on the same background factor is defined as standardized difference and determined according to the following formula: ； In the formula, For the first The vineyard management unit in the first Standardized differences on the background factors; For the first The reference for the balance test of the vineyard management unit is the first... The background factor value; when the reference to be tested for balance is a single candidate reference record. For the single candidate reference record corresponding to the first The background factor value; when the reference to be tested for equilibrium is a synthetic reference. The first candidate reference record is obtained by weighting multiple candidate reference records according to their corresponding synthesis weights. Background factor value; and The meaning is the same as the definition in the aforementioned matching distance calculation; When the standardized differences of each continuous background factor are not greater than the corresponding preset balance threshold, the reference to be tested for balance is determined to have passed the continuous background factor balance test. Among the candidate reference records that simultaneously meet the constraints of soil type, landform type and climate zone category, matching distance threshold and preset balance threshold, the candidate reference record with the smallest matching distance is selected as the single counterfactual reference record. S3. Determine the ecosystem service function mass of the vineyard management unit and its corresponding counterfactual reference service function, wherein the counterfactual reference service function mass is the ecosystem service function mass corresponding to a single counterfactual reference record, or the ecosystem service function mass calculated by multiple candidate reference records according to the composite weight; form the service function vector of the vineyard management unit and the counterfactual reference service function vector, and obtain the net ecosystem service effect vector after reducing the impact of environmental background differences based on the mass difference between the two. S4. Establish a set of influencing factors consisting of natural background factors, ecological state factors, and cultivation management factors. Replace each individual influencing factor in the vineyard management unit with its corresponding counterfactual reference influencing factor value, and calculate the first-order contribution of each influencing factor to different ecosystem service functions. Then, replace both influencing factors with their corresponding counterfactual reference influencing factor values ​​simultaneously, and calculate the second-order non-additive coupling contribution of the two influencing factors. The counterfactual reference influencing factor value is derived from a single counterfactual reference record or generated from multiple candidate reference records according to a composite weight. S5. Convert the net effect of the ecosystem services into a non-negative evaluation value, and determine the basic weight of each ecosystem service function based on the non-negative evaluation values ​​of each vineyard management unit within the same phenological time window. The reliability coefficients of each ecosystem service function are determined according to the following formula: ； In the formula, For the first Within the first phenological time window Credibility coefficient of ecosystem service function; This corresponds to the missing data rate; The root mean square error between the model-calculated values ​​and the ground-measured values; The standard deviation of the measured values ​​on the ground; The scale consistency coefficient is the result of summarizing the calculation results at different spatial resolutions into the same vineyard management unit. In order to be with the first Within the first phenological time window The quantities of ecosystem service functions have the same dimensions and are used to prevent positive numbers with a denominator of zero; the data missing rate satisfy The scale consistency coefficient satisfy When the first Ecosystem service functions in the first When no valid measured samples are available within a phenological time window, the confidence coefficient is... Set to 0; The basic weights are then corrected using the credibility coefficient. S6. The first-order contribution and the second-order non-additive coupling contribution of the model are dimensionless. The factor coupling correction of each ecosystem service function is determined based on the dimensionless first-order contribution and the second-order non-additive coupling contribution. The ecological pressure constraint is determined based on the measured values ​​of the ecological pressure factors and the corresponding pressure thresholds and upper limits. The net effect vector of ecosystem services, the factor coupling correction, the corrected weights, and the ecological pressure constraint are then fused to obtain the multi-factor coupling evaluation value of each vineyard management unit in different phenological time windows. S7. Based on the corrected weights, determine the comprehensive first-order contribution intensity of each influencing factor and the comprehensive coupling contribution intensity of each combination of influencing factors. Then, determine the controlling influencing factor and the combination of controlling factors based on the comprehensive first-order contribution intensity and the comprehensive coupling contribution intensity. Determine the allowable adjustment set of cultivation management factors based on grape cultivation standards, the control range of management equipment, and historical effective management records. Perform influencing factor replacement simulations within the allowable adjustment set. Determine the ecological recoverability index based on the increment of the multi-factor coupling evaluation value before and after adjustment. Determine the evaluation zoning threshold and the recovery zoning threshold based on the multi-factor coupling evaluation value and the ecological recoverability index of each vineyard management unit within the same phenological time window. Output the spatial zoning results based on the high-low relationship between the multi-factor coupling evaluation value and the ecological recoverability index relative to the evaluation zoning threshold and the recovery zoning threshold.

2. The multi-factor coupling evaluation method according to claim 1, characterized in that: The phenological time window includes at least the dormancy period, the budding to flowering period, the young fruit development period, the color-changing period, and the mature harvest period; The remote sensing data includes one or more of the following: vegetation cover, leaf area index, net primary productivity, surface temperature, evapotranspiration, and remote sensing ecological indices. The cultivation management data includes one or more of the following: irrigation amount, fertilizer amount, pesticide application amount, inter-row coverage ratio, soil tillage method, pruning intensity, windbreak belt parameters, and orchard establishment year.

3. The multi-factor coupling evaluation method according to claim 1, characterized in that, When no single candidate reference record simultaneously satisfies the constraints of soil type, landform type, and climate zone category, as well as the matching distance threshold and the preset balance threshold, multiple candidate reference records are selected from the candidate reference records that satisfy the constraints of soil type, landform type, and climate zone category in ascending order of matching distance, and a synthetic counterfactual reference value is generated according to the following formula: ； ； In the formula, To avoid matching dimensionless positive numbers with a distance of zero; For candidate reference record index, For the first A set of candidate reference records corresponding to each vineyard management unit; For the first The composite weight of each candidate reference record; For the first The vineyard management unit in the first Within the first phenological time window The synthetic counterfactual reference quality of ecosystem service functions; For the first The candidate reference record in the first Within the first phenological time window The quality of ecosystem service functions; Ecosystem service function sequence number. This refers to the sequence number of the phenological time window; For continuous impact factors, the composite counterfactual reference impact factor value is generated according to the following formula: ； In the formula, For the first The vineyard management unit in the first Within the first phenological time window A composite counterfactual reference value for each influencing factor; For the first The candidate reference record in the first Within the first phenological time window The values ​​of the impact factors, The impact factor number; For categorical impact factors, the category state corresponding to the candidate reference record with the largest composite weight is used as the composite counterfactual reference state. When there is no management activity corresponding to the cultivation management factor to be replaced in the candidate reference record, the counterfactual reference state corresponding to the cultivation management factor is determined as the preset no management baseline state. For continuous background factors of synthetic counterfactual references, the corresponding synthetic background factor value is obtained by weighted summation according to the synthetic weight of each candidate reference record; for categorical background factors of soil type, landform type and climate zone, the category state corresponding to the candidate reference record that satisfies the category constraint is adopted. Perform a balance test on the generated synthetic counterfactual reference; When the standardized difference of any continuous background factor exceeds the corresponding preset equilibrium threshold, a new candidate reference record is selected and the synthetic counterfactual reference value is recalculated; when the standardized difference of each continuous background factor is not greater than the corresponding preset equilibrium threshold, the corresponding synthetic counterfactual reference value is determined as the counterfactual reference of the vineyard management unit.

4. The multi-factor coupling evaluation method according to claim 1, characterized in that: The ecosystem service function vector includes at least product supply function, water regulation function, soil conservation function, windbreak and sand fixation function, net carbon sequestration function, and habitat maintenance function; Among them, the soil conservation function is determined based on the difference between potential soil erosion and actual soil erosion; the windbreak and sand fixation function is determined based on the difference between potential wind erosion and current wind erosion under bare land conditions; and the net carbon sequestration function is determined jointly based on changes in vegetation carbon storage, changes in soil organic carbon, and carbon emissions from cultivation and management.

5. The multi-factor coupling evaluation method according to claim 1, characterized in that, First, determine the difference in ecosystem service functions using the following formula: ； The difference in ecosystem service functions is then robustly standardized using the following formula: ； In the formula, For the first The vineyard management unit in the first The difference in quantity of the i-th ecosystem service function relative to the counterfactual reference within a phenological time window; The quality of service functions for vineyard management units; The corresponding counterfactual reference service function quality is either the service function quality corresponding to a single counterfactual reference record or the service function quality synthesized from multiple candidate reference records. and Each of the vineyard management units in the 1st 95th and 5th quantiles of the difference in the i-th ecosystem service function within a phenological time window; In order to be with the first Within the first phenological time window The quantities of each ecosystem service function have the same dimensions and are used to prevent positive numbers with a denominator of zero; in, The truncation function has the following expression: 。 6. The multi-factor coupling evaluation method according to claim 5, characterized in that: No. The impact factor on the first The vineyard management unit in the first Within the first phenological time window The first-order contribution of an ecosystem service function model is determined according to the following formula: ； No. The first influence factor and the first The second-order non-additive coupling contribution of each influencing factor is determined according to the following formula: ； In the formula, For the first Within the first phenological time window Computational models for ecosystem service functions; For the first The vineyard management unit in the first Vector of influencing factors within a phenological time window; , The index of each impact factor in the set of impact factors; Indicates the first The impact factor vector is formed by replacing each impact factor with its corresponding counterfactual reference value; Indicates the first The impact factor vector is formed by replacing each impact factor with its corresponding counterfactual reference value; Indicates that at the same time the first The first impact factor and the first The impact factor vector is formed by replacing each impact factor with its corresponding counterfactual reference value; when When, the corresponding combination of influencing factors is determined as a synergistic enhancement combination; when When this happens, the corresponding combination of influencing factors is determined to be a mutually restrictive combination; To eliminate the influence of dimensional differences between different ecosystem service functions on the comparison results of contribution intensity, the first-order contribution and the second-order non-additive coupling contribution of the model were dimensionless according to the quantile of the difference in the corresponding ecosystem service functions: ； ； In the formula, This represents the first-order contribution of the dimensionless model. This represents the contribution of the second-order non-additive coupling to the dimensionless model. and The first Within the first phenological time window The 95th and 5th percentiles of the differences in ecosystem service functions; In order to be with the first Within the first phenological time window The quantities of each ecosystem service function have the same dimensions and are used to prevent positive numbers with a denominator of zero.

7. The multi-factor coupling evaluation method according to claim 6, characterized in that, Net ecosystem service effects corresponding to each vineyard management unit The values ​​were converted to non-negative evaluation values, and the basic weights of each ecosystem service function were determined using the entropy weight method based on the non-negative evaluation values ​​of each vineyard management unit within the same phenological time window. ; The corrected weights are determined according to the following formula: ； In the formula, For the first Within the first phenological time window Adjusted weights for ecosystem service functions; These serve as the basic weights for corresponding ecosystem service functions; Summation index for ecosystem service functions ; The number of ecosystem service functions; when ； At that time, the first If no ecosystem service function with valid credibility exists within a phenological time window, the multi-factor coupling assessment value will not be calculated for that phenological time window.

8. The multi-factor coupling evaluation method according to claim 7, characterized in that: First, calculate the following formula: The vineyard management unit in the first Factor coupling correction for the i-th ecosystem service function within a phenological time window: ； in, , Let be the sequence number of each impact factor in the set of impact factors, and , A dimensionless positive number used to prevent the denominator from being zero during the calculation of factor coupling corrections; Then calculate the multi-factor coupling evaluation value according to the following formula: ； In the formula, For the first The vineyard management unit in the first Multi-factor coupled evaluation values ​​within a phenological time window; The coefficient contributes to the coupling adjustment. This represents the ecological pressure constraint coefficient. This is a constraint on ecological pressure.

9. The multi-factor coupling evaluation method according to claim 8, characterized in that: Adjusting weights In the ecosystem service functions, the first The vineyard management unit in the first Weight adjustment within each phenological time window Net effect of ecosystem services The ecosystem service function corresponding to the minimum value of the product is determined as the limited ecosystem service function of the vineyard management unit within the phenological time window; Calculate the first according to the following formula. The combined first-order contribution strength of the influencing factors: ； In the formula, For the first The vineyard management unit in the first Within the first phenological time window The combined first-order contribution intensity of each influencing factor; The number of ecosystem service functions; For the first Within the first phenological time window Adjusted weights for ecosystem service functions; This represents the first-order contribution of the dimensionless model. Calculate the first according to the following formula. The first impact factor and the first The combined coupling contribution strength of the influencing factors: ； In the formula, For the first The vineyard management unit in the first Within the first phenological time window The first impact factor and the first The combined coupling contribution strength of each influencing factor; This represents the contribution of the second-order non-additive coupling to the dimensionless model. The overall first-order contribution intensity will be considered. The first corresponding to the maximum value One influencing factor was identified as the controlling influencing factor, and the combined contribution intensity was determined. The first corresponding to the maximum value The first impact factor and the first One influencing factor is determined to be the combination of the main controlling factors; Within the allowable adjustment set of irrigation amount, fertilizer amount, inter-row cover ratio, soil tillage method, and windbreak belt parameters, an influencing factor replacement simulation was conducted, and the ecological restoreability index of the vineyard management unit was calculated according to the following formula: ； In the formula, For the first The vineyard management unit in the first Ecological recoverability indicators within a phenological time window; The set of permissible adjustments is defined by grape cultivation standards, the scope of control of management equipment, and historically valid management records. To allow adjustment of the vector of influence factors after replacement within the set; The multi-factor coupling evaluation value is calculated using the replaced influence factor vector; The values ​​are the multi-factor coupling assessment values ​​before adjustment; In calculating the At that time, only the selected cultivation management factors within the allowed adjustment set are replaced with the adjusted values, while the unselected influencing factors and the corresponding ecological pressure constraints of the vineyard management unit within that phenological time window are not considered. The ecological pressure constraint is not recalculated with each single impact factor replacement simulation, keeping the original values ​​unchanged. Based on the multi-factor coupled assessment values ​​and ecological restoreability indicators of each vineyard management unit within the same phenological time window, the [number]th [period] was determined. The evaluation zone threshold corresponding to each phenological time window and restorative partition threshold ; when and At that time, the corresponding vineyard management unit will be divided into an ecological conservation area; when and At that time, the corresponding vineyard management unit will be divided into a maintenance and improvement zone; when and At that time, the corresponding vineyard management unit will be divided into a regulation and improvement zone; when and At that time, the corresponding vineyard management unit will be divided into key restoration areas.

10. The multi-factor coupling evaluation method according to claim 1, characterized in that: The vineyard management unit is divided into independent calibration and verification units according to a preset ratio. The parameters of the corresponding ecosystem service function calculation model are corrected by using the grape yield and quality test data, soil erosion monitoring data, wind and sand transport monitoring data, soil moisture content and irrigation and drainage records, vegetation and soil carbon storage, cultivation management energy consumption and input usage data, and biological survey data from the calibration unit. A match balance test was performed on the topographic, soil, and climate background factors between vineyard management units and single or composite counterfactual reference records or values. For continuous background factors, when the standardized variance of any continuous background factor exceeds a preset equilibrium threshold, a single counterfactual reference record is re-screened or a synthetic counterfactual reference value is recalculated. For categorical background factors, the corresponding categorical constraints are used for testing. When the standardized variance of each continuous background factor is not greater than the equilibrium threshold and each categorical background factor satisfies the corresponding categorical constraint, the model error of each ecosystem service function is calculated using the validation unit, and the corresponding confidence coefficient is determined accordingly. The validation unit is only used for testing the stability of model errors and evaluation results, does not participate in the correction of ecosystem service function calculation model parameters, and keeps the counterfactual reference relationships that have passed the equilibrium test and the ecosystem service function calculation model parameters determined by the calibration unit unchanged.

Images