A Method for Constructing an Assessment Model for the Impacts of Climate Change on Biodiversity

Abstract

This invention discloses a method for constructing an assessment model of the impact of climate change on biodiversity, specifically relating to the fields of ecological modeling and climate response analysis. Based on microclimate disturbance data and remote sensing imagery, a high-resolution feature raster is generated, ecological patches are extracted, and a heterogeneous spatial index model is established. A connectivity map and a coupled response path library are constructed. Nonlinear dimensionality reduction and spectral clustering are used to extract response patterns of typical species. Dynamic fitting is performed using biological survey data, and an adaptive adjustment function is constructed to optimize the model's prediction accuracy. Finally, a biodiversity decline trend prediction curve and a risk heatmap are generated, identifying the critical point of biodiversity collapse and outputting intervention priorities. This method can be widely applied in regional ecological early warning and protection planning.

Description

Technical Field

[0001] This invention relates to the field of ecological modeling and climate response analysis technology, specifically to a method for constructing an assessment model of the impact of climate change on biodiversity. Background Technology

[0002] With the intensification of global climate change, extreme weather events, abnormal seasonal patterns, and habitat degradation are becoming increasingly serious, especially in ecologically sensitive areas such as high-altitude mountains, tropical rainforests, and polar tundra, where biodiversity loss has increased significantly. Traditional biodiversity assessment methods mainly rely on species surveys, habitat change monitoring, and expert scoring, which suffer from problems such as low monitoring frequency, poor timeliness, insufficient regional coverage, and a lack of dynamic analysis capabilities for the coupling relationship between climate drivers and biological responses. These methods are insufficient for achieving precise predictions of micro-scale ecological disturbances and species response chains.

[0003] Some existing methods attempt to make macroscopic assessments based on remote sensing imagery or climate model data, but they often fail to characterize the impact of species within a region on local climate variables (such as diurnal temperature variations, intermittent extreme rainfall, and lower-level atmospheric conditions). The successive response characteristics (such as sudden changes in concentration) lead to significant scale mismatch and feature distortion in the model output. Especially in areas with severe ecological fragmentation and highly heterogeneous species diversity, traditional methods lack the ability to analyze the multi-factor coupled system of "micro-topography + microclimate + micro-habitat + sensitive species", making it difficult to conduct forward-looking ecological risk management and species protection early warning.

[0004] A biodiversity impact assessment method integrating microclimate variable perturbation, multidimensional habitat degradation feature extraction, and species response mechanism modeling is proposed. By constructing a refined coupled feature matrix and a nonlinear response prediction path, it enables proactive prediction and targeted intervention support for biodiversity change trends driven by climate change. Summary of the Invention

[0005] The purpose of this invention is to provide a method for constructing an assessment model of the impact of climate change on biodiversity, in order to address the shortcomings of the prior art.

[0006] To achieve the above objectives, the present invention provides the following technical solution: a method for constructing an assessment model for the impact of climate change on biodiversity, comprising:

[0007] Acquire microclimate time-series data and local habitat remote sensing image data within the target ecological area, and generate high-resolution micro-perturbation feature raster through wavelet compression and terrain convolution enhancement algorithms;

[0008] Based on feature raster extraction of local ecological patch units, a patch-transition zone-boundary heterogeneity spatial index model is established, climate-sensitive weight factors are assigned to each local ecological patch unit, and a micro-region environmental adaptability distribution map is constructed.

[0009] The adjacent patch nodes in the map are linked by the ecological connectivity probability set, and the dynamic evolution map of habitat connection is derived based on the linkage intensity distribution.

[0010] Based on the microclimate anomaly trigger points of each node in the graph, extreme variable variation sequences are extracted and climate factor response tensor matrices are constructed to form a biological response coupling path library.

[0011] Low-dimensional embedding and complex network spectral clustering are performed on response sequences of typical species in the coupling path library to generate diversity decline trend prediction curves for local regions.

[0012] Heterogeneous fitting of the predicted curves with existing biological survey databases was performed to calculate the dynamic bias factor and construct an adaptive microscale adjustment function to optimize the adaptive accuracy of the biodiversity impact assessment model.

[0013] The optimized biodiversity impact assessment model is applied to the generation of risk heat maps for target ecological areas, identifying potential biodiversity collapse thresholds and providing visualized early warning pathway sequences and intervention priorities.

[0014] Preferably, the generation of the high-resolution micro-perturbation feature grid includes:

[0015] Acquire time-series microclimate data for a fixed period within the target ecological region. The microclimate data includes regional diurnal temperature range, rainfall intensity, wind speed fluctuations, relative humidity, and surface temperature. Concentration change sequence, and a multivariate time synchronization matrix was constructed by time normalization;

[0016] High-resolution remote sensing image data consistent with the regional spatial distribution were acquired. The remote sensing data included multispectral reflectance, thermal infrared radiation index and vegetation index, and was pixel-level aligned with microclimate data through an image registration algorithm.

[0017] Discrete wavelet transform is applied to the registered multi-source data for compression and noise reduction to extract key temporal disturbance feature signals; and based on the terrain elevation model, terrain-aware convolutional kernels in convolutional neural networks are used to perform reinforcement learning on areas with severe terrain undulations, outputting a micro-perturbation feature raster layer that reflects the degree of spatial heterogeneous disturbance.

[0018] Preferably, the establishment of the patch-transition zone-boundary heterogeneity spatial index model includes:

[0019] Based on the micro-perturbation feature raster layer, local ecological patch units are extracted using a region growing algorithm and a minimum heterogeneity merging criterion. These local ecological patch units are based on the spatial continuity of the perturbation index distribution. Spatially adjacent pixel groups with perturbation gradients below a set threshold are merged to generate an ecologically significant patch layer. Profile analysis is performed on the edge regions of each extracted ecological patch to calculate the rate of change of the perturbation index gradient between it and adjacent patches. A transition zone layer is generated based on the rate of change. Furthermore, a spatial multi-scale boundary roughness evaluation method is used to construct a patch-transition zone-boundary heterogeneity spatial index model.

[0020] Preferably, the construction of the micro-region environmental adaptability distribution map includes:

[0021] Based on the mean internal disturbance, boundary heterogeneity index, area size, and spatial morphological compactness of each local ecological patch unit, the climate sensitivity weighting factor of the patch is calculated; the weighted patch layer is superimposed with the high-frequency disturbance response map of microclimate variables, and the Kriging spatial interpolation method is used to generate a micro-region environmental adaptability distribution map.

[0022] Preferably, the derivation of the habitat connectivity dynamic evolution diagram based on the linkage intensity distribution includes:

[0023] Based on the micro-region environmental adaptability distribution map, the spatial centroid coordinates of all local ecological patch units were extracted, and the regional ecological node set was constructed using the patch climate sensitivity weight factor as the node attribute. The threshold adjacency method was used to preliminarily identify adjacent patch node pairs.

[0024] For each pair of adjacent patch nodes, the probability value of ecological connectivity is calculated using an exponential decay function, taking into account their spatial distance, differences in patch sensitivity, mean disturbance index of the interval region, and width of the transition zone.

[0025] All patch nodes and their connectivity probability relationships are represented as a weighted undirected graph structure to generate an ecological connectivity probability set graph. Beta-connectivity coefficients are introduced to calculate the connectivity strength and network centrality of each node in the graph structure.

[0026] By simulating connectivity changes driven by time-series perturbation data, and extracting the trend of connectivity probability changes based on the sliding window analysis method, a dynamic evolution map of habitat connectivity is further constructed, with time as the horizontal axis and linkage strength as the vertical axis, to describe the connectivity evolution trajectory of the ecosystem under the background of perturbation.

[0027] Preferably, the formation of the biological response coupling pathway library includes:

[0028] In the ecological connectivity probability set map, based on the historical perturbation data of each patch node, a sliding window extreme value detection algorithm is used to identify microclimate anomaly trigger points. These anomaly trigger points include temperature jumps, sudden drops in precipitation, extreme wind speed changes, and... The timing of the concentration spike event;

[0029] For each type of abnormal trigger point, extract the microclimate variable sequence within a fixed time window before and after the trigger to construct an extreme variable variation sequence set;

[0030] The extreme variable variation sequences of all nodes are integrated across variables to construct a three-dimensional climate factor response tensor matrix. The three-dimensional climate factor response tensor matrix uses "node number-variable type-time segment" as the dimension to represent the climate factor response behavior of each ecological node under the background of extreme events.

[0031] Pattern clustering and path similarity analysis methods are used to couple and model the subset of nodes with similar response trajectories in the three-dimensional climate factor response tensor matrix, generating a biological response coupling path library.

[0032] Preferably, generating the diversity decline trend prediction curve for a local area includes:

[0033] Species response sequences were selected from the biological response coupling path library, including species distribution changes, behavioral migration, reproductive interruption and population density change characteristics, and samples were screened based on response intensity, frequency and spatial coverage indicators.

[0034] After processing species response sequences through time alignment and variable normalization, low-dimensional embedding is performed using t-SNE or UMAP nonlinear dimensionality reduction algorithms to map high-dimensional ecological response behaviors to two-dimensional or three-dimensional space, preserving the local and global structural features of the response trajectory.

[0035] We construct a weighted graph structure based on response sequence similarity in a low-dimensional embedding space, and use a spectral clustering algorithm to divide the nodes in the graph structure into modules to identify species response clusters with common ecological response mechanisms.

[0036] For each response cluster, diversity decay modeling is performed. Based on historical response trends and perturbation driving factors, a local regional diversity change fitting function is constructed, and a diversity decay trend prediction curve is output.

[0037] Preferably, the step of heterogeneous fitting based on the predicted curve and an existing biological survey database to calculate the dynamic bias factor includes:

[0038] Obtain historical biological survey databases corresponding to the target ecological region, including species occurrence frequency, population density and habitat utilization indicators at different periods, and organize them into a structured dataset according to time series.

[0039] The predicted diversity decline trend curve was heterogeneously fitted with the actual observation data in the biological survey database. The prediction error at different time points was calculated using the dynamic residual evaluation method, and a time series of dynamic deviation factors was constructed.

[0040] Based on the spatial distribution characteristics and error evolution trend of the deviation factor, a microscale adaptive adjustment function is constructed. The function takes the node location, the magnitude of the prediction deviation and the degree of fluctuation of the climate factor as input variables and outputs the correction gain used to adjust the prediction curve.

[0041] The adjustment function is embedded into the original biodiversity impact assessment model.

[0042] Preferably, the identification of potential biodiversity collapse thresholds includes:

[0043] The optimized biodiversity impact assessment model is spatially mapped within the target ecological region. Raster-by-raster calculations are performed on each ecological unit within the region to output a raster layer of biodiversity vulnerability index in continuous space.

[0044] A biodiversity risk heat map is generated based on the vulnerability index layer. The kernel density estimation method is used to spatially aggregate and express the heat gradient of high vulnerability value areas, and risk level areas are set according to heat threshold classification.

[0045] Based on heat maps and combined with dynamic trend prediction curves, spatial areas with significant and continuous downward trends in vulnerable regions are identified, and their diversity decline slope and fluctuation stability index are calculated to identify and mark them as potential biodiversity collapse critical points.

[0046] The technical effects and advantages provided by the present invention in the above technical solution are as follows:

[0047] 1. This invention constructs a multi-source fusion biodiversity impact assessment model, systematically integrating microclimate disturbance characteristics, ecological patch structure information, biological response pathway behavior, and dynamic trend fitting methods to achieve a refined, multi-scale, and multi-dependent variable-driven dynamic assessment of ecosystems under the background of climate change. Compared with traditional static or single-factor prediction models, this invention has stronger spatiotemporal resolution, higher ecological response capture capability, and stronger mechanism explanatory power, and can effectively characterize the complex coupled change processes of regional ecosystems.

[0048] 2. This invention further enhances the model's adaptability to regional differences and data heterogeneity by introducing a dynamic bias factor correction mechanism and an adaptive microscale adjustment function. Furthermore, it constructs a heatmap and collapse threshold identification method in the risk identification stage, ultimately generating a visualized early warning path sequence and intervention priority layer, achieving a complete closed loop from prediction to decision support. This technical solution has broad application prospects in biodiversity conservation planning, early warning of ecologically fragile areas, and refined management of natural resources. Attached Figure Description

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

[0050] Figure 1 This is a flowchart of the method of the present invention. Detailed Implementation

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

[0052] For examples, please refer to Figure 1 As shown in this embodiment, a method for constructing an assessment model for the impact of climate change on biodiversity includes:

[0053] Acquire microclimate time-series data and local habitat remote sensing image data within the target ecological area, and generate high-resolution micro-perturbation feature raster through wavelet compression and terrain convolution enhancement algorithms;

[0054] Based on feature raster extraction of local ecological patch units, a patch-transition zone-boundary heterogeneity spatial index model is established, climate-sensitive weight factors are assigned to each local ecological patch unit, and a micro-region environmental adaptability distribution map is constructed.

[0055] The adjacent patch nodes in the map are linked by the ecological connectivity probability set, and the dynamic evolution map of habitat connection is derived based on the linkage intensity distribution.

[0056] Based on the microclimate anomaly trigger points of each node in the graph, extreme variable variation sequences are extracted and climate factor response tensor matrices are constructed to form a biological response coupling path library.

[0057] Low-dimensional embedding and complex network spectral clustering are performed on response sequences of typical species in the coupling path library to generate diversity decline trend prediction curves for local regions.

[0058] Heterogeneous fitting of the predicted curves with existing biological survey databases was performed to calculate the dynamic bias factor and construct an adaptive microscale adjustment function to optimize the adaptive accuracy of the biodiversity impact assessment model.

[0059] The optimized biodiversity impact assessment model is applied to the generation of risk heat maps for target ecological areas, identifying potential biodiversity collapse thresholds and providing visualized early warning pathway sequences and intervention priorities.

[0060] This method preferentially acquires multi-source microclimate time-series data of the target ecological region within a fixed time period (e.g., the past 10 years, sampled daily or hourly) to reflect the historical dynamics of climate disturbances in the region. The microclimate data specifically includes the following five categories of key variables:

[0061] Regional daily temperature range (i.e., the difference between the highest and lowest daily temperatures);

[0062] Rainfall intensity (precipitation amount per unit time);

[0063] Wind speed fluctuation (based on the deviation between instantaneous wind speed and daily average wind speed);

[0064] Relative humidity (indicates the water vapor content in the air);

[0065] surface Concentration change sequence (unit: ppm, using near-surface air quality monitoring data or satellite inversion data).

[0066] The above five types of microclimate data are first uniformly converted to the same time resolution (such as daily or hourly), and then time normalized using Z-score standardization, that is, the value of each variable at each time step is subtracted from the mean of the variable and then divided by its standard deviation, thereby constructing a multivariate time synchronization matrix. This matrix can be represented as a two-dimensional matrix with dimension N×T, where N represents the number of variables and T represents the total number of sampling time points.

[0067] At the spatial distribution level, to ensure the geographical consistency of raster output and further acquire remote sensing image data consistent with the spatial distribution of the target area, medium-to-high resolution (e.g., 10-meter or 30-meter) remote sensing sources are preferred. Data sources include public satellite platforms such as Landsat-8, Sentinel-2, and MODIS. The remote sensing data contains the following three typical ecological indicator variables:

[0068] Multispectral reflectance (including visible light, near-infrared, and short-wave infrared bands);

[0069] Thermal infrared radiation index (used to identify abnormal surface temperature);

[0070] The vegetation index NDVI (Normalized Difference Vegetation Index) is used to monitor changes in vegetation cover.

[0071] The remote sensing data is aligned pixel-level with the microclimate data using a geographic coordinate system-based image registration algorithm. A preferred method is "mutual information maximization + pyramid scale matching" to ensure complete spatial overlap of the multi-source data, with an error controlled within one pixel. The registered results are then uniformly converted to raster format with the same resolution as the remote sensing image.

[0072] The registered microclimate time series matrix and remote sensing image raster are feature-fused, and the Discrete Wavelet Transform (DWT) algorithm is used for compression and noise reduction. The DWT algorithm decomposes each variable sequence into two parts: approximate coefficients (low-frequency information) and detail coefficients (high-frequency perturbations), which effectively suppresses background noise while preserving key perturbation features.

[0073] In the two-dimensional image domain, the remote sensing multi-band image is decomposed into two levels using Haar wavelets;

[0074] In the time series, the microclimate variables are decomposed into 1D DWT, and the components with periodic fluctuations in the high-frequency components are selected as perturbation features.

[0075] The aforementioned perturbation features are recombined to construct a perturbation fusion tensor with dimensions X×Y×K, where X and Y are spatial resolutions and K represents the dimension of the fused features, which serves as the input layer of the neural network.

[0076] Since topographic relief has a significant impact on microclimate disturbances (such as the temperature and humidity differences between windward and leeward slopes), this method further introduces a topographic elevation model (DEM) as an auxiliary input to enhance features in areas with drastic topographic changes. By calling publicly available DEM databases (such as SRTM, ASTER GDEM, etc.), elevation information layers consistent with the spatial distribution of remote sensing images are obtained.

[0077] A terrain-aware convolutional neural network (Terrain-Aware CNN) is constructed, and a terrain-aware convolutional kernel is introduced into the standard convolutional structure. This kernel can dynamically adjust its weights based on local elevation differences as it slides.

[0078] When the standard deviation of elevation within the sliding window exceeds a set threshold (e.g., 10 meters), the high-frequency characteristic response is enhanced.

[0079] When the elevation change is gradual, the conventional convolution features are preserved.

[0080] In its implementation, the weight adjustment function of the terrain-aware convolutional kernel is an adaptive mapping function based on the ReLU activation function, where each convolutional weight is multiplied by the elevation gradient gain coefficient of its corresponding region. The network adopts a dual-channel input structure, inputting both the remote sensing perturbation tensor and the DEM elevation raster, and extracting spatial heterogeneous perturbation features through multi-scale convolution.

[0081] After the neural network completes training and feature extraction, the network output is decoded and restored to a micro-perturbation feature raster layer with the same resolution as the original image. The value of each pixel unit in this layer represents the intensity level of microclimate disturbance at that spatial location within a given time range. The unit can be the dimensionless disturbance index (DI), with a value range of [0,1]. The larger the value, the stronger the disturbance.

[0082] The output micro-perturbation feature raster layer serves as the basic input for subsequent models (such as habitat linkage probability map construction, biological response path modeling, etc.), and has technical advantages such as high temporal resolution, strong spatial heterogeneity expression ability, and good physical-ecological consistency.

[0083] First, based on the previously constructed micro-perturbation feature raster layer, local ecological patch units are extracted. The patch extraction is achieved using a region growing algorithm combined with a minimum heterogeneity merging criterion.

[0084] Region growing uses a pixel with a high perturbation value as a seed point and determines whether its 8 neighboring pixels meet the heterogeneity merging condition.

[0085] The heterogeneity merging condition is defined as follows: the difference between the pixel perturbation index and the current patch mean is less than a specified perturbation gradient threshold.

[0086] The threshold is preferably set to 0.15, meaning that if the difference between the perturbation index of a pixel and the average perturbation value of the current patch is within ±0.15, it is considered that the pixel can be merged.

[0087] After iteration, several local ecological patch unit layers that meet the spatial continuity condition are output, with non-overlapping closed regions between patches.

[0088] Subsequently, profile disturbance analysis was performed on the boundary regions of each extracted ecological patch. The specific method was as follows:

[0089] Extend 1 to 2 pixels outward along the patch boundary as an "edge buffer";

[0090] Calculate the rate of change between the pixel perturbation index within the buffer and the mean perturbation value within the patch. The rate of change is the difference in perturbation indices divided by the spatial distance.

[0091] This generates a perturbation exponential gradient change rate layer, whose value range represents the severity of the boundary transition.

[0092] Based on the aforementioned rate of change layer, the transition zone region is further identified through boundary directional scanning, and the rules are defined as follows:

[0093] If the boundary perturbation gradient between adjacent patches is greater than 0.25 and the buffer width exceeds 3 pixels, it is defined as a significant transition zone.

[0094] The transition zone area is saved as a separate layer, and the starting patch number and the target patch number are recorded.

[0095] To quantify the heterogeneity of the overall spatial structure, a patch-transition zone-boundary heterogeneity spatial index model is constructed, the core of which includes:

[0096] Boundary Roughness Index: Calculated by measuring the ratio of the perimeter of the boundary of each patch to the circumference of its equivalent circle;

[0097] Heterogeneity Intensity Index: a weighted statistical analysis based on the average rate of change of perturbation gradients across all transition zones;

[0098] Spatial structure index matrix: Using patch numbers as row and column indices, it is filled with the heterogeneity intensity and directionality between adjacent patches (direction can be encoded as...). wait).

[0099] This spatial indexing model provides a quantitative basis for subsequent assessment of patch response sensitivity to microclimate disturbances.

[0100] After obtaining the complete patch layer and boundary index model, this invention further calculates the Climate Sensitivity Weight Factor for each patch, which comprehensively considers the following four spatial attributes of the patch:

[0101] Internal perturbation mean (μ1): The average of the perturbation indices of all pixels within the patch;

[0102] Boundary heterogeneity index (μ2): the average rate of change of perturbation in the transition zone between adjacent patches;

[0103] Area size (μ3): Number of pixels within a patch;

[0104] Spatial morphological compactness (μ4): The ratio of patch area to its boundary length. Compact shapes have a higher weight than stretched or flaky shapes.

[0105] Constructing climate sensitivity weights: After standardizing the above four indicators to the [0,1] interval, a weighted sum is performed. The boundary heterogeneity index weight is set to 0.4, the mean internal disturbance weight to 0.3, the area weight to 0.2, and the morphological compactness weight to 0.1, resulting in a sensitivity score for each patch, denoted as . .

[0106] After patch weighting was completed, this layer was overlaid and analyzed with the previously generated microclimate high-frequency disturbance response map. The response map is a spatial distribution raster of multivariate disturbance frequencies, representing the frequency or intensity of extreme climate disturbances experienced at each location over a past time period.

[0107] The Kriging spatial interpolation algorithm is used to extend the discrete patch sensitivity factor to the entire region. The interpolation process takes into account the spatial correlation and directionality of the disturbance effect.

[0108] The interpolation radius is preferably 1.5 times the average side length of the patch to ensure sufficient coverage of the affected area;

[0109] Spatial interpolation weights are calculated using the covariance function, and spatial fitting is performed using the spherical variogram.

[0110] The final generated micro-region environmental adaptability distribution map is a high-resolution raster layer, where each pixel value represents the climate change adaptability or biodiversity vulnerability at that location. Higher values ​​indicate that the location is more sensitive to climate change and has lower adaptability, and should be given priority for intervention in biodiversity conservation or restoration.

[0111] Based on the generated micro-region environmental adaptability distribution map, the geometrical centroid coordinates of each local ecological patch unit are first extracted. The centroid is the weighted average coordinates of all pixels constituting the patch in two-dimensional space, calculated as follows:

[0112] The centroid X coordinate is obtained by averaging the X coordinates of all pixels within the patch, and the Y coordinate is processed in the same way.

[0113] Based on this, each patch centroid is treated as a node in the graph structure, and the node's attribute fields include:

[0114] Node ID (patch number);

[0115] Spatial location (X,Y coordinates);

[0116] The climate sensitivity weighting factor, which is the sensitivity score calculated by weighting the perturbation mean, boundary heterogeneity, area and morphological compactness, ranges from 0 to 1.

[0117] To identify patch combinations that can form ecological corridors, a threshold adjacency method was used for preliminary screening:

[0118] Calculate the pairwise Euclidean distances between the centroids of all patches;

[0119] If the center distance between two patches is less than a preset connectivity distance threshold (e.g., 500 meters), they are identified as a potential connectivity pair.

[0120] The threshold can be determined based on regional scale, species migration capacity, or remote sensing resolution, and is adjustable.

[0121] The output is a set of candidate connectivity edges consisting of adjacent patches, and the connectivity probability values ​​will be calculated later.

[0122] For each pair of candidate edges, the probability of ecological connectivity is calculated using a comprehensive exponential decay function model, considering the following four key influencing factors:

[0123] Spatial distance D: refers to the Euclidean distance between the centroids of two patches;

[0124] Sensitivity difference S: The absolute difference in climate sensitivity weighting factors between the two patches;

[0125] Disturbance drag R: The mean of the disturbance index in the region along the connecting path (i.e., the region between two patches);

[0126] Transition zone width W: If two patches share a common edge region, measure the width of their shared boundary buffer zone.

[0127] Combining the above variables, the following connectivity probability function is constructed: Connectivity Probability Where A is the normalization coefficient, ensuring The weighted parameters controlling the exponential decay rate represent the sensitivity to distance, sensitivity differences, and disturbance drag, respectively; δ is the positive enhancement factor of the transition zone, reflecting the positive effect of the ecological channel buffer zone. Recommended values: α=0.005, β=1.0, γ=0.8, δ=0.02, A=1 (adjusted after normalization).

[0128] The above probabilities are calculated for all patch pairs that satisfy the initial connectivity conditions, resulting in a relational table consisting of node pairs and their connectivity probability values.

[0129] The nodes and connections described above are organized into a weighted undirected graph structure G=(V,E,W), where: V is the set of ecological nodes; E is the potential connected edges; and W is the connectivity probability value P_ij corresponding to each edge. This graph structure is called the ecological connectivity probability set graph.

[0130] Further ecological network analysis was performed on this graph structure, focusing on introducing the following two network indicators:

[0131] Link Strength: The sum of the connectivity probabilities of all connected edges to each node;

[0132] Network centrality (Beta-connectivity coefficient): measures the degree of influence of a node on the overall connectivity of the graph structure.

[0133] Beta connectivity coefficient is defined as the ratio of the sum of the connectivity between a node and its neighboring nodes to the maximum possible connectivity weight of its local subgraph. It ranges from 0 to 1, with higher values ​​indicating more significant ecological hub node attributes.

[0134] The analysis results are output in tabular and graphical form, providing a structural basis for dynamic evolution.

[0135] Considering that climate disturbances have time-evolutionary characteristics, this invention introduces a time-series disturbance data-driven mechanism to dynamically simulate changes in connectivity probability.

[0136] Divide the past continuous time period (e.g., 10 years) into several sliding time windows (e.g., one window every 2 years);

[0137] Regenerate maps of microclimate disturbance indices and remote sensing variables within each window;

[0138] Repeat the aforementioned node identification and connectivity probability calculation process to obtain the ecological network structure under multiple time slices.

[0139] By comparing the changing trends of patch connectivity strength and connectivity probability under different time windows, a dynamic evolution diagram of habitat connectivity is plotted, with the horizontal axis representing time points and the vertical axis representing any one or more of the following ecological connectivity indicators:

[0140] Regional average connectivity;

[0141] The number of nodes in the largest connected component;

[0142] Network density;

[0143] Evolution curve of Beta connection coefficients of core nodes.

[0144] This evolution map is used to identify key breakpoints, vulnerable time windows, and structural change trends in ecological networks under disturbance conditions, and has important forward-looking ecological early warning value. For example, if the average connectivity declines by more than 20% in a certain period, it can be regarded as a "critical period of decline" in regional ecological connectivity, and should be given priority for intervention and management.

[0145] First, based on the historical microclimate disturbance data of each patch node in the ecological connectivity probability set map, anomalous climate events that may lead to sensitive responses in the ecosystem are identified. The historical disturbance data comes from more than 10 years of continuous meteorological observation records, with a time resolution of daily or hourly.

[0146] A sliding window extremum detection algorithm is used to identify outliers for each variable. The algorithm is defined as follows:

[0147] Given a time series X(t), set the sliding window length to w (e.g., 7 days);

[0148] Calculate the mean μ and standard deviation σ of the data within each sliding window;

[0149] Determine whether the current value X(t) satisfies the extreme value condition:

[0150] If X(t) > μ + 2σ (positive extreme) or (Negative extreme) is recorded as an abnormal trigger point.

[0151] This invention independently detects the following four types of extreme microclimate variables:

[0152] Temperature surge event: The continuous rate of temperature increase exceeds a set threshold (e.g., daily increase > 5°C).

[0153] Sudden rainfall event: Total 24-hour rainfall exceeds twice the average for the same period in previous years;

[0154] Extreme wind speed events: daily wind speed variations exceeding 10 m / s;

[0155] Concentration spike: Daily average concentration increased by more than 100 ppm compared to the previous week's average.

[0156] All identified anomaly trigger points are recorded for their occurrence time, variable type, and anomaly intensity, and used as time anchors for subsequent modeling.

[0157] For each anomalous trigger point, extract the microclimate variable sequence within a fixed time window before and after it. The preferred time window length is 5 days before and after the event, for a total of 11 days.

[0158] The extracted variables include temperature, precipitation, wind speed, humidity, and The concentration and other main driving factors were used to normalize all variable sequences, as follows:

[0159] For each type of variable, calculate the maximum and minimum values ​​in the historical sequence of this region;

[0160] Map the original value x to the standard value. This ensures that all variable sequences are uniformly within the range of 0 to 1.

[0161] The final sequence set of each abnormal triggering event includes:

[0162] Timeline: 5 days before the event to 5 days after the event;

[0163] Variable axis: 5 microclimate drivers;

[0164] Value axis: Normalized perturbation intensity.

[0165] Each event can be represented as a two-dimensional matrix, which will then be used to construct higher-dimensional data structures.

[0166] The variable mutation sequences corresponding to the abnormal trigger points of all the above-mentioned patch nodes are merged and organized into a three-dimensional tensor structure according to the format of "node number-variable type-time segment". The structure is defined as follows:

[0167] There are N nodes, V variables, and T time slices;

[0168] Constructing tensors , where R(i,j,k) represents the response value of the i-th node on the j-th variable at the k-th time slice.

[0169] This three-dimensional tensor reflects the perturbation trajectories of various climate factors in different ecological spatial units under the background of extreme events.

[0170] To improve the efficiency of subsequent analysis, tensor compression methods (such as principal component tensor dimensionality reduction or Tucker decomposition) are used to process the original data.

[0171] Based on the aforementioned response tensor matrix, pattern recognition and trajectory clustering techniques are used to identify subsets of climate response trajectories that are similar among different patch nodes. The analysis method includes:

[0172] Dynamic Time Warping (DTW): Used to calculate the similarity between different time series, suitable for situations where there is a time misalignment in the peak position of the response;

[0173] Multivariate K-Shape clustering: Combining DTW distance with the synchronous variation trend of variables, it forms a response type for clustering subsequences in tensors.

[0174] After clustering, each response path type represents a group of ecological nodes' common response behavior to a certain type of climate anomaly, for example:

[0175] Vegetation response is delayed after rapid warming;

[0176] Sudden drop in precipitation accompanied by a sharp drop in temperature triggers Concentration decreased;

[0177] Humidity recovers rapidly after a sudden change in wind speed, etc.

[0178] Based on this, a biological response coupling pathway library is constructed. This library uses the response pathway number as an index and records the following information:

[0179] Response pattern ID;

[0180] The set of nodes covered;

[0181] Initial exception event type;

[0182] Multivariate response trajectory;

[0183] Potential ecological impacts (such as vegetation loss, species avoidance, population loss, etc.).

[0184] The coupling path library can be used as input for subsequent diversity prediction models, or independently for ecological management applications such as climate risk zoning and species conservation priority ranking.

[0185] First, representative species response sequences are selected from the previously constructed biological response coupling pathway library as the input basis for the prediction model. These species response sequences may include, but are not limited to, the following ecological behavioral characteristics:

[0186] Changes in species distribution: shrinkage, displacement, or fluctuation of spatial distribution range;

[0187] Behavioral migration patterns: changes in daily activity range and foraging behavior routes;

[0188] Reproductive interruption signals: Increased frequency of delayed, interrupted, or failed breeding seasons;

[0189] Population density change: The upward or downward trend of population size per unit area.

[0190] To ensure data representativeness and modeling stability, species response sequence samples were screened according to the following three criteria:

[0191] Response intensity: such as the variation amplitude or extreme response frequency being greater than a set threshold (e.g., annual change rate exceeding 20%).

[0192] Response frequency: The number of times an abnormal event triggers a response is no less than 3 times during the observed period;

[0193] Spatial coverage: The response record covers at least 3 or more ecological patch nodes.

[0194] The selection results serve as the core sample set for subsequent dimensionality reduction and modeling.

[0195] Since typical species response sequences involve multivariate, long-term, and spatially cross-dimensional high-dimensional structures (with dimensions reaching hundreds), this invention performs low-dimensional embedding processing to facilitate the analysis of differences and similarities in response patterns among species.

[0196] First, time alignment and variable normalization are performed on all sample response sequences:

[0197] Time alignment employs a sliding window unification strategy, mapping all response sequences to a fixed-length time step.

[0198] Variable normalization uses the Min-Max standardization method to map the value of each variable to the interval [0,1] to eliminate dimensional differences.

[0199] Then, a nonlinear dimensionality reduction algorithm is used for embedding mapping:

[0200] Algorithms such as t-SNE (t-distributed stochastic neighbor embedding) or UMAP (Uniform Manifold Approximation and Projection) are preferred.

[0201] The algorithm maps high-dimensional response behavior to a two-dimensional or three-dimensional embedding space while keeping the local neighborhood structure in the high-dimensional space unchanged.

[0202] After mapping, each point represents a species' response sequence, and the coordinate position reflects the comprehensive performance of its ecological response characteristics.

[0203] Example of parameter settings:

[0204] For t-SNE, a learning rate of 200 and a neighborhood dimension of 30 were chosen.

[0205] For UMAP, a minimum distance of 0.1 and a neighborhood number of 15 are selected.

[0206] In a low-dimensional embedding space, a similarity-weighted graph structure with species response sequences as nodes is constructed, defined as follows:

[0207] Each node represents a response sequence;

[0208] The weight of the connecting edge is the similarity between the node and its k nearest neighbors, preferably based on Euclidean distance or DTW distance, taking the similarity as 1 divided by the distance value;

[0209] The resulting graph structure is G=(V,E,W), where V is the set of response sequences, E is the adjacent connection pairs, and W is the edge weight matrix.

[0210] Spectral clustering is applied to this graph structure to identify the hidden modular structure among the response sequences. Specifically, this includes:

[0211] Construct the Laplacian matrix of the graph;

[0212] Calculate the eigenvectors corresponding to the first k smallest eigenvalues;

[0213] The feature vectors are clustered using K-means to divide them into response clusters;

[0214] In the clustering results, each cluster represents a set of species that share common ecological response mechanisms driven by climate anomalies.

[0215] For each species response cluster, the trend of species number changes during the observation period was statistically analyzed, and trend modeling was performed by combining climate factor perturbation data, using the following regression modeling method:

[0216] Multinomial regression model: used to fit nonlinear trends, the curve form is: diversity index Where, 'a' is the quadratic coefficient, used to control the direction and degree of curvature of the curve; 'b' is the linear coefficient, used to control the slope trend of the curve, representing the linear change of the diversity index over time; 'c' is a constant; and 't' is the time node.

[0217] Exponential decay model: Applicable to abrupt decline scenarios, the model takes the form of a diversity index. ;

[0218] Perturbation factor driven model: Introduce the frequency of extreme climate events (such as the number of temperature jumps and the total amount of precipitation) as independent variables, and construct linear or nonlinear fitting functions.

[0219] After fitting, the model is evaluated through cross-validation, and the best model is selected. The curve is used as the reliable prediction output.

[0220] Finally, the modeling results for each response cluster are output as a diversity decay trend prediction curve, which is visualized as follows:

[0221] The horizontal axis represents the future forecast period (e.g., divided by quarter or year);

[0222] The vertical axis represents the number of species, diversity index, or Shannon index;

[0223] The rate of decline of the curve, the location of the inflection point, and the plateau period are used to identify the risk of local ecosystem decline.

[0224] First, a historical biological survey database of the target ecological area is obtained. This database is usually collected by ecological monitoring units, research institutions, or nature reserve management departments. Data sources include transect surveys, remote sensing species identification, automatic camera identification systems, and manual surveys.

[0225] The biological survey database should meet the following conditions:

[0226] The time span should be no less than 5 years, preferably more than 10 years;

[0227] The spatial coverage area is consistent with or has a mapping relationship with the patch node space in the evaluation model;

[0228] Data types include, but are not limited to:

[0229] Species occurrence frequency (the number of times a species is observed per unit of time);

[0230] Population density (number of individuals per unit area);

[0231] Habitat utilization (the frequency or degree of utilization of a species in a region).

[0232] The above data is organized into a structured dataset according to time series, and a uniform time resolution (such as quarterly or annual) is set to construct the following data structure:

[0233] Biological survey data Where t represents a time point, n is the total number of data points, and x represents a vector containing multiple biological indicators, for example: . For at a certain point in time Or, within the survey period, the number of times the target species was observed or the number of days it was recorded in the designated area. This refers to the number of individuals of the target species observed per unit area. The habitat utilization ratio of the target species in the region is the ratio of the actual distribution area of ​​the species to the area of ​​usable ecological patches in the surveyed area.

[0234] The predicted biodiversity decline trend curve from the model was heterogeneously fitted to the aforementioned measured data. Since the prediction results are derived from simulated response paths, while the survey data originate from actual ecological observations, there are differences in accuracy, temporal distribution, and sampling mechanisms. Therefore, a fitting method suitable for heterogeneous data fusion is required.

[0235] The preferred method is Dynamic Residual Evaluation (DRE), which is based on the following logic:

[0236] Predicted curve Compared with actual observed values Perform point-to-point comparisons;

[0237] Calculate time series residuals ;

[0238] Considering the inherent volatility of the data, a perturbation confidence interval δ is introduced into the residual calculation. The corrected deviation is defined as: dynamic deviation factor. ;in, σ_obs(t) represents the standard deviation of the predicted values ​​over the specified period, and B(t) is the standard deviation of the observed values. B(t) is a dimensionless coefficient used to reflect the degree of prediction bias during the specified period.

[0239] The final output is the time series of the dynamic bias factor. This is used for subsequent spatial correction and model adjustment.

[0240] Based on the distribution characteristics of the dynamic deviation factor, a micro-region response correction function is constructed on a spatial scale. This function uses patch nodes in the model as units, and an adaptive micro-scale adjustment function A(x) is constructed for each node, defined as follows: Where: x is an ecological patch node; p represents the spatial location of the node (such as centroid coordinates or grid number); b is the deviation factor B(t) corresponding to the node in the current time window; v is the degree of fluctuation of microclimate variables in the region corresponding to the node (such as temperature variance, precipitation frequency index, etc.); the function form f can be an adaptive form such as linear combination, local weighted regression (LOESS) or radial basis function network (RBFN).

[0241] The preferred construction method is:

[0242] Collect the deviation factor and climate variable fluctuation values ​​for the past N periods at each node;

[0243] Construct the input matrix ;

[0244] The corresponding output is the corrected gain value ΔY for each period;

[0245] The function f is trained using machine learning regression methods (such as gradient boosting trees, SVR, etc.) so that the input M can accurately output ΔY;

[0246] A space-time gain correction mapping model for the predicted curve is formed.

[0247] The aforementioned adjustment function A(x) is embedded into the existing biodiversity impact assessment model as a post-processing correction module or a closed-loop feedback channel. Integration methods include:

[0248] For each time point in the prediction curve The corresponding node's correction gain value is superimposed. The final corrected prediction value is obtained as follows: ;in The optimized predicted value is obtained by introducing an adjustment function as a regularization term in model training to dynamically correct the initial parameters at the micro-region level, thereby achieving adaptive evolutionary updates of the model under spatial heterogeneity.

[0249] During the model evaluation phase, the output of the adjustment function is introduced as a third type of input variable to improve the error interpretation capability and model generalization capability.

[0250] First, the optimized biodiversity impact assessment model is deployed across the entire target ecological region, performing grid-by-grid calculations on each ecological unit (such as local ecological patches or basic grid units). Model inputs include:

[0251] Regional microclimate disturbance index grid;

[0252] Response path coupling strength;

[0253] Dynamically adjust the function gain value;

[0254] Biodiversity response probability output.

[0255] Based on the response strength and sensitivity indices output by the model, the biodiversity vulnerability index for each spatial unit is calculated, defined as follows:

[0256] Vulnerability Index ;

[0257] V(x,t) represents the vulnerability index at location x and time t;

[0258] S(x,t) represents the sensitivity weight factor of this unit (from the response model);

[0259] R(x,t) represents the current disturbance intensity value (such as extreme temperature fluctuations, drastic changes in precipitation, etc.).

[0260] A(x,t) is the gain value of the adaptive adjustment function; w1, w2, and w3 are weighting coefficients, preferably set as w1=0.4, w2=0.5, and w3=0.1.

[0261] All calculation results are output as a vulnerability index raster layer in a continuous space, with each cell value ranging from 0 to 1. Higher values ​​indicate greater sensitivity to disturbances and a greater susceptibility to diversity loss.

[0262] After obtaining the vulnerability index raster layer, a regional biodiversity risk heat map is generated using the kernel density estimation (KDE) method. This method estimates the "risk density" of each point in space by aggregating cell values ​​within a spatial sliding window.

[0263] The kernel density function is defined as follows: risk density Where x is the current position; For the positions of other pixels within the sliding window; The vulnerability value at that location is given by K; K is the Gaussian kernel function, and h is the bandwidth (sliding radius). Preferably, h is between 200 meters and 500 meters to cover the patch scale.

[0264] Based on the density of heat distribution, the region is divided into different risk levels (e.g., extremely high, high, medium, low). The specific classification threshold can be set based on the quantiles (e.g., P90, P75, P50) of the global V(x) distribution. The output heatmap is a visual raster layer with color gradients.

[0265] By combining the heat map with the previously generated diversity decline trend prediction curve, the changing trends in each risk area are analyzed by spatial overlay, and spatial areas with significant continuous decline characteristics are screened out.

[0266] The criteria for judging a continuous decline include:

[0267] The diversity index declined over three consecutive time windows (e.g., three consecutive years of negative growth).

[0268] The annual average decay slope exceeds the set threshold (e.g., a decrease of more than 0.05 per year).

[0269] If the standard deviation σ of the predicted curve fluctuation is less than a given value (e.g., σ < 0.02), it indicates a stable trend and non-random fluctuations.

[0270] Based on the above judgment, areas that meet the criteria are marked as potential biodiversity collapse thresholds and recorded:

[0271] Start time of change;

[0272] Current diversity index;

[0273] Historical minimum and maximum values;

[0274] Peripheral response path strength and coupling degree.

[0275] Meanwhile, the identified critical points are clustered and grouped to facilitate the formation of intervention blocks, and a vector layer is output in the form of boundaries for spatial decision support.

[0276] Time series backtracking analysis is performed on the identified critical point regions to construct risk evolution paths. The specific steps are as follows:

[0277] Track the historical intensity sequence of disturbances at this point (temperature, precipitation, etc.).

[0278] Track its coupling path response trajectory (e.g., whether it undergoes multi-species synchronous response);

[0279] Record its structural position in the network (whether it is a core node, a bridge node, etc.).

[0280] Based on the above information, a risk evolution path sequence for each critical point is constructed, with the time axis representing the order of events, node attributes representing the state of biodiversity, and edges representing the relationship between response paths and disturbances.

[0281] The following three types of indicators are further introduced to prioritize interventions:

[0282] Patch connectivity: the strength of connection with other ecological units. The higher the connectivity, the wider the impact after intervention.

[0283] Response coupling path strength: The more frequently this region appears in the coupling path library, the more likely it is to be a critical node in the response;

[0284] Ecological importance weight: Weights are assigned based on external assessment data such as national key protection level and ecosystem service value.

[0285] A multi-factor weighted scoring model is used to rank all critical points and output an intervention priority layer for decision support in ecological restoration, species migration corridor planning, and disaster early warning.

[0286] The above description is merely a specific embodiment of this application, but the scope of protection of this application is not limited thereto. Any changes or substitutions that can be easily conceived by those skilled in the art within the scope of the technology disclosed in this application should be included within the scope of protection of this application.

Claims

1. A method of constructing an assessment model of the impact of climate change on biodiversity, characterized by comprising: The method comprises the following steps: ​ Obtain microclimate time series data and local habitat remote sensing image data in the target ecological region, and generate high-resolution micro-perturbation feature grid through wavelet compression and terrain convolution enhancement algorithm; Extract local ecological patch unit based on the feature grid, establish patch-transition zone-boundary heterogeneity spatial index model, assign climate sensitivity weight factor to each local ecological patch unit, and construct micro-region environmental adaptability distribution atlas; Link adjacent patch nodes in the atlas through ecological connectivity probability set, and derive habitat connection dynamic evolution map based on the connectivity strength distribution; Extract extreme variable variation sequence and construct climate factor response tensor matrix based on the microclimate anomaly trigger points of each node in the atlas, and form a biological response coupling path library, which comprises: In the ecological connectivity probability set atlas, identify the microclimate anomaly trigger points based on the disturbance history data of each patch node using the sliding window extreme value detection algorithm. The anomaly trigger points include temperature jump, sudden rainfall, wind speed extreme change and CO2 concentration surge event time nodes; For each type of anomaly trigger point, extract the microclimate variable sequence within a fixed time window before and after the trigger point, and construct an extreme variable variation sequence set; Integrate all node extreme variable variation sequences across variables to construct a three-dimensional climate factor response tensor matrix. The three-dimensional climate factor response tensor matrix has "node number-variable type-time segment" as the dimension, representing the climate factor response behavior of each ecological node under extreme event background; Use pattern clustering and path similarity analysis method to couple the nodes in the three-dimensional climate factor response tensor matrix with similar response trajectories, and generate a biological response coupling path library; Perform low-dimensional embedding and complex network spectral clustering on the typical species response sequence in the coupling path library to generate a diversity decay trend prediction curve for the local region; Based on the prediction curve and the existing biological survey database, calculate the dynamic deviation factor, and construct an adaptive micro-scale adjustment function to optimize the adaptive accuracy of the biological diversity impact assessment model; Based on the spatial distribution characteristics and error evolution trend of the deviation factor, construct a micro-scale adaptive adjustment function. The function takes node position, prediction deviation size and climate factor fluctuation degree as input variables, and outputs a correction gain for adjusting the prediction curve. Apply the optimized biological diversity impact assessment model to the risk heat map generation of the target ecological region, identify the potential biological diversity collapse critical point, and provide visual warning path sequence and intervention priority.

2. The method of claim 1, wherein the method is characterized by: The method for generating high-resolution micro-perturbation feature grid comprises the following steps: Obtain microclimate time series data in the target ecological region within a fixed time period. The microclimate data includes regional daily temperature difference, rainfall intensity, wind speed fluctuation, relative humidity and surface CO2 concentration change sequence, and a multivariate time synchronization matrix is constructed through time normalization processing; Obtain high-resolution remote sensing image data consistent with the spatial distribution of the region. The remote sensing data includes multispectral reflectance, thermal infrared radiation index and vegetation index, and the pixel-level alignment is performed through image registration algorithm and microclimate data; The registered multi-source data is compressed and denoised by applying discrete wavelet transform to extract key time-series disturbance feature signals; and a terrain perception convolution kernel in a convolutional neural network is used based on a terrain elevation model to perform enhanced learning on areas with dramatic terrain fluctuations, and output a micro-disturbance feature raster layer reflecting the degree of spatial heterogeneity disturbance.

3. The method of claim 2, wherein the method further comprises: determining a relationship between the climate change and the biological diversity. The model of the patch-transition zone-boundary heterogeneity space index is established, and includes: Based on the micro-disturbance feature raster layer, a region growing algorithm and a minimum heterogeneity merging criterion are used to extract a local ecological patch unit, the local ecological patch unit is based on the spatial continuity of the disturbance index distribution, and a group of pixels that are spatially adjacent and have a disturbance gradient lower than a set threshold are merged to generate a patch layer with ecological significance; a profile analysis is performed on the edge area of each extracted ecological patch to calculate the disturbance index gradient change rate between the patch and the adjacent patch, and a transition zone layer is generated based on the change rate, and a spatial multi-scale boundary roughness evaluation method is further used to construct the patch-transition zone-boundary heterogeneity space index model.

4. The method of claim 3, wherein the method further comprises: determining the impact of climate change on biodiversity based on the determined relationship between the climate change and the biodiversity. The micro-environment adaptability distribution map is constructed, and includes: According to the internal disturbance mean, boundary heterogeneity index, area size and spatial form compactness of each local ecological patch unit, a climate sensitivity weight factor of the patch is calculated; the weighted patch layer is superimposed on a microclimate variable high-frequency disturbance response map, and a Kriging spatial interpolation method is used to generate a micro-environment adaptability distribution map.

5. The method of claim 4, wherein the method further comprises: The habitat connection dynamic evolution map is derived based on the linkage strength distribution, and includes: ​ Based on the micro-environment adaptability distribution map, the spatial centroid coordinates of all local ecological patch units are extracted, and the patch climate sensitivity weight factor is used as a node attribute to construct a regional ecological node set, and a threshold adjacency method is used to preliminarily identify adjacent patch nodes; For each adjacent patch node, the spatial distance, patch sensitivity difference, interval region disturbance index mean and transition zone width are comprehensively considered, and an exponential decay function is used to calculate the ecological connectivity probability value; All patch nodes and the connectivity probability relationship between them are represented as a weighted undirected graph structure to generate an ecological connectivity probability set map, and a beta-connection coefficient is introduced to calculate the connection strength and network centrality of each node in the graph structure; Through time-series disturbance data-driven connectivity change simulation, a sliding window analysis method is used to extract the connectivity probability change trend, and a habitat connection dynamic evolution map is further constructed, with time as the horizontal axis and linkage strength as the vertical axis, to describe the connectivity evolution trajectory of the ecological system under the disturbance background. 6.The method of claim 1, wherein the method further comprises: determining a change in the biodiversity index; and determining a change in the biodiversity index based on the change in the biodiversity index. The diversity attenuation trend prediction curve for a local region is generated, and includes: Species response sequences are selected from a biological response coupling path library, including species distribution changes, behavior migration, reproduction interruption and population density change characteristics, and sample screening is performed according to response intensity, frequency and spatial coverage; After time alignment and variable normalization of the species response sequences, a t-SNE or UMAP nonlinear dimension reduction algorithm is used for low-dimensional embedding, and high-dimensional ecological response behaviors are mapped to two-dimensional or three-dimensional space to retain the local and global structural features of the response trajectory; A weighted graph structure is constructed based on the similarity of response sequences in a low-dimensional embedding space, and a spectral clustering algorithm is used to divide the nodes in the graph structure into modules to identify species response clusters with common ecological response mechanisms; For each response cluster, a diversity decay model is established, a local regional diversity change fitting function is constructed based on historical response trends and disturbance driving factors, and a diversity decay trend prediction curve is output.

7. The method of claim 6, wherein the method further comprises: determining the impact of climate change on biodiversity based on the determined relationship between the climate change and the biodiversity. The prediction curve is fitted with the existing biological survey database to calculate the dynamic deviation factor, including: Obtaining the historical biological survey database corresponding to the target ecological region, including the species occurrence frequency, population density and habitat utilization rate index at different periods, and arranging them into a structured data set according to the time sequence; The predicted diversity decay trend curve is fitted with the actual observation data in the biological survey database, the prediction error at different time points is calculated using a dynamic residual evaluation method, and a dynamic deviation factor time series is constructed.

8. The method of claim 7, wherein the method further comprises: determining a relationship between the climate change and the biological diversity. The potential biological diversity collapse critical point is identified, including: Mapping the optimized biological diversity impact evaluation model in the target ecological region, performing raster-by-raster operation on each ecological unit in the region, and outputting a biological diversity vulnerability index raster layer on a continuous space; Based on the vulnerability index layer, a biological diversity risk heat map is generated, the kernel density estimation method is used to aggregate the high vulnerability value regions in space and express the heat gradient, and the risk level regions are set according to the heat threshold; Based on the heat map and combined with the dynamic trend prediction curve, the spatial area with a significant continuous downward trend in the vulnerable region is identified, the diversity decay slope and fluctuation stability index of the spatial area are calculated, and the potential biological diversity collapse critical point is identified and marked.

Images