A population stability evaluation method and device based on Lyapunov function and medium

Abstract

The application discloses a kind of population stability evaluation method, equipment and medium based on Lyapunov function, it is related to modeling simulation technical field, including, the multiple-source record of target population is collected, pre-processing is executed, and time-varying growth rate candidate sequence, constant survival rate candidate parameter and various phase identifiers are extracted, to obtain phase candidate population basic data set;Single-step model simulation and multi-step model simulation are carried out on phase condition Lyapunov atlas, to obtain population state evolution result, identify stable core area, critical hinge area, escape boundary and return path, form attractor domain boundary chain verification set.The application realizes the unified expression of local contraction trend, cross-region transition characteristics and boundary continuous relationship in population state space by constructing phase condition Lyapunov atlas.

Description

Technical Field

[0001] This invention relates to the field of modeling and simulation technology, and in particular to a method, device and medium for evaluating population stability based on Lyapunov functions. Background Technology

[0002] In recent years, discrete population dynamics, Lyapunov stability analysis, and multi-source time series data modeling have continued to develop. Parameter identification, state representation, and evolutionary simulation based on population size, survival capacity, resource supply, migration activities, and external disturbance records have become important technical approaches in population system research. At the same time, computer-aided modeling, verification, and simulation methods for complex objects have matured, providing a new foundation for population state space characterization, local behavior analysis, and stability assessment.

[0003] While existing technologies can perform parameter fitting or stability assessment of population evolution processes, most methods still rely on a single equilibrium point, a single criterion, or holistic analysis. They are unable to uniformly model and continuously verify the local evolution, inter-regional transition, and stability boundary of the population state under the combined effects of time-varying growth rate, constant survival rate, and multi-phase perturbation. As a result, it is difficult to form a stability assessment result that combines boundary margin representation, instability trend identification, and early warning suggestion output capabilities. Summary of the Invention

[0004] In view of the aforementioned existing problems, the present invention is proposed.

[0005] Therefore, this invention provides a population stability assessment method based on Lyapunov functions to solve the problem of difficulty in modeling, simulating, and verifying the stability boundary of time-varying populations under complex perturbations.

[0006] To solve the above-mentioned technical problems, the present invention provides the following technical solution: In a first aspect, this invention provides a population stability assessment method based on Lyapunov functions, comprising: collecting multi-source records of a target population, performing preprocessing, and extracting candidate sequences of time-varying growth rates, candidate parameters of constant survival rates, and various phase identifiers to obtain a phase candidate population base dataset; writing the phase candidate population base dataset into population state entries at each time step, and performing equilibrium deviation normalization and phase binding to form a phase-conditional population state sequence; performing population evolution relationship rolling fitting and local partitioning of the state space on the phase-conditional population state sequence to determine the reference center of each local partition, Adjacency relationships and transition directions are identified. Local Lyapunov evaluation entries are constructed and assembled around the reference center, adjacency relationships, and transition directions of each local partition to form a phase-conditional Lyapunov map. Single-step and multi-step model simulations are performed on the phase-conditional Lyapunov map to obtain population state evolution results, identifying stable core regions, critical hinge regions, escape boundaries, and regression paths to form an attraction domain boundary chain validation set. Based on the attraction domain boundary chain validation set, the current population's stability level, boundary margin, instability trend, and key driving phases are calculated, generating early warning results and control suggestions.

[0007] As a preferred embodiment of the population stability assessment method based on the Lyapunov function described in this invention, the steps for obtaining the basic dataset of the phase candidate population are as follows: The multi-source records include population size records, viability records, resource availability records, climate fluctuation records, migration activity records, and historical intervention records for the current assessment period and historical assessment periods. The preprocessing includes time unification, missing data completion, anomaly removal, dimension normalization, and confidence labeling. Based on the preprocessed multi-source records, according to the non-autonomous discrete population evolution relationship characterized by population size change, time-varying growth rate and constant survival rate at adjacent time points, candidate parameters of constant survival rate are extracted, and candidate sequences of time-varying growth rate are extracted by combining population size records and candidate parameters of constant survival rate. By combining the candidate sequences of time-varying growth rates, candidate parameters of constant survival rates, and preprocessed multi-source records, resource, perturbation, migration, and reproduction phase identification is determined, and a basic dataset of phase candidate populations is obtained.

[0008] As a preferred embodiment of the population stability assessment method based on the Lyapunov function described in this invention, the steps for forming the phase-conditional population state sequence are as follows: The data entries at each time step in the phase candidate population base dataset are written into the corresponding initial population state entries in time index order to obtain the initial population state sequence. Based on the degree of deviation and local fluctuation amplitude of the population size value in the initial population state sequence relative to the unique positive equilibrium point corresponding to the non-autonomous discrete population evolution relationship, the equilibrium deviation normalization value at each time moment is calculated to generate the equilibrium deviation state sequence. The equilibrium deviation normalization value, equilibrium deviation change, time-varying growth rate candidate value, constant survival rate candidate parameter, and various phase identifiers in the equilibrium deviation state sequence are bound to each other and arranged in time index order to obtain the phase conditional population state sequence.

[0009] As a preferred embodiment of the population stability assessment method based on the Lyapunov function described in this invention, the steps for determining the reference center, adjacency relationship, and transition direction of each local partition are as follows: Confidence-weighted rolling fitting and aggregate partitioning are performed on the phase-conditional population state sequence to determine the reference center of each local partition; Based on the reference center of each local partition, the state transition between different local partitions is statistically analyzed, and the adjacency relationship and transition direction between each local partition are determined.

[0010] As a preferred embodiment of the population stability assessment method based on the Lyapunov function described in this invention, the steps for forming the phase-conditional Lyapunov map are as follows: Around the reference center of each local partition, the deviation trajectory, advancement trajectory and cross-regional transition trajectory of the phase condition population state entries in the corresponding local partition are extracted to obtain the directional interface record; Based on the direction interface records, determine the degree of change of each local partition shrinking towards the reference center, the degree of change of advancement along the jump direction, and the smoothness of the boundary transition, and write the corresponding local Lyapunov evaluation entries to obtain the sequence of local Lyapunov evaluation entries to be spliced. Based on the sequence of local Lyapunov evaluation entries to be spliced, the shared boundary consistency check is performed on adjacent local partitions with adjacency relationship. When the difference between local Lyapunov evaluation entries on both sides of the shared boundary exceeds the preset boundary consistency threshold, the corresponding cross-region transition trajectory, advancement trajectory and boundary correspondence are backtracked and adjusted to obtain a sequence of local Lyapunov evaluation entries with consistent boundaries. The local Lyapunov evaluation entries along the boundary are sequentially spliced ​​together according to the stability of the adjacency relationship and the degree of conflict of the transition direction, forming the phase condition Lyapunov map.

[0011] As a preferred embodiment of the population stability assessment method based on the Lyapunov function described in this invention, the steps for obtaining the population state evolution results are as follows: The simulation starting point is set around the reference center of each local partition, the location of the shared boundary, and the location of the main transition direction in the phase condition Lyapunov plot; Using the simulation starting point as input, single-step state progression is performed along the phase condition Lyapunov map, and the single-step landing point position, local partition changes, and local Lyapunov evaluation value changes are recorded. Using the single-step landing point location, local partition changes, and local Lyapunov evaluation value changes as inputs, multi-step state advancement is performed, recording the continuous landing point trajectory, cross-region traversal order, and local Lyapunov evaluation value evolution order to form the population state evolution results.

[0012] As a preferred embodiment of the population stability assessment method based on the Lyapunov function described in this invention, the steps for forming the attraction domain boundary chain verification set are as follows: Based on the population state evolution results, regions where continuous landing trajectories remain in the same local partition or stably adjacent local partitions and where the local Lyapunov evaluation value continuously decreases are identified as stable core regions. Based on the population state evolution results outside the stable core region, the region with continuous landing trajectory switching and alternating local Lyapunov evaluation value is identified as the critical articulation region, and the shared boundary of repeatedly crossing out of the current local partition and its stable adjacent local partition is identified as the escape boundary. Based on the population state evolution results corresponding to the critical hinge region and the escape boundary, the cross-regional trajectory for re-entering the stable core region is extracted as the regression path. The stable core region, critical hinge region, escape boundary and regression path are associated and written into the corresponding records to form the attraction domain boundary chain validation set.

[0013] As a preferred embodiment of the population stability assessment method based on the Lyapunov function described in this invention, the steps for generating early warning results and control suggestions are as follows: The stability of the core region, the distribution of the critical articulation region, the escape boundary triggering, the recovery of the regression path, and the frequency of the corresponding phases in the validation set of the attraction domain boundary chain are statistically analyzed, and the stability level, boundary margin, and instability trend of the current population are calculated. The system jointly determines the stability level, boundary margin, and instability trend. It also identifies key driving phases by combining the distribution of critical articulation regions, escape boundary triggering conditions, and the frequency of occurrence of various phases, and generates early warning results and control recommendations.

[0014] In a second aspect, the present invention provides a computer device, including a memory and a processor, wherein the memory stores a computer program, wherein: when the computer program is executed by the processor, it implements any step of the population stability assessment method based on the Lyapunov function as described in the first aspect of the present invention.

[0015] Thirdly, the present invention provides a computer-readable storage medium having a computer program stored thereon, wherein: when the computer program is executed by a processor, it implements any step of the population stability assessment method based on the Lyapunov function as described in the first aspect of the present invention.

[0016] The beneficial effects of this invention are as follows: by constructing a phase-conditional Lyapunov map, a unified expression of local contraction trends, cross-regional transition characteristics, and boundary continuity relationships in the population state space is achieved; single-step model simulation and multi-step model simulation can be carried out continuously along the same map framework, and further support fine-grained identification of stable core regions, critical hinge regions, escape boundaries, and regression paths. Attached Figure Description

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

[0018] Figure 1 This is a flowchart of a population stability assessment method based on the Lyapunov function.

[0019] Figure 2 A flowchart for obtaining the phase-conditional population state sequence.

[0020] Figure 3 A flowchart for generating phase-conditional Lyapunov plots.

[0021] Figure 4 A flowchart for generating early warning results and control recommendations.

[0022] Figure 5 A plot of numerical evolution data for local Lyapunov evaluation entries.

[0023] Figure 6 This is a diagram of local partitions of the state space and the reference center data. Detailed Implementation

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

[0025] Many specific details are set forth in the following description in order to provide a full understanding of the invention. However, the invention may also be practiced in other ways different from those described herein, and those skilled in the art can make similar extensions without departing from the spirit of the invention. Therefore, the invention is not limited to the specific embodiments disclosed below.

[0026] Secondly, the term "one embodiment" or "embodiment" as used herein refers to a specific feature, structure, or characteristic that may be included in at least one implementation of the present invention. The phrase "in one embodiment" appearing in different places in this specification does not necessarily refer to the same embodiment, nor is it a single or selective embodiment that is mutually exclusive with other embodiments.

[0027] Reference Figures 1-6 This is one embodiment of the present invention, which provides a population stability assessment method based on the Lyapunov function, including the following steps: S1. Collect multi-source records of the target population, perform preprocessing, and extract candidate sequences of time-varying growth rate, candidate parameters of constant survival rate, and various phase identifiers to obtain the basic dataset of phase candidate populations.

[0028] S1.1: Multi-source records include population size records, viability records, resource availability records, climate fluctuation records, migration activity records, and historical intervention records for the current assessment period and historical assessment periods; It should be noted that the population size record comes from the population monitoring, survey statistics or automatic counting results of the target population in the current assessment period and historical assessment periods. Its purpose is to characterize the population size change of the target population in adjacent time periods and to serve as the basis for extracting candidate sequences of time-varying growth rate and state evolution analysis. Survival records are derived from the survival rate statistics, retention observations, or life cycle tracking results of the target population in the current assessment period and historical assessment periods. Their purpose is to extract candidate parameters for constant survival rate and characterize the continuous survival level of the target population. Resource availability records are derived from monitoring results of the supply of food, space, water, nutrients or other survival resources in the living environment of the target population. Their function is to characterize changes in resource conditions and support the determination of resource phase markers. Climate fluctuation records are derived from continuous observations of temperature, rainfall, humidity, wind speed, light, or other meteorological factors in the environment in which the target population is located. Their function is to characterize the intensity of external environmental disturbances and support the determination of disturbance phase markers. Migration activity records are derived from the movement observation, migration statistics, trajectory collection, or passage records of the target population between different regions. Their purpose is to characterize the spatial flow status of the target population and support the determination of migration phase markers. Historical intervention records are derived from records of feeding, supplementation, isolation, migration guidance, protection control, or other management measures implemented on the target population during the current assessment period and historical assessment periods. Their purpose is to characterize the impact of external regulation on population evolution and to support the reproductive and perturbation phases.

[0029] S1.2: Preprocessing includes time unification, missing data completion, anomaly removal, dimension normalization, and confidence labeling; Specifically, for multi-source records of the target population, the timestamps of various records are converted into a unified time index format, and records from different sources are aligned to the same time step based on the unified time index. For missing data in the aligned multi-source records, interpolation is performed using effective observations at adjacent times to fill in the gaps. Based on the preset numerical range thresholds corresponding to different types of multi-source records, outlier values ​​that exceed normal physiological limits, environmental physical limits, monitoring equipment ranges, or recording logic constraints are identified and removed to achieve anomaly removal. The range transformation method is used to map the values ​​of various records after anomaly removal to a unified dimensionless interval to eliminate the influence of different physical dimensions and achieve dimension normalization. Based on the measurement error percentage of the monitoring equipment corresponding to each record, each normalized multi-source record is assigned a corresponding weight value, and multiple weight values ​​participating in the formation of the comprehensive state entry under the same time index are aggregated to complete the confidence labeling, resulting in the preprocessed multi-source records.

[0030] It should be noted that the numerical range thresholds are set based on the biological characteristics limits of the target population, the physical limits of its environment, the range of monitoring equipment, and the logical constraints of the recorded fields. Specifically, the temperature values ​​in the climate fluctuation records are exemplarily set between -50 degrees Celsius and 50 degrees Celsius; the population size records are set to a range of not less than 0 and not exceeding the historical maximum number of the target population or the carrying capacity limit of its habitat; the survival capacity records are set to a range of 0 to 1 or 0% to 100%; the resource availability records are set to a range of not less than 0 and not exceeding the range limit of the corresponding resource monitoring equipment or the historical resource supply limit; the migration activity records are set to a range of not less than 0 and not exceeding the maximum count value of the channel monitoring per unit time; and the value ranges of the historical intervention records are set according to the intervention type field, intervention duration field, and intervention intensity field, with the intervention duration not exceeding the length of the current assessment period and the intervention intensity not exceeding the maximum implementation intensity of the corresponding management measures.

[0031] For multi-source records after dimension normalization, the weight value is determined based on the percentage of measurement error of the monitoring equipment corresponding to each record. When the percentage of measurement error is not higher than 5%, a weight value between 0.8 and 1.0 is assigned; when the percentage of measurement error is higher than 5% but not higher than 10%, a weight value between 0.6 and 0.8 is assigned; when the percentage of measurement error is higher than 10% but not higher than 20%, a weight value between 0.4 and 0.6 is assigned; when the percentage of measurement error is higher than 20%, the corresponding record is marked as a low-confidence record and assigned a weight value of 0.4. The specific weight value is determined by interpolation within the linear interval of the percentage of measurement error; the lower the percentage of measurement error, the higher the corresponding weight value. When multiple records under the same time index jointly form a comprehensive status entry, the arithmetic mean of the weight values ​​of the participating records is taken as the confidence label value of the comprehensive status entry. For records obtained by missing data completion, the average of the weight values ​​of adjacent valid observation records is taken and multiplied by the completion reduction coefficient as the corresponding weight value. The example value of the completion reduction coefficient is 0.9.

[0032] S1.3: Based on the preprocessed multi-source records, according to the non-autonomous discrete population evolution relationship characterized by population size change, time-varying growth rate and constant survival rate at adjacent time points, candidate parameters of constant survival rate are extracted, and candidate sequences of time-varying growth rate are extracted by combining population size records and candidate parameters of constant survival rate. Specifically, the preprocessed survival records are read time by time, and all survival records in the current assessment period and the historical assessment period are arranged according to the unified time index to form a survival sample sequence. The constant survival rate is used to characterize the relatively stable survival level of the target population during the evolution of adjacent time periods. Therefore, the overall representative value is extracted from the survival sample sequence as a candidate parameter for the constant survival rate. After extracting the candidate parameters for constant survival rate, the preprocessed population size records are read time-by-time, and paired according to the unified time index to form population size change pairs for adjacent time moments. For each pair of population size change pairs for adjacent time moments, the ratio of the population size change at the next time moment to the population size at the previous time moment is calculated. The candidate parameters for constant survival rate are substituted into the non-autonomous discrete population evolution relationship, and the time-varying growth rate at the corresponding time moment is inverted and solved to obtain the candidate values ​​for the time-varying growth rate at each time moment. All the candidate values ​​for the time-varying growth rate are arranged in order of time index to form the candidate sequence of time-varying growth rate.

[0033] It should be noted that the non-autonomous discrete population evolution relationship refers to the quantitative relationship used to describe the discrete evolution law of the target population over time, characterized by the population size change, time-varying growth rate and constant survival rate at adjacent time points.

[0034] The calculation expressions for the constant survival rate candidate parameter and the time-varying growth rate candidate sequence are as follows: The inversion of candidate values ​​for the time-varying growth rate is based on the positive non-autonomous discrete population evolution relationship:

[0035] In the formula, Candidate parameters for constant survival rate; Indicates median value operation; This shows a sequence of normalized survival values ​​arranged by time index, where... Indicates the first The normalized value of survival ability at a given moment; This represents the total number of time points after combining the current assessment period and historical assessment periods; Indicates the first The moment to the The ratio of population size change at each time point; Indicates the first The normalized value of the population size at each moment; Indicates the first The normalized value of the population size at each moment; Indicates the first Candidate values ​​for the time-varying growth rate at each time point.

[0036] It should be noted that the above calculations apply to The adjacent time data; when the adjacent time data does not meet the above conditions, the corresponding adjacent time data pair is removed, or the adjacent valid time data pair is used to invert the candidate value of the time-varying growth rate.

[0037] S1.4: The candidate sequences of time-varying growth rates, candidate parameters of constant survival rates, and preprocessed multi-source records are used to determine the resource, perturbation, migration, and reproduction phase identifiers to obtain the basic dataset of phase candidate populations.

[0038] Specifically, for the preprocessed multi-source records, the moments when the food supply value in the resource availability record is higher than the preset resource threshold are marked as resource phase identifiers; the moments when the temperature change amplitude in the climate fluctuation record is higher than the preset perturbation threshold are marked as perturbation phase identifiers; the moments when the cross-regional migration frequency value in the migration activity record changes abruptly are marked as migration phase identifiers; and the time periods in the historical intervention records that contain reproduction observation markers, or that contain feeding operation markers, replenishment operation markers, and have reproduction observation markers within the corresponding time period are marked as reproduction phase identifiers. The resource phase identifiers, perturbation phase identifiers, migration phase identifiers, and reproduction phase identifiers are aligned with the time-varying growth rate candidate sequence and the constant survival rate candidate parameter to obtain the phase candidate population basic dataset.

[0039] It should be noted that the resource threshold is set based on the minimum food supply required for the target population to maintain basic survival metabolism. It is selected within the range of 0 to 1 after the resource availability record is normalized. An example value range is 0.4 to 0.8. The value is based on the observed minimum food supply that the target population can maintain a stable population size without decline in the historical environment. For different species or under different environments, the resource threshold is adjusted according to the historical upper limit of resource supply, the minimum survival supply, and the stable population size range. The disturbance threshold is set based on the tolerance limit of the target population to climate and environmental changes. It is selected in the range of 0 to 1 after the temperature change amplitude is dimensionally normalized. An example value range is 0.3 to 0.7. The moment of mutation refers to the moment when the normalized difference between the cross-regional movement frequency values ​​in the migration activity record within two adjacent unified time indices exceeds the preset mutation threshold. The mutation threshold is set based on the normal movement fluctuation range of the target population in the natural environment and is selected in the range of 0 to 1 after the migration activity record is normalized. An example value range is 0.2 to 0.5.

[0040] S2. Write the phase candidate population basic dataset into the population state entries at each time step, and perform balance deviation normalization and phase binding to form a phase conditional population state sequence.

[0041] S2.1: Write the data entries of each time step in the phase candidate population basic dataset into the corresponding initial population state entries in the order of time index to obtain the initial population state sequence; Specifically, based on the unified time index contained in the phase candidate population basic dataset, the data entries at each time point in the phase candidate population basic dataset are arranged in chronological order, and the arranged data entries at each time point are filled into the time field, population size normalization value field, time-varying growth rate field, constant survival rate field, resource phase field, disturbance phase field, migration phase field, and reproduction phase field of the corresponding initial population state entry, forming an initial population state sequence distributed according to time series.

[0042] S2.2: Based on the degree of deviation and local fluctuation amplitude of the population size value in the initial population state sequence relative to the unique positive equilibrium point corresponding to the non-autonomous discrete population evolution relationship, calculate the equilibrium deviation normalization value at each time point and generate the equilibrium deviation state sequence. Specifically, the population size normalization value in the initial population state sequence is read time by time, and the unique positive equilibrium point corresponding to the non-autonomous discrete population evolution relationship is used as a unified reference center. The deviation of the population size normalization value at each time time from the unique positive equilibrium point is calculated. A local time window is constructed around the current time, and the median scale of the squared deviation value within the local time window is used as the local fluctuation amplitude. The equilibrium deviation normalization value at each time time is calculated, and the difference calculation is performed on the equilibrium deviation normalization values ​​of two adjacent times to obtain the change in equilibrium deviation, thus forming the equilibrium deviation state sequence.

[0043] The calculation expressions for the unique positive equilibrium point, local fluctuation amplitude, equilibrium deviation normalization value, and equilibrium deviation change are as follows: Among them, the unique positive equilibrium point is based on the positive non-autonomous discrete population evolution relationship, in and Determined under equilibrium conditions;

[0044] In the formula, It represents the unique positive equilibrium point corresponding to the evolutionary relationship of a non-autonomous discrete population; Indicates the first Each moment corresponds to the local fluctuation amplitude within a local time window; Indicates the first time within the local time window Normalized population size for a local time window; Indicates the local time window index; Indicates the length of the local time window; Indicates the first The equilibrium deviation from the normalized value at each moment; This indicates a very small number that prevents the denominator from being zero; Indicates the first The change in equilibrium deviation at each moment; Indicates the first The equilibrium deviation from the normalized value at each moment.

[0045] It should be noted that, Based on the data precision setting of the population size normalization value, an exemplary value range is provided. ; The value is set based on the sampling frequency of the multi-source records and the current evaluation period length, with an exemplary range of 3 to 10 time steps.

[0046] S2.3: Bind the equilibrium deviation normalization value, equilibrium deviation change, time-varying growth rate candidate value, constant survival rate candidate parameter and various phase identifiers in the equilibrium deviation state sequence, and arrange them in time index order to obtain the phase conditional population state sequence.

[0047] Specifically, based on the unified time index contained in the equilibrium deviation state sequence and the initial population state sequence, the equilibrium deviation normalization value and equilibrium deviation change at each time point in the equilibrium deviation state sequence are matched with the time-varying growth rate candidate value, constant survival rate candidate parameter, resource phase identifier, disturbance phase identifier, migration phase identifier, and reproduction phase identifier corresponding to the same time index in the initial population state sequence. The matched equilibrium deviation normalization value, equilibrium deviation change, time-varying growth rate candidate value, constant survival rate candidate parameter, resource phase identifier, disturbance phase identifier, migration phase identifier, and reproduction phase identifier are combined into a single comprehensive state entry. All comprehensive state entries are arranged in the order of their time indices to obtain the phase-conditional population state sequence.

[0048] S3. Perform population evolution relationship rolling fitting and state space local partitioning on the phase condition population state sequence to determine the reference center, adjacency relationship and transition direction of each local partition.

[0049] S3.1: Perform confidence-weighted rolling fitting and aggregate partitioning on the phase-conditional population state sequence to determine the reference center of each local partition; Specifically, the comprehensive state entries in the phase condition population state sequence are read time-by-time, and the weight values ​​corresponding to the confidence labeling results are used as the confidence basis for each comprehensive state entry. A sliding time window is constructed around the current time. Within each sliding time window, a confidence-weighted rolling fit is performed using a linear fitting method with the time index as the independent variable and the equilibrium deviation normalization value as the dependent variable, where the weight values ​​corresponding to the confidence labeling results are used as the fitting weights. The absolute value of the difference between the fitting slopes of adjacent sliding time windows is used as the change difference index. When the change difference index exceeds the preset state mutation threshold, the corresponding time index is marked as a state mutation point. The phase condition population state sequence is aggregated and partitioned according to the identified state mutation points, and the comprehensive state entries continuously distributed between state mutation points are divided into the same local partition. After completing the aggregation and partitioning, the median value of the equilibrium deviation normalization value and the equilibrium deviation change amount of the comprehensive state entries in each local partition are taken to determine the reference center of the corresponding local partition.

[0050] It should be noted that the state change threshold is set based on the distribution of the absolute value of the difference in the fitting slope of adjacent sliding time windows within the historical stability evaluation period. The exemplary value range is 0.1 to 0.3, and the exemplary value is 0.15.

[0051] S3.2: Based on the reference center of each local partition, statistically analyze the state transition between different local partitions, and determine the adjacency relationship and transition direction between each local partition.

[0052] Specifically, based on the order of the reference centers of each local partition in the time index, the frequency of transitions of the comprehensive state entries in the phase condition population state sequence from the current local partition to the adjacent local partition is counted. Two local partitions with a transition frequency greater than zero are determined to be adjacent, and the flow direction corresponding to the transition frequency is determined as the transition direction between each local partition.

[0053] Figure 6 The horizontal axis represents the equilibrium deviation normalization value, and the vertical axis represents the change in equilibrium deviation. Different colored scatter points represent different local partitions, and the yellow cross marks represent the reference centers of each local partition. As can be seen in the figure, the state points in the phase conditional population state sequence are not uniformly distributed, but form multiple local clusters within the state space composed of the equilibrium deviation normalization value and the change in equilibrium deviation. The reference center of each local partition is located at the representative position of the corresponding state point cluster and is used to calculate the degree of contraction of the state point towards the reference center, the distance to the reference center, and the cross-regional state transition relationship.

[0054] S4. Construct and stitch together local Lyapunov evaluation entries around the reference center, adjacency relationship and transition direction of each local partition to form a phase condition Lyapunov map.

[0055] S4.1: Around the reference center of each local partition, extract the deviation trajectory, advancement trajectory and cross-regional transition trajectory of the phase condition population state entry in the corresponding local partition to obtain the direction interface record; Specifically, based on the reference center of each local partition, the distance between each integrated state entry in the corresponding local partition and the reference center is calculated, and the distance changes are read in the order of the time index. Continuous integrated state entries whose distance from the reference center gradually decreases after entering the corresponding local partition are classified as the entry stage; continuous integrated state entries whose distance from the reference center remains within a preset center neighborhood threshold are classified as the stable stage; and continuous integrated state entries whose distance from the reference center gradually increases and moves closer to the shared boundary of adjacent local partitions are classified as the exit stage. The path of the change of the balance deviation normalization value with the time index in the entry stage is extracted as the deviation trajectory; the path of the fluctuation of the balance deviation normalization value around the reference center in the stable stage is extracted as the advancement trajectory; and the path of the transition of the integrated state entries to adjacent local partitions in the exit stage is extracted as the cross-regional transition trajectory. The deviation trajectory, advancement trajectory, and cross-regional transition trajectory are combined in the order of the time index to form a direction interface record.

[0056] It should be noted that the center neighborhood threshold is set based on the historical distribution of the distance from the comprehensive state entry to the reference center within the corresponding local partition, with an exemplary value range of 0.05 to 0.20.

[0057] S4.2: Based on the direction interface record, determine the degree of change of each local partition shrinking towards the reference center, the degree of change of advancement along the transition direction, and the degree of boundary transition smoothness, and write the corresponding local Lyapunov evaluation entries to obtain the sequence of local Lyapunov evaluation entries to be spliced. Specifically, the deviation trajectory, advancement trajectory, and cross-regional transition trajectory in the direction interface records are read one by one; for the deviation trajectory, the ratio of the distance change from the integrated state entry to the reference center at adjacent time moments is calculated as the value of the degree of change of the contraction towards the reference center. When the degree of change value is greater than 0, it indicates that the contraction is enhanced, and when it is less than 0, it indicates that the contraction is weakened. For the propulsion trajectory, the change in state coordinates is determined based on the normalized value of balance deviation and the change in balance deviation in the comprehensive state entries at adjacent time points. The same-direction change ratio is calculated in combination with the jump direction as a numerical value of the degree of change in propulsion along the jump direction. When the same-direction change ratio is greater than 0, it indicates that the propulsion is enhanced, and when it is less than 0, it indicates that the propulsion is weakened. For cross-regional transition trajectories, the Euclidean distance between the final integrated state entry of the exit stage and the initial integrated state entry of the adjacent local partition in terms of the equilibrium deviation normalization value and the equilibrium deviation change amount is calculated as the change difference value. When the change difference value is not higher than the preset transition continuity threshold, the boundary transition state is determined to be continuous, and the difference between 1 and the ratio of the change difference value to the preset transition continuity threshold is recorded as the boundary transition smoothness value. When the change difference value is higher than the preset transition continuity threshold, the boundary transition state is determined to be discontinuous, and the boundary transition smoothness value is recorded as 0. The numerical values ​​of the degree of contraction towards the reference center, the degree of change along the transition direction, and the degree of boundary transition smoothness, along with the equilibrium deviation normalization value and the change in equilibrium deviation at the corresponding time, are written into the local Lyapunov evaluation entries and arranged in time index order to obtain the sequence of local Lyapunov evaluation entries to be spliced.

[0058] It should be noted that the transition continuity threshold is set based on the distribution of differences in the changes of state entries in the normal transition phase of adjacent local partitions in the historical phase condition population state sequence, with an exemplary value range of 0.1 to 0.3.

[0059] S4.3: Based on the sequence of local Lyapunov evaluation entries to be spliced, perform a shared boundary consistency check on adjacent local partitions with adjacency relationships. When the difference between local Lyapunov evaluation entries on both sides of the shared boundary exceeds the preset boundary consistency threshold, perform a backtracking adjustment on the corresponding cross-region transition trajectory, advancement trajectory and boundary correspondence to obtain a sequence of local Lyapunov evaluation entries with consistent boundaries. Specifically, based on adjacent local partitions with adjacency relationships in the sequence of local Lyapunov evaluation entries to be spliced, the corresponding values ​​of the degree of contraction towards the reference center, the degree of advancement along the transition direction, and the degree of boundary transition smoothness in the local Lyapunov evaluation entries on both sides of the shared boundary are read, and the average of the differences of the three values ​​is calculated as the difference of local Lyapunov evaluation entries on both sides of the shared boundary. When the difference of local Lyapunov evaluation entries exceeds the preset boundary consistency threshold, the corresponding position of the candidate boundary is traced back according to the time index in the original direction interface record, and the cross-regional transition trajectory, advancement trajectory, and boundary correspondence are rematched to ensure that the difference between the final integrated state entry of the exit stage and the initial integrated state entry of the adjacent local partition is not higher than the preset transition continuity threshold. If the preset transition continuity threshold is still not met after several consecutive backtrackings (e.g., 3 times), the corresponding shared boundary is marked as the boundary to be checked, and the original cross-regional transition trajectory and advancement trajectory are retained. The rematched cross-regional transition trajectory, advancement trajectory, and boundary correspondence are written into the local Lyapunov evaluation entries to obtain a sequence of local Lyapunov evaluation entries with consistent boundaries.

[0060] It should be noted that the boundary consistency threshold is set based on the difference distribution of local Lyapunov evaluation entries on both sides of the shared boundary of adjacent local partitions under historical stable conditions, with an exemplary value range of 0.10 to 0.25; the rollback adjustment is only used to rematch the corresponding position of the boundary in the original direction interface record, and does not change the original integrated state entries.

[0061] S4.4: The local Lyapunov evaluation entries along the boundary are sequentially spliced ​​together according to the stability of the adjacency relationship and the degree of conflict of the transition direction, forming the phase condition Lyapunov map.

[0062] Specifically, based on the consistent local Lyapunov evaluation item sequence, the variances of the degree of contraction towards the reference center and the degree of advancement along the transition direction of adjacent local partitions at the shared boundary are calculated. Adjacent local partitions with variance values ​​lower than a preset adjacency stability threshold are classified as having high adjacency stability. For adjacent local partitions with high adjacency stability, the tangential direction of the shared boundary is determined based on the two endpoints of the shared boundary, and the direction perpendicular to the tangential direction of the shared boundary and pointing towards the reference center of the adjacent local partition is determined as the boundary normal direction. The angle between the transition direction and the boundary normal direction is calculated, and adjacent local partitions with an angle less than a preset direction conflict threshold are classified as having low transition direction conflict. According to the priority order of high adjacency stability and low transition direction conflict, the local Lyapunov evaluation items of adjacent local partitions are numerically fused at the shared boundary to form a phase conditional Lyapunov map.

[0063] For example, the adjacency stability threshold is set based on the variance distribution of the degree of contraction and advancement changes at the boundary of historical stable local partitions, with an exemplary value range of 0.05 to 0.20; the directional conflict threshold is set based on the angle distribution between the transition direction and the boundary normal direction during historical normal transitions, with an exemplary value range of 15 degrees to 45 degrees.

[0064] S5. Perform single-step and multi-step model simulations on the phase-conditional Lyapunov map to obtain the population state evolution results, identify the stable core region, critical hinge region, escape boundary and regression path, and form an attraction domain boundary chain validation set.

[0065] S5.1: Set the simulation starting point around the reference center of each local partition, the location of the shared boundary, and the location of the main transition direction in the phase condition Lyapunov plot; Specifically, based on the reference center coordinates of each local partition in the phase-condition Lyapunov map, a first type of simulation starting point is generated within the reference center's neighborhood radius, where the reference center's neighborhood radius is set according to the median value of the distance from the synthetic state entry to the reference center within the corresponding local partition. Based on the shared boundary position between adjacent local partitions in the phase-condition Lyapunov map, the mean value of the coordinates of the two endpoints of the shared boundary is used as the center point of the shared boundary, and a second type of simulation starting point is generated at the position offset from the boundary spacing along the main transition direction from the center point of the shared boundary, where the boundary spacing is set according to a preset ratio of the distance between the reference centers of adjacent local partitions. Based on the main transition direction of each local partition in the phase-condition Lyapunov map, the coordinates of the first and last synthetic state entries arranged along the main transition direction within the corresponding local partition are used as the starting and ending endpoints of the main transition direction, respectively, and a third type of simulation starting point is generated. The first, second, and third types of simulation starting points are merged into a set of simulation starting points.

[0066] It should be noted that the coordinates of the comprehensive status entries consist of the equilibrium deviation normalization value and the equilibrium deviation change; the preset ratio is set based on the local partition scale within the historical stability assessment period, and the exemplary value range is 5% to 20%.

[0067] S5.2: Using the simulation starting point as input, perform single-step state advancement along the phase condition Lyapunov map, and record the single-step landing point position, local partition changes, and local Lyapunov evaluation value changes; Specifically, based on the position of the simulation starting point in the phase condition Lyapunov map from the set of simulation starting points, the corresponding local Lyapunov evaluation entry is read, and the population state coordinates at the next moment are calculated using the non-autonomous discrete population evolution relationship as the single-step landing point position. The single-step landing point position is compared with the local partition to which the simulation starting point belongs to determine the local partition change. The degree of change towards the reference center in the local Lyapunov evaluation entry corresponding to the single-step landing point position is read and compared with the degree of change towards the reference center corresponding to the simulation starting point to obtain the change of the local Lyapunov evaluation value.

[0068] It should be noted that the change in local partitioning refers to the change in the relationship between the single-step landing point and the simulation starting point in the region of the phase condition Lyapunov map spatial structure. Specifically: Staying within the original local partition: The single-step landing point and the simulation starting point are located within the same local partition and do not cross any shared boundary positions; Crossing shared boundaries and entering adjacent local partitions: The single-step landing point moves from the current local partition to an adjacent local partition with an adjacency relationship based on the main transition direction or adjacency relationship; Leaving the current local partition and entering an undefined region: The single-step landing point exceeds the coverage area of ​​all spliced ​​local partitions in the phase condition Lyapunov map.

[0069] The change in local Lyapunov evaluation value refers to the numerical evolution characteristics of the stability index in the local Lyapunov evaluation item corresponding to the single-step landing position relative to the index corresponding to the simulation starting point. Among them, the local Lyapunov evaluation item is used to record the degree of change of contraction towards the reference center, the degree of change of advancement along the transition direction, and the smoothness of the boundary transition. The local Lyapunov evaluation value takes the value of the change of contraction towards the reference center as the main judgment index, and the distance to the reference center is used to calculate the value of the change of contraction towards the reference center. Specifically, the numerical values ​​of the degree of contraction towards the reference center corresponding to the single-step landing position and the simulation starting point are compared. If the numerical value of the degree of contraction towards the reference center corresponding to the single-step landing position is greater than the numerical value of the degree of contraction towards the reference center corresponding to the simulation starting point, it is determined to be stable convergence; otherwise, it is determined to be divergence or fluctuation. At the same time, the numerical values ​​of the degree of advancement along the transition direction corresponding to the single-step landing position and the simulation starting point are compared to characterize the rate or intensity change of the state moving along the main transition direction.

[0070] Figure 5 The graph shows the changes in local Lyapunov evaluation item values ​​across the entire time domain and within a local magnified window. The upper part is an overview, and the lower part is the magnified interval corresponding to the red dashed box. The horizontal axis represents the time index, and the vertical axis represents the local Lyapunov evaluation item values. The blue curve represents the degree of contraction towards the reference center, and the orange curve represents the distance to the reference center. As can be seen from the graph, the orange curve fluctuates significantly overall, indicating that the distance from the population state point to the local partition reference center is more sensitive to changes in perturbation. The blue curve is generally at a low level, but shows a significant rise at the end of the local magnified window, indicating that the local contraction trend of the population state has changed during this period.

[0071] S5.3: Using the single-step landing point position, local partition changes, and local Lyapunov evaluation value changes as inputs, perform multi-step state advancement, record the continuous landing point trajectory, cross-region traversal order, and local Lyapunov evaluation value evolution order to form the population state evolution results.

[0072] Specifically, the single-step landing point is used as the propagation base point for the next moment. The state propagation along the phase-conditional Lyapunov map is repeatedly executed to generate a continuous landing point coordinate sequence as a continuous landing point trajectory. Based on the local partition identifier to which each coordinate point in the continuous landing point trajectory belongs, the sequence of local partition switching events that occur between coordinate points at adjacent moments is extracted as the cross-region traversal order. Based on the local Lyapunov evaluation entry corresponding to each coordinate point in the continuous landing point trajectory, the degree of change in contraction towards the reference center and the degree of change in propagation along the transition direction are read to form the local Lyapunov evaluation value evolution order arranged with time index. The continuous landing point trajectory, cross-region traversal order and local Lyapunov evaluation value evolution order are combined to form the population state evolution result.

[0073] S5.4: Based on the population state evolution results, identify regions where continuous landing trajectories remain in the same local partition or stably adjacent local partitions and the local Lyapunov evaluation value continues to decrease, as stable core regions; Specifically, based on the continuous landing trajectories in the population state evolution results, trajectory segments in which all landing coordinates are located within the same local partition or move back and forth only between adjacent local partitions with an adjacency relationship are selected. For the selected trajectory segments, the evolution order of the corresponding local Lyapunov evaluation values ​​is read. The local Lyapunov evaluation value is represented by the square of the distance from the continuous landing trajectory to the reference center of the corresponding local partition. The degree of change in shrinkage towards the reference center is used to confirm the decreasing direction of the local Lyapunov evaluation value. When the local Lyapunov evaluation value decreases continuously with the time index and the degree of change in shrinkage towards the reference center is greater than 0, the spatial range covered by the continuous landing trajectory that satisfies the positional constraint and the numerical descent constraint is marked as the stable core region.

[0074] S5.5: Based on the population state evolution results outside the stable core region, identify the region where the continuous landing trajectory switches back and forth and the local Lyapunov evaluation value fluctuates alternately as the critical articulation region, and take the shared boundary of the current local partition and its stable adjacent local partitions that have been continuously crossed multiple times as the escape boundary. Specifically, based on the population state evolution results, continuous landing trajectories located outside the stable core region are selected, and trajectory segments exhibiting repeated back-and-forth switching characteristics between adjacent local partitions in the cross-regional traversal sequence are identified. For the selected trajectory segments, the evolution sequence of the corresponding local Lyapunov evaluation values ​​is read, and segments with non-monotonic alternating fluctuation characteristics in the degree of contraction towards the reference center are extracted. The spatial range covered by continuous landing trajectories that satisfy the back-and-forth switching constraint and numerical fluctuation constraint is marked as the critical articulation region. Based on the continuous landing trajectories within the critical articulation region, the number of consecutive times the same continuous landing trajectory crosses the shared boundary position is counted, and the shared boundary position where the number of consecutive crossings exceeds the preset escape frequency threshold is marked as the escape boundary.

[0075] It should be noted that the escape frequency threshold is set based on the maximum frequency distribution of the normal back-and-forth fluctuations of the population state within the critical articulation region between adjacent local partitions. An exemplary value range is 2 to 5 times, and the value is based on the upper quartile boundary of the number of times the continuous landing trajectory in the critical articulation region switches back and forth between adjacent local partitions in the historical population state evolution data.

[0076] S5.6: Extract the cross-regional trajectory for re-entering the stable core region as a regression path from the population state evolution results corresponding to the critical hinge region and escape boundary, and write the stable core region, critical hinge region, escape boundary and regression path into the corresponding record to form the attraction domain boundary chain validation set.

[0077] Specifically, based on the population state evolution results, continuous landing trajectories starting from the critical hinge region or escape boundary and ending at the stable core region are selected, and these selected continuous landing trajectories are marked as regression paths. The spatial range identifiers, trajectory coordinate sequences, and local Lyapunov evaluation value evolution order of the stable core region, critical hinge region, escape boundary, and regression path are correlated and written into the same validation record entry to form an attraction domain boundary chain validation set.

[0078] S6. Calculate the stability level, boundary margin, instability trend and key driving phase of the current population based on the validation set of the attraction domain boundary chain, and generate early warning results and regulation suggestions.

[0079] S6.1: Statistically analyze the stability of the core region, the distribution of the critical articulation region, the escape boundary triggering, the recovery of the regression path, and the frequency of the corresponding phase in the validation set of the attraction domain boundary chain, and calculate the stability level, boundary margin, and instability trend of the current population. Specifically, each record in the attraction domain boundary chain validation set is expanded by time index, and the number of consecutive stays in the stable core region, the number of triggers in the critical articulation region, the number of triggers in the escape boundary, and the number of successful recoverys of the regression path at the current evaluation time are counted. The median distances from the continuous landing trajectories corresponding to the current evaluation time to the nearest escape boundary and the nearest critical articulation region are calculated. The stability index, boundary margin, and instability trend at the current evaluation time are then calculated. Based on the stability index, boundary margin, and instability trend, when the stability index is higher than the high stability threshold, the boundary margin is higher than the high buffer threshold, and instability occurs... When the stability trend is not greater than 0, the current population is classified as high stability. When the stability index is higher than the medium stability threshold but lower than the high stability threshold, and the boundary margin is higher than the medium buffer threshold but lower than the high buffer threshold, the current population is classified as medium stability. When the stability index is higher than the low stability threshold but lower than the medium stability threshold, and the boundary margin is higher than the low buffer threshold but lower than the medium buffer threshold, the current population is classified as low stability. When the stability index is lower than the low stability threshold, the boundary margin is lower than the low buffer threshold, or the instability trend is continuously higher than the preset upper limit of the instability trend, the current population is classified as unstable.

[0080] The stability index, boundary margin, and instability trend at the current assessment moment are calculated using the following expression:

[0081] In the formula, Indicates the first The stability index at any given moment; Indicates the first The number of consecutive stays in the stable core area at any given moment; Indicates the first The number of times the regression path was successfully recovered at each time point; Indicates the first The number of times the critical hinge region is triggered at any given moment; Indicates the first The number of times the escape boundary is triggered at any given moment; Indicates the first Boundary margin at each moment; Indicates the first At the [time]th moment The normalized distance from the continuous landing trajectory to the nearest escape boundary; Indicates the index of continuous landing point trajectories; Indicates the first At the [time]th moment The normalized distance from the continuous landing point trajectory to the nearest critical hinge region; Indicates the first The instability trend at each assessment moment; Indicates the first The stability index at any given moment.

[0082] It should be noted that the high stability threshold, medium stability threshold, and low stability threshold are set based on the historical stability index distribution, with exemplary values ​​of 0.70, 0.45, and 0.25, respectively; the high buffer threshold, medium buffer threshold, and low buffer threshold are set based on the historical boundary margin distribution, with exemplary values ​​of 0.50, 0.30, and 0.15, respectively; and the preset upper limit of the instability trend is set based on the historical decline in the stability index before instability, with an exemplary value of 0.05. when hour, Alternatively, use 0, or use the stability index at the end of the previous evaluation period as... , A value greater than 0 indicates a decline in the stability index and an increasing trend towards instability. ≤0 indicates that the stability index has not decreased.

[0083] S6.2: Jointly determine the stability level, boundary margin, and instability trend. Identify key driving phases by combining the distribution of critical articulated regions, escape boundary triggering conditions, and the frequency of occurrence of various phases, and generate early warning results and control suggestions.

[0084] Specifically, based on the stability level, boundary margin, and instability trend at the current assessment time, the frequency of occurrence of resource phase markers, disturbance phase markers, migration phase markers, and reproductive phase markers corresponding to each record in the attraction domain boundary chain validation set is counted within a preset statistical window. Phase markers whose frequency exceeds a preset driving frequency threshold and whose synchronous occurrence reaches a preset co-correlation threshold during the assessment period of rising instability trend are identified as key driving phases. Combining the stability level, boundary margin, instability trend, and key driving phases, early warning results are generated, including the current population stability level, spatial margin from the escape boundary, future evolutionary direction prediction, and main driving phase types. Corresponding control suggestions are generated based on the key driving phases. Specifically, when the key driving phase is a resource phase marker, suggestions for resource replenishment or habitat resource adjustment are generated; when the key driving phase is a disturbance phase marker, suggestions for disturbance source reduction or environmental buffering are generated; when the key driving phase is a migration phase marker, suggestions for migration channel guidance or habitat area connection are generated; and when the key driving phase is a reproductive phase marker, suggestions for reproductive protection or reduction of reproductive period disturbance are generated.

[0085] It should be noted that the statistical window is set according to the length of the current evaluation period, with an exemplary range of 3 to 6 consecutive evaluation periods; the driving frequency threshold is set based on the frequency distribution of various phase markers during the historical population state evolution process in the rising stage of instability trend, with an exemplary range of 2 to 6 times; the co-directional correlation threshold is set based on the number of times key driving phases and instability trends occur simultaneously during the historical rising stage of instability trend, with an exemplary range of 2 to 4 times.

[0086] This embodiment also provides a computer device applicable to the population stability assessment method based on the Lyapunov function, including: a memory and a processor; the memory is used to store computer-executable instructions, and the processor is used to execute the computer-executable instructions to implement the population stability assessment method based on the Lyapunov function as proposed in the above embodiment.

[0087] The computer device can be a terminal, comprising a processor, memory, communication interface, display screen, and input devices connected via a system bus. The processor provides computing and control capabilities. The memory includes non-volatile storage media and internal memory. The non-volatile storage media stores the operating system and computer programs. The internal memory provides an environment for the operation of the operating system and computer programs stored in the non-volatile storage media. The communication interface is used for wired or wireless communication with external terminals; wireless communication can be achieved through Wi-Fi, carrier networks, NFC (Near Field Communication), or other technologies. The display screen can be an LCD screen or an e-ink screen. The input devices can be a touch layer covering the display screen, buttons, a trackball, or a touchpad on the computer device's casing, or an external keyboard, touchpad, or mouse.

[0088] This embodiment also provides a storage medium storing a computer program that, when executed by a processor, implements the population stability assessment method based on the Lyapunov function as proposed in the above embodiments. The storage medium can be implemented by any type of volatile or non-volatile storage device or a combination thereof, such as Static Random Access Memory (SRAM), Electrically Erasable Programmable Read-Only Memory (EEPROM), Erasable Programmable Read Only Memory (EPROM), Programmable Red-Only Memory (PROM), Read-Only Memory (ROM), magnetic storage, flash memory, magnetic disk, or optical disk.

[0089] In summary, this invention achieves a unified expression of local contraction trends, cross-regional transition characteristics, and boundary continuity relationships in the population state space by constructing a phase-conditional Lyapunov map; enables single-step and multi-step model simulations to be carried out continuously along the same map framework; and further supports fine-grained identification of stable core regions, critical hinge regions, escape boundaries, and regression paths.

[0090] It should be noted that the above embodiments are only used to illustrate the technical solutions of the present invention and are not intended to limit it. Although the present invention has been described in detail with reference to preferred embodiments, those skilled in the art should understand that modifications or equivalent substitutions can be made to the technical solutions of the present invention without departing from the spirit and scope of the technical solutions of the present invention, and all such modifications or substitutions should be covered within the scope of the claims of the present invention.

Claims

1. A population stability assessment method based on Lyapunov functions, characterized in that, include: Multi-source records of the target population are collected, preprocessed, and candidate sequences of time-varying growth rate, candidate parameters of constant survival rate, and various phase identifiers are extracted to obtain the basic dataset of phase candidate populations. The phase candidate population base dataset is written into the population state entries at each time step, and balance deviation normalization and phase binding are performed to form a phase conditional population state sequence. Perform population evolutionary relationship rolling fitting and state space local partitioning on the phase condition population state sequence to determine the reference center, adjacency relationship and transition direction of each local partition; Based on the reference center, adjacency relationship and transition direction of each local partition, local Lyapunov evaluation entries are constructed and spliced ​​to form a phase conditional Lyapunov map; Single-step and multi-step model simulations were performed on the phase-conditional Lyapunov map to obtain the population state evolution results, identify the stable core region, critical hinge region, escape boundary and regression path, and form an attraction domain boundary chain validation set; Based on the validation set of the attraction domain boundary chain, the stability level, boundary margin, instability trend and key driving phase of the current population are calculated, and early warning results and control suggestions are generated.

2. The population stability assessment method based on the Lyapunov function as described in claim 1, characterized in that, The steps to obtain the basic dataset of the phase candidate population are as follows: The multi-source records include population size records, viability records, resource availability records, climate fluctuation records, migration activity records, and historical intervention records for the current assessment period and historical assessment periods. The preprocessing includes time unification, missing data completion, anomaly removal, dimension normalization, and confidence labeling. Based on the preprocessed multi-source records, according to the non-autonomous discrete population evolution relationship characterized by population size change, time-varying growth rate and constant survival rate at adjacent time points, candidate parameters of constant survival rate are extracted, and candidate sequences of time-varying growth rate are extracted by combining population size records and candidate parameters of constant survival rate. By combining the candidate sequences of time-varying growth rates, candidate parameters of constant survival rates, and preprocessed multi-source records, resource, perturbation, migration, and reproduction phase identification is determined, and a basic dataset of phase candidate populations is obtained.

3. The population stability assessment method based on the Lyapunov function as described in claim 1 or 2, characterized in that, The steps for forming the phase conditional population state sequence are as follows: The data entries at each time step in the phase candidate population base dataset are written into the corresponding initial population state entries in time index order to obtain the initial population state sequence. Based on the degree of deviation and local fluctuation amplitude of the population size value in the initial population state sequence relative to the unique positive equilibrium point corresponding to the non-autonomous discrete population evolution relationship, the equilibrium deviation normalization value at each time moment is calculated to generate the equilibrium deviation state sequence. The equilibrium deviation normalization value, equilibrium deviation change, time-varying growth rate candidate value, constant survival rate candidate parameter, and various phase identifiers in the equilibrium deviation state sequence are bound to each other and arranged in time index order to obtain the phase conditional population state sequence.

4. The population stability assessment method based on the Lyapunov function as described in claim 3, characterized in that, The steps for determining the reference center, adjacency relationship, and transition direction of each local partition are as follows: Confidence-weighted rolling fitting and aggregate partitioning are performed on the phase-conditional population state sequence to determine the reference center of each local partition; Based on the reference center of each local partition, the state transition between different local partitions is statistically analyzed, and the adjacency relationship and transition direction between each local partition are determined.

5. The population stability assessment method based on the Lyapunov function as described in claim 4, characterized in that, The steps for forming the phase-conditional Lyapunov map are as follows: Around the reference center of each local partition, the deviation trajectory, advancement trajectory and cross-regional transition trajectory of the phase condition population state entries in the corresponding local partition are extracted to obtain the directional interface record; Based on the direction interface records, determine the degree of change of each local partition shrinking towards the reference center, the degree of change of advancement along the jump direction, and the smoothness of the boundary transition, and write the corresponding local Lyapunov evaluation entries to obtain the sequence of local Lyapunov evaluation entries to be spliced. Based on the sequence of local Lyapunov evaluation entries to be spliced, the shared boundary consistency check is performed on adjacent local partitions with adjacency relationship. When the difference between local Lyapunov evaluation entries on both sides of the shared boundary exceeds the preset boundary consistency threshold, the corresponding cross-region transition trajectory, advancement trajectory and boundary correspondence are backtracked and adjusted to obtain a sequence of local Lyapunov evaluation entries with consistent boundaries. The local Lyapunov evaluation entries along the boundary are sequentially spliced ​​together according to the stability of the adjacency relationship and the degree of conflict of the transition direction, forming the phase condition Lyapunov map.

6. The population stability assessment method based on the Lyapunov function as described in claim 5, characterized in that, The steps to obtain the population state evolution results are as follows: The simulation starting point is set around the reference center of each local partition, the location of the shared boundary, and the location of the main transition direction in the phase condition Lyapunov plot; Using the simulation starting point as input, single-step state progression is performed along the phase condition Lyapunov map, and the single-step landing point position, local partition changes, and local Lyapunov evaluation value changes are recorded. Using the single-step landing point location, local partition changes, and local Lyapunov evaluation value changes as inputs, multi-step state advancement is performed, recording the continuous landing point trajectory, cross-region traversal order, and local Lyapunov evaluation value evolution order to form the population state evolution results.

7. The population stability assessment method based on the Lyapunov function as described in claim 1 or 6, characterized in that, The steps for forming the attraction domain boundary chain verification set are as follows: Based on the population state evolution results, regions where continuous landing trajectories remain in the same local partition or stably adjacent local partitions and where the local Lyapunov evaluation value continuously decreases are identified as stable core regions. Based on the population state evolution results outside the stable core region, the region with continuous landing trajectory switching and alternating local Lyapunov evaluation value is identified as the critical articulation region, and the shared boundary of repeatedly crossing out of the current local partition and its stable adjacent local partition is identified as the escape boundary. Based on the population state evolution results corresponding to the critical hinge region and the escape boundary, the cross-regional trajectory for re-entering the stable core region is extracted as the regression path. The stable core region, critical hinge region, escape boundary and regression path are associated and written into the corresponding records to form the attraction domain boundary chain validation set.

8. The population stability assessment method based on the Lyapunov function as described in claim 1 or 7, characterized in that, The steps for generating early warning results and control recommendations are as follows: The stability of the core region, the distribution of the critical articulation region, the escape boundary triggering, the recovery of the regression path, and the frequency of the corresponding phases in the validation set of the attraction domain boundary chain are statistically analyzed, and the stability level, boundary margin, and instability trend of the current population are calculated. The system jointly determines the stability level, boundary margin, and instability trend. It also identifies key driving phases by combining the distribution of critical articulation regions, escape boundary triggering conditions, and the frequency of occurrence of various phases, and generates early warning results and control recommendations.

9. A computer device comprising a memory and a processor, wherein the memory stores a computer program, characterized in that: When the processor executes the computer program, it implements the steps of the population stability assessment method based on the Lyapunov function as described in any one of claims 1 to 8.

10. A computer-readable storage medium having a computer program stored thereon, characterized in that: When the computer program is executed by the processor, it implements the steps of the population stability assessment method based on the Lyapunov function as described in any one of claims 1 to 8.

Images