A method and system for identifying and dynamically evaluating a source of a glacial disaster chain

Abstract

The application discloses a glacier disaster chain source identification and dynamic source evaluation method and system, the method comprises the following steps: preprocessing multi-source basic data, constructing a source state transition matrix to quantify the time sequence evolution law of the source, and calculating a plurality of rule indexes and a rule source potential index SPI by counting a plurality of dynamic indexes of the source. Weak supervision samples are constructed based on the time sequence leakage prevention rule, and a unit high-activity source probability PML is predicted by a tree ensemble model to obtain a high-precision machine learning calibration source potential index ML-SPI. Finally, combined with multi-level quality inspection rules, abnormal identification and potential grading are completed, and various spatial results and evaluation reports are output. The application realizes full-process automatic operation, solves the problems of strong subjectivity, poor time sequence adaptability and low result precision of traditional glacier source evaluation, and can support glacier disaster chain disaster prevention monitoring and risk judgment work.

Description

Technical Field

[0001] This invention belongs to the technical fields of geological disaster risk assessment, remote sensing monitoring of glacier disaster chains, geographic information spatial analysis and data-driven prediction and evaluation, and specifically relates to a method and system for identifying the source of materials and evaluating the dynamic supply of materials in glacier disaster chains. Background Technology

[0002] Within high-altitude canyon glacier development zones, glacial retreat, glacial lake expansion, moraine exposure, rock wall collapse, landslide expansion, gully deposition, and glacial fan growth collectively constitute the material supply system for the glacier disaster chain. Whether the material source enters the gully and forms a continuous supply directly affects the probability of chain disasters such as debris flows, glacial lake outburst floods, river blockages, and landslide dams.

[0003] Existing remote sensing intelligent identification methods mostly revolve around automatic image segmentation or target detection, typically requiring a large number of labeled samples. They often focus on identifying landslides, glacial lakes, or glaciers in single-period images, and use semantic segmentation or target detection models to output land cover boundaries. These methods primarily identify single hazard types, making it difficult to simultaneously represent the multiple source systems within a glacial hazard chain. Furthermore, the model outputs often remain at the land cover distribution layer, lacking standardized databases, confidence markers, and inter-period consistency corrections for subsequent source trend prediction.

[0004] Existing methods for assessing source material in disaster chains are mostly based on single-period remote sensing interpretation results or static susceptibility factors. They emphasize the spatial correlation between current source material distribution and influencing factors such as topography, tectonics, and channels, and use analytic hierarchy process (AHP), entropy weighting, or machine learning methods for evaluation. While these methods can illustrate the correspondence between source material and environmental conditions at a certain period, they struggle to characterize the processes of new additions, declines, persistence, migrations, and type transformations of multiple types of source material across consecutive periods. They also find it difficult to use historical evolution patterns to predict highly active source units in the next stage.

[0005] Existing methods typically treat remote sensing identification, change detection, susceptibility assessment, and results compilation as independent processes. They lack an integrated technical solution that integrates multi-period image or source patch input, source type standardization, inter-period state transition analysis, future source trend prediction, and evaluation results generation. This results in repetitive data processing, untraceable evaluation indicators, and difficulty in automatically calibrating weights based on historical samples and data quality. Consequently, the evaluation results cannot be used in a timely manner to serve disaster prevention and monitoring in high-altitude glacier watersheds.

[0006] Therefore, there is a need for a dynamic source supply evaluation method and system for glacier disaster chains that is compatible with human interpretation results, existing source labeling data, and semantic segmentation model outputs. The method should be able to construct time-series weakly supervised samples by utilizing multi-period source state transition relationships without relying on artificial disaster labels and large-scale relabeling, learn future highly active source supply trends, and combine regular geoscientific indices, machine learning prediction probabilities, quality control, and result verification mechanisms to generate source supply potential classification results, anomaly unit lists, and regional evaluation result reports. Summary of the Invention

[0007] The purpose of this invention is to provide a method and system for identifying and dynamically evaluating the source of glaciers in hazard chains, in order to solve the problems of disconnect between the identification results of multiple types of source and the evaluation of source, the difficulty of single-period static evaluation in depicting the evolution process of multiple sources, the lack of prediction of future high-activity source trends, the difficulty of automatically calibrating evaluation weights, and the lack of quality inspection and reporting output of results.

[0008] To achieve the above objectives, the technical solution of this invention is as follows: A method for identifying the source material of glacier disaster chains and evaluating their dynamic supply includes the following steps: S1. Acquire multi-source data of remote sensing, source patches, and topographic structure of the study area from multiple periods, complete preprocessing, and generate a basic spatial dataset; S2. Read the basic spatial dataset, and complete the unified category mapping and confidence labeling of the three types of source data, which are compatible with manual interpretation, annotation and semantic segmentation, to generate a standardized source patch library; S3. Using the standardized source patch library as the data source, construct the source state transition matrix by spatially superimposing adjacent patches, and output the time-series evolution statistics for dynamic index calculation. S4. Overlay the source state transition matrix with the evaluation unit layer and statistically analyze the dynamic indicators of multiple sources. S5. Using the dynamic indicators and topographic factors from step S4, the four indicators—Dynamic Transformation Intensity Index (DTI), Dynamic Source Accumulation Index (DASCI), Topographic Structure Coupling Index (TTCI), and Channel Connectivity Index (CCI)—are solved by normalization. The weighted sum of these four indicators yields the Regular Source Potential Index (SPI). S6. Extract the source dynamic index, DTI, DASCI, TTCI, CCI and SPI time series quantization data obtained from steps S4 and S5 as model input features, and use the source change results of the next period to generate unlabeled weakly supervised samples, which are then sent to machine learning training. S7. Based on the unannotated weakly supervised sample training tree ensemble model constructed in step S6, the high-activity source probability PML of each evaluation unit is output. PML and SPI are used together for weight fusion calculation. S8. Based on rolling time verification, the internal weights, SPI combined weights, and fusion weights of DTI, DASCI, TTCI, and CCI are automatically calibrated, and the rule-based source potential index SPI and the high-activity source probability PML of each evaluation unit are substituted into the fusion formula to obtain ML-SPI. S9. Based on the ML-SPI, perform multi-level source potential classification, combine the patch, index, and probability data generated in each step of the whole process to perform multi-level quality inspection, and output spatial map, anomaly list and disaster prevention evaluation report.

[0009] Furthermore, step S3 constructs a 9×9 source state transition matrix, with rows and columns corresponding to the source states of two adjacent time periods. State code 0 represents a region without source, and codes 1 to 8 correspond to eight types of standardized glacial disaster chain sources. The value of each element in the matrix is ​​the proportion of the area of ​​the corresponding state transition type to the total area of ​​the evaluation unit.

[0010] Furthermore, the calculation method for the Dynamic Conversion Intensity Index (DTI) is as follows: ; in, The index represents the dynamic conversion intensity of evaluation unit i within time interval t; N represents the normalization function. to These are the weighting coefficients; This indicates the proportion of newly added material sources in evaluation unit i within time interval t; Indicates the category conversion ratio; Indicates the proportion of disaster chain direction transformation; Represents the state transition entropy; Indicates the percentage of continuous area.

[0011] The method for calculating the Dynamic Source Accumulation Index (DASCI) is as follows: ; in: This represents the dynamic material source accumulation index of evaluation unit i within time interval t; Indicates net cumulative density; Indicates the total density of change; Indicates the frequency density of change; Indicates the contribution rate of the dominant material source; to These are the weighting coefficients; The method for calculating the terrain structure coupling index (TTCI) is as follows: ; in: This represents the terrain structure coupling index of evaluation unit i; Indicates the average elevation of the evaluation unit; Indicates the average slope; Indicates the degree of terrain relief; Indicates the distance from the active fault; Indicates peak ground acceleration; This indicates the assignment of slope aspect category; Indicates the lithology category assignment; to These are the weighting coefficients. This represents the negative normalization function; The method for calculating the Channel Connectivity Index (CCI) is as follows: ;

[0012] in, This represents the channel connectivity index of evaluation unit i within time interval t; This represents the arithmetic mean of the effective pixels in the raster grid within the evaluation cell, representing the distance from the nearest channel. This represents the average gradient of the main channel within the evaluation unit; This indicates the cumulative flow at the outlet of the evaluation unit or the maximum cumulative flow in the main channel; This indicates the proportion of the source area intersecting with the channel buffer zone to the total source area of ​​the evaluation unit; to These are the weighting coefficients.

[0013] Furthermore, the calculation method for the rule-based source potential index (SPI) is as follows: ;

[0014] in: This represents the rule-based supply potential index of evaluation unit i within time interval t; This represents the dynamic material source accumulation index of evaluation unit i within time interval t; This represents the terrain structure coupling index of evaluation unit i; This represents the channel connectivity index of evaluation unit i within time interval t; This represents the dynamic conversion intensity index of evaluation unit i within time interval t; to These are the weighting coefficients.

[0015] Furthermore, the calculation method for the high-activity source probability (PML) of each evaluation unit is as follows: Based on the time-series unlabeled weakly supervised samples constructed in step S6, a tree ensemble classification model is selected as the prediction model, and the time-series quantized data of the material source dynamic indicators, terrain environmental factors, DTI, DASCI, TTCI, CCI, and SPI corresponding to each evaluation unit are used as input features; the model divides the training set and validation set according to the chronological order, and obtains the PML according to the probability output mechanism of the selected tree ensemble classification model; among them, random forest outputs the PML based on the positive class voting ratio or the average positive class probability of each decision tree; gradient boosting decision trees, extreme gradient boosting trees, histogram gradient boosting trees, and other gradient boosting models can map the original prediction score to the 0-1 probability interval through the sigmoid function; the model loss function is determined according to the selected model. For probabilistic binary classification models using the binary cross-entropy loss function, the binary cross-entropy loss function can be used for training.

[0016] Furthermore, the ML-SPI calculation method is as follows: ; in, This represents the ML-SPI value of evaluation unit i within time interval t; This represents the probabilistic fusion weights for machine learning, with values ​​ranging from 0 to 1; This represents the rule-based supply potential index of evaluation unit i within time interval t; The evaluation unit i represents the probability of being a highly active source within time interval t; the unified parameter vector used for automatic calibration is... group r includes to Group d includes to , Group includes to Group C includes to Group w includes to , For fusion weights; parameter space Ω constraints r, d, The weights of groups c and w are all not less than 0, and the sum of the weights within each group is 1. The value range is 0 to 1; the optimal parameter is according to Sure, This indicates taking the maximum value.

[0017] A system for identifying and dynamically evaluating the source of glaciers' disaster chains includes a data access and preprocessing module, an output unified spatial dataset, and a downstream module for accessing and standardizing source identification results. The source identification results are accessed and standardized by the module, which reads the spatial dataset to generate a standard source patch library and outputs it to the intertemporal consistency correction and state transition analysis module. The intertemporal consistency correction and state transition analysis module constructs a state transition matrix based on the patch library, and the results are input into the evaluation unit construction and dynamic index statistics module. The evaluation unit construction and dynamic index statistics module calculates a complete set of material source dynamic indicators and transmits them to the rule index calculation module. The rule index calculation module generates the SPI index and simultaneously supplies it to the time-series weakly supervised sample construction module and the weight automatic calibration and index fusion module. The time-series weakly supervised sample construction module generates training samples based on time-series indicators and outputs them to the machine learning supply trend prediction module. The machine learning supply trend prediction module trains the model to output PML, and then passes the PML to the weight automatic calibration and exponential fusion module. The automatic weight calibration and index fusion module receives SPI and PML, calibrates the internal weights, SPI combined weights, and fusion weights of DTI, DASCI, TTCI, and CCI based on rolling time verification, calculates ML-SPI, and pushes them to the quality inspection and review module and the result output module, respectively. The quality inspection and review module is used to read intermediate data from all the aforementioned modules to identify various anomalies; The output module is used to integrate all data such as ML-SPI and quality inspection marks to output various results.

[0018] Furthermore, the data access preprocessing module outputs a unified spatial dataset to the material source standardization module, which then transmits the state matrix, dynamic indicators, and rule indices layer by layer to the machine learning module. The SPI and PML are simultaneously input into the weight fusion calibration module to generate ML-SPI.

[0019] Furthermore, the quality inspection and review module reads source patches, state matrix, SPI, and PML data, identifies missing periods, low-confidence patches, extreme values ​​of indicators, and index grading divergences, classifies them into three levels of quality inspection, and outputs a standardized list of anomalies.

[0020] Furthermore, the output module integrates ML-SPI classification, PML probability, dominant source material, transfer path, and quality inspection markers to generate multiple types of spatial classification maps, a list of high-potential units, and a regional assessment report containing disaster prevention and monitoring recommendations.

[0021] The method of the present invention has the following technical effects: 1. This invention integrates multiple remote sensing images, manual interpretation results, existing source labeling data, and semantic segmentation model outputs into a unified database, and converts them into a standardized multi-type source patch library. This solves the problem of difficulty in uniformly managing source identification results from different sources and participating in subsequent source evaluation.

[0022] 2. This invention, through intertemporal consistency correction and a 9×9 source state transition matrix, can characterize the processes of new addition, decline, persistence, migration, and category transformation of sources in consecutive periods of the glacier disaster chain. Compared with single-period static source evaluation, it can better reflect the dynamic evolution characteristics of sources.

[0023] 3. This invention constructs a dynamic conversion intensity index, a dynamic source accumulation index, a topographic structure coupling index, a channel connectivity index, and a regular source potential index, enabling source evolution, topographic structure, and channel connectivity conditions to participate in source potential assessment in a calculable and traceable manner.

[0024] 4. This invention uses the actual changes in material sources in the next time interval to generate time-series weakly supervised labels, without relying on artificial disaster point labels, and can make full use of multi-period material source evolution data to learn the spatial distribution patterns of future highly active material source units.

[0025] 5. This invention automatically calibrates the internal weights of the rule index and the weights of the machine learning probability fusion based on the rolling time verification results, which solves the problem that the weights in the traditional source potential evaluation rely on human experience and are difficult to adjust with historical samples and model performance.

[0026] 6. This invention integrates the regular source potential index with the machine learning prediction probability to obtain the machine learning-calibrated source potential index, so that the evaluation results simultaneously possess geoscientific interpretability, historical evolution learning ability, and future source trend prediction ability.

[0027] 7. This invention, through quality inspection and verification mechanisms such as missing periods, low-confidence patches, small-area anomalies, extreme change density, state transition anomalies, and model rule divergences, can automatically identify evaluation units that need to be verified, thereby improving the reliability of evaluation results and the safety of engineering applications.

[0028] 8. This invention can output source potential classification maps, machine learning probability maps, dominant source type maps, dominant transfer path maps, dynamic conversion intensity index spatial distribution maps, model interpretation maps, anomaly unit quality inspection maps, and regional evaluation results reports, which facilitates disaster prevention monitoring, engineering geological safety evaluation, and results archiving in high-altitude glacier watersheds. Attached Figure Description

[0029] Figure 1 This is a schematic diagram of the overall process of the method of the present invention.

[0030] Figure 2This is a structural diagram of the system modules of the present invention.

[0031] Figure 3 This is a schematic diagram illustrating the access and eight-category standardization of multi-source material identification results in this invention.

[0032] Figure 4 This is a schematic diagram illustrating the intertemporal source state transition and dynamic index calculation of the present invention.

[0033] Figure 5 This is a schematic diagram illustrating the construction of time-series weakly supervised samples and leakage prevention training in this invention.

[0034] Figure 6 This is a flowchart of the SPI, PML and ML-SPI fusion evaluation process of the present invention.

[0035] Figure 7 This is a schematic diagram illustrating the output and quality inspection of the results of this invention.

[0036] Figure 8 This is a schematic diagram illustrating the evaluation results and model verification of an embodiment of the present invention.

[0037] Figure 8 (a) in the figure is a schematic diagram of the number of ML-SPI power supply potential levels in an embodiment of the present invention.

[0038] Figure 8 (b) in the figure is a schematic diagram of the PML probability distribution in an embodiment of the present invention.

[0039] Figure 8 (c) in the figure is a schematic diagram of the top ten features of the embodiments of the present invention.

[0040] Figure 8 (d) in the figure is a schematic diagram of the rolling verification performance of the model in the embodiment of the present invention.

[0041] Figure 8 (e) in the figure is a hierarchical diagram of the ML-SPI power supply potential space in an embodiment of the present invention. Detailed Implementation

[0042] To make the objectives, technical solutions, and advantages of this invention clearer, the invention will be further described in detail below with reference to the accompanying drawings and embodiments. It should be understood that the specific embodiments described herein are only some, not all, of the embodiments of this invention, and are not intended to limit the invention. All other embodiments obtained by those skilled in the art based on the embodiments of this invention without inventive effort are within the scope of protection of this invention.

[0043] This invention discloses a method for identifying and evaluating the dynamic supply potential of glacial disaster chains. It is used to integrate multi-type glacial source identification results in high-altitude glacial basins, perform intertemporal evolution analysis, predict dynamic supply potential, and generate evaluation results. The dynamic supply potential of multi-type glacial sources in glacial disaster chains refers to the dynamic capacity and future potential activity of various loose materials (sources) within high-altitude glacial basins as they continuously evolve over time and replenish channels. This method assesses whether, with what probability and scale, glacial materials such as glaciers, moraines, landslides, and rock debris will continuously replenish channels in a region, potentially triggering chain disasters such as glacial lake outburst floods, debris flows, and river blockages. It is a comprehensive quantitative evaluation of the dynamic supply capacity and future disaster-causing potential of glacial sources. Figure 1 This is a schematic diagram of the overall process of the method of the present invention, showing the complete process of multi-source data acquisition, standardization of source identification results, inter-period consistency correction, state transition analysis, dynamic index calculation, SPI rule index calculation, time-series weakly supervised learning, PML prediction, ML-SPI fusion evaluation, quality inspection and verification, and output of results. The method of the present invention includes the following steps: S1: Acquisition of multi-source data and basic spatial preprocessing of the study area.

[0044] The study area boundary of the target alpine glacier basin is acquired, along with multiple periods of remote sensing imagery, and at least one source identification result from manually interpreted source patches, existing source-labeled data, or semantic segmentation model outputs. Auxiliary geospatial data is also acquired. This auxiliary geospatial data includes digital elevation models, channel networks, basin boundaries, active faults, lithology, peak ground acceleration, and other environmental data related to the source supply of the glacier hazard chain.

[0045] A unified target coordinate system, spatial extent, raster resolution, and time period index were determined based on the study area boundaries. Coordinate transformation, spatial cropping, and image registration were performed on multi-period remote sensing images. Coordinate transformation, spatial cropping, resolution matching, and integrity checks were performed on auxiliary geospatial data. Environmental factors such as slope, aspect, topographic relief, and gully distance were calculated based on digital elevation models, gully networks, and other basic data. The data source, file format, spatial coverage, and time period of the source identification results were recorded to form a basic spatial dataset for the study area.

[0046] This invention takes a high-mountain canyon glacier disaster chain development area that has completed multiple phases of remote sensing interpretation and watershed delineation as the object. It uses five phases of source patch data from 2016, 2018, 2020, 2022 and 2024, evaluation watershed boundary data and auxiliary geospatial data to implement the glacier disaster chain source identification and dynamic source evaluation method described in this invention.

[0047] The five-phase source patch data were compiled from multiple phases of remote sensing interpretation results and existing source boundary data, and are areal vector data. Each patch includes at least the patch boundary, year, source type, and area attribute. The evaluation watershed boundary data includes 1192 evaluation watersheds, each with a unique number. The auxiliary geospatial data includes environmental factors such as digital elevation model, slope, aspect, topographic relief, gully network, active faults, lithology, peak ground acceleration, and gully distance. Before execution, the above data undergoes unified coordinate projection, spatial clipping, geometric restoration, field standardization, category coding, and area recalculation to ensure that source patches, evaluation watersheds, and environmental factors are under the same spatial reference and field caliber. In this embodiment of the invention, the number of source patches, the number of evaluation watersheds, the state transition area, the model validation indicators, and the quality inspection records all originate from the same study area data processing process. Spatial overlay, area statistics, sample construction, model validation, and quality inspection review are all completed using a unified data caliber.

[0048] S2: Access and standardization of source identification results for multi-type glacier disaster chains.

[0049] It receives one or more source identification results from manual interpretation results, existing source-labeled data, and semantic segmentation model outputs. For manual interpretation results and existing source-labeled data, it reads the vector boundaries of source patches and identifies the original category field, time field, data source field, and confidence field. For semantic segmentation model outputs, it receives category probability maps, raster classification maps, or vector candidate patches, converts the category probability maps or raster classification maps into candidate source patches, and records their category, average identification confidence, boundary confidence, and data source. Figure 3 This diagram illustrates the process of integrating multi-source provenance identification results with eight standardized categories. It shows how, after field cleaning, geometric repair, category mapping, and confidence recording, the results of manual interpretation, existing provenance labeling data, and semantic segmentation model output are uniformly formed into a standardized provenance patch library of eight categories: glaciers, glacial lakes, glacial till, rock wall deposits, glacial hazardous rock masses, landslides, gully deposits, and depositional fans.

[0050] The original material sources are uniformly mapped to eight categories of glacial hazard chain material sources, namely glaciers, glacial lakes, glacial moraines, rock wall deposits, glacial hazardous rock masses, landslides, channel deposits, and depositional fans. Category codes 1 to 8 correspond to the above eight material sources, and status code 0 indicates an area without material sources. The time intervals between adjacent original material sources include 2016 to 2018, 2018 to 2020, 2020 to 2022, and 2022 to 2024. The source identification results from various sources are converted to a unified target coordinate system and study area spatial range. Invalid geometric repair, topological checking, overlapping conflict checks between different categories within the same year, time field normalization, attribute field normalization, area recalculation, and small-area patch verification marking are performed on each source patch. A unique number is assigned to each patch, recording the standard category name, category code, year, area, data source, identification confidence level, boundary confidence level, and quality control mark. For patches with missing years, unmappable categories, empty geometry, zero or negative areas, or overlapping conflicts between different categories within the same year, the original records are retained and corresponding quality control marks are set. After verification, it is determined whether to proceed to subsequent analysis. Patches with areas smaller than a preset threshold are marked as small-area verification patches. The area threshold is determined based on the image spatial resolution and study area scale; the preferred small-area patch verification threshold is 0.001 square kilometers. Patches below this threshold are marked as small-area verification patches as needed for research. The small-area patch threshold is used for verification marking rather than batch deletion, facilitating subsequent quality control and manual verification. After completing the above processing, a standardized source patch library with unified fields, clear categories, and quality tags is formed.

[0051] The standardized source patch library includes at least the following fields: unique patch number, year, standard category name, category code, area in square meters, and area in square kilometers. In a specific embodiment of the present invention, a total of 45,020 standard source patches were obtained after standardization of the five phases of source patches. The field integrity check results show that the number of patches with unidentified categories, abnormal years, empty geometry, and areas with zero or negative values ​​is zero; overlapping patches of different categories in the same year, spatial overlay tail differences, and local area closure errors are not directly deleted in this step, but are entered into the subsequent state transition quality inspection and comprehensive quality inspection checklist.

[0052] The weak supervision label threshold is determined by the 80th percentile of four indicators: total change density, change frequency density, dynamic conversion intensity index, and net cumulative density within the next adjacent time interval. High activity is defined when any one of these conditions is met, with the net cumulative density condition requiring a value greater than 0. This embodiment uses the optimal fusion weight determined through rolling time validation. =0.6; the potential for resource generation is divided into four levels: low, medium, relatively high, and high; the top 20% of high-potential resource generation units are identified according to the machine learning-calibrated resource generation potential index (ML-SPI). All the above parameters are written into the processing record and saved along with the results table. For different study areas, the above thresholds can be adjusted within the limits specified in the instructions, but a uniform parameter caliber is used in the same evaluation task.

[0053] The integration of multi-source prototyping results and the standardization of eight categories specifically include five steps: integration of manual interpretation results, integration of existing prototyping labeled data, integration of semantic segmentation output results, multi-source conflict checking and standard category mapping, and generation of a standard prototyping patch library. The multi-source prototyping results include three categories: manual interpretation results, existing prototyping labeled data, and semantic segmentation model output results.

[0054] The manually interpreted results are source patch vector data obtained by manually delineating multiple periods of remote sensing imagery. Each patch includes at least the patch boundary, interpretation year, original land cover name, and area attribute. The existing source labeling data are source boundary data formed from historical surveys, thematic layers, or previous research results. Each patch includes at least the labeling category, spatial boundary, data source, and survey time. The semantic segmentation model output is a category probability map, raster classification map, or candidate source patches generated by the remote sensing image segmentation model. The semantic segmentation model is a pixel-level segmentation network for remote sensing images, including but not limited to U-Net, DeepLabv3+, and Mask2Former.

[0055] In this embodiment of the invention, the system uniformly converts three types of inputs—manual interpretation results, existing source labeling data, and semantic segmentation model output results—into area source patch data and establishes a unified field structure. The unified fields include at least the patch's unique ID, year, original category name, standard category name, category code, area, data source, identification confidence level, boundary confidence level, and quality inspection mark.

[0056] The data source field is set to three categories: manual interpretation, existing annotation, and semantic segmentation, used to record the source type of each source patch. The identification confidence and boundary confidence values ​​both range from 0 to 1. The low confidence threshold can be set according to image quality and model performance, preferably between 0.6 and 0.8. Patches below the threshold are not directly deleted, but are marked as low-confidence patches or patches requiring boundary verification, and added to the quality inspection checklist.

[0057] The integration of manually interpreted results involves the system reading multi-period provenance patch vector files generated through a multi-source identification result integration module. The system first checks the validity of the patch geometry; if self-intersections, empty geometry, or topological anomalies exist, geometric repair is performed. Then, the original category, year, and area fields are identified and mapped to unified fields. For manually interpreted results, the system matches the original category names with eight standard provenance names. For example, "glacier," "glacier," or synonyms in the original fields are mapped to "glacier"; "glacial lake," "glacial lake," or synonyms are mapped to "glacial lake"; the remaining categories are mapped sequentially to glacial till, rock wall deposits, glacial hazardous rock masses, landslides, channel deposits, and depositional fans. After mapping, the system generates a unique patch number for each patch and records its data source as manually interpreted.

[0058] Accessing existing source-labeled data involves the system reading source patch data from historical source boundaries, thematic survey layers, or existing research findings through the existing labeled data access module. Because existing labeled data may have inconsistencies in field names, category systems, or incomplete year information, the system first identifies and standardizes the fields. When existing labeled data contains a specific year, the system directly retains that year information; when existing labeled data only contains the survey time or image acquisition time, the system converts it to the corresponding year field; when existing labeled data lacks year information, the system marks the patch as a patch requiring subsequent review. Before the year information is completed, the patch can be retained as auxiliary source identification evidence but will not participate in the calculation of inter-period state transition matrices, time-series weakly supervised samples, or dynamic indicators. For existing labeled data with inconsistent categories, the system uniformly merges them into eight standard source categories according to the category mapping table, and records the original category name, standard category name, and mapping result in the quality inspection field.

[0059] The semantic segmentation output result access module allows the system to read the category probability map, raster classification map, or candidate source patches output by the semantic segmentation model. When the input is a raster classification map, the system merges adjacent pixels of the same category into candidate patches and performs vectorization, hole filling, fragmentation merging, and small area labeling. When the input is a category probability map, the system determines the category of the candidate patch according to the highest category probability and calculates the average category probability within the patch area as the recognition confidence. The system further calculates the boundary confidence based on the category probability gradient near the patch boundary, the degree of boundary fragmentation, the proportion of holes, or the results of manual review. When the recognition confidence or boundary confidence is lower than a preset threshold, the system does not directly remove the patch but marks it as a low-confidence patch or a patch requiring boundary review in the quality inspection field.

[0060] Multi-source conflict checking and standard category mapping involve the system not directly deleting any source result when duplicate patches or conflicting patches from different data sources exist within the same spatial range. Instead, conflict marking is performed based on data source priority, identification confidence, boundary confidence, and time period consistency. In a preferred embodiment, manually reviewed interpretation results have higher priority than unreviewed semantic segmentation results; high-confidence semantic segmentation patches have higher priority than low-confidence already labeled patches. For conflict patches whose priority cannot be automatically determined, the system retains the original source information and marks them as multi-source conflict patches requiring review in the quality inspection field. After completing the conflict check, the system uniformly maps all patches to eight standard sources, generating a unified category code. For categories that cannot be automatically mapped, the system marks them as unrecognized categories and outputs them to the quality inspection checklist.

[0061] The generation of the standard provenance patch library involves the system uniformly writing data from three sources into the library through a provenance standardization module. Each patch in the library has a unified year field, standard category name, category code, area field, data source field, identification confidence field, boundary confidence field, and quality inspection field. During the standardization process, the system recalculates the patch area, unifying the area units to square meters and square kilometers; unifies the coordinate projection to ensure all patches are in the same spatial reference; unifies the time period to allow provenance patches within the same year or adjacent time intervals to participate in subsequent state transition analysis; and unifies the category code so that all patches are mapped to standard provenance categories 1 to 8.

[0062] Through the five steps described above, the results of manual interpretation, existing source labeling data, and semantic segmentation output are all converted into a standardized source patch library in a unified format. Data from different sources all have unified year, category, area, source, and confidence fields before entering the subsequent dynamic source evaluation process. The standardized source patch library can be directly used for evaluation unit intersection, 9×9 state transition matrix construction, dynamic index statistics, time-series weakly supervised sample construction, and source potential evaluation. For patches with different sources, confidence levels, or boundary reliability, the system can distinguish them through data source, identification confidence, and boundary confidence fields, thereby marking and verifying low-quality patches in the subsequent evaluation process. The integration of multi-source source identification results demonstrates that this invention does not rely on a single source identification method. Regardless of whether the source identification results come from manual interpretation, existing labeled data, or semantic segmentation model output, they can all be standardized into eight categories of source patch libraries through unified field cleaning, category mapping, confidence recording, conflict marking, and quality checks, thus ensuring consistent data entry for the subsequent dynamic source evaluation process.

[0063] S3: Intertemporal Consistency Correction and Source State Transition Analysis.

[0064] Spatial overlay of source patches from two adjacent periods is used to identify newly added patches, fading patches, persistent patches, migrating patches, and category-transformation patches. This spatial overlay is a GIS vector-polygon feature overlay analysis. The layers of source patches from two adjacent periods, after unifying coordinates and restoring geometry, are sequentially subjected to topological operations such as intersection, erasure, and spatial distance calculation. Intersection yields the overlapping area between the two periods; erasure extracts the unique areas of each period; and the Euclidean distance between the patch centroids is calculated to determine migrating patches, thus identifying newly added, fading, persistent, migrating, and category-transformation patches. Newly added patches are source patches that did not exist in the previous period but appear in the current period; fading patches are source patches that existed in the previous period but have disappeared in the current period; persistent patches are source patches that exist in both adjacent periods and overlap spatially; migrating patches are source patches that exist in both adjacent periods but whose centroid or boundary displacement exceeds a preset distance threshold; and category-transformation patches are patches that spatially overlap in both adjacent periods but whose source category changes. Figure 4 It is a schematic diagram of inter-period source state transition and dynamic index calculation. It shows the process of identifying newly added, regressing, persistent, migrating and category-transformed patches after the superposition of source patches in adjacent periods, and constructing a 9×9 source state transition matrix to further calculate dynamic indicators such as net cumulative area, total changed area, change frequency density, persistent area ratio, category transformation ratio, disaster chain direction transformation ratio and state transition entropy.

[0065] The semantic segmentation model output, including source patches (raster, vector patches, category identification results, patch boundaries, etc.), is further consistent with the logic of topography, gullies, faults, slope, aspect, and intertemporal changes. If the same area exhibits frequent category jumps that do not conform to the evolution of the glacial hazard chain in adjacent periods, or if there are cases with low identification confidence, fragmented patch boundaries, or abnormal category transformations, the system marks these patches as requiring review and records the anomaly type.

[0066] Furthermore, the source state within the evaluation unit is defined as nine states: 0 represents no source, 1 represents glacier, 2 represents glacial lake, 3 represents glacial moraine, 4 represents rock wall deposits, 5 represents glacial dangerous rock mass, 6 represents landslide, 7 represents gully deposits, and 8 represents depositional fan.

[0067] For each evaluation unit and each adjacent time interval, a 9×9 source state transition matrix is ​​constructed. The matrix elements represent the proportion of the transition area from state a in the previous period to state b in the next period within the evaluation unit, calculated using the following formula: ;

[0068] in, This represents the area ratio of evaluation unit i that transitions from state a to state b within time interval t. Each group (a, b) corresponds to a transition path, and each transition path corresponds to a transition area or a transition area ratio. Within the same evaluation unit and the same time interval, the a→b state transition path with the largest transition area or transition area ratio is the dominant transition path. This represents the spatial range of evaluation unit i in the previous period when its state was a; This represents the spatial range in which evaluation unit i is in state b in the next period; This represents the area of ​​evaluation unit i; the values ​​of a and b both range from 0 to 8. This indicates finding the intersection.

[0069] When a=0 and b≠0, it indicates that the region without material source is transformed into the region with material source, corresponding to the addition of a material source; when a≠0 and b=0, it indicates that the region with material source is transformed into the region without material source, corresponding to the disappearance of material source; when a=b and a≠0, it indicates that the same type of material source continues to exist; when a≠b and a≠0 and b≠0, it indicates that a category transformation occurs between different types of material sources.

[0070] To ensure the state transition matrix can be calculated, the system retains four types of records during spatial overlay: overlapping records where there are sources in both the previous and next periods, newly added records where there were no sources in the previous period, fading records where there were sources in the previous period but no sources in the next period, and records where no sources are retained within the evaluation unit. This process also addresses topological conflicts between overlapping patches of different categories in the same year, small negative areas caused by floating-point precision in spatial overlay operations, and discrepancies between the sum of the watershed source transfer areas and the total watershed area (less than 1 × 10⁻⁶). -8 The system judges the closure tail error as a calculation accuracy defect and does not use it as a valid basis for the actual source transfer. It is uniformly written into the quality inspection list and marked in the source potential classification map and regional results report to prompt manual review of the patch boundary and area calculation results.

[0071] The inter-period consistency correction and source state transition analysis module is invoked to spatially overlay the standard source patches of two adjacent periods with the boundaries of 1192 evaluation watersheds, and construct a 9×9 source state transition matrix for each evaluation watershed and each adjacent time interval.

[0072] The matrix rows represent the state of the previous period, and the matrix columns represent the state of the next period. For any evaluated watershed, if the previous period's value is 0 and the next period's value is non-zero, it is determined to be a newly added sediment source; if the previous period's value is non-zero and the next period's value is 0, it is determined to be a sediment source decline; if the previous period and the next period have the same non-zero state, it is determined to be a continuation of the same category; if the previous period and the next period have different non-zero states, it is determined to be a category transformation; if both the previous period and the next period's value are 0, it is determined to be no sediment source retention. The no sediment source retention term is retained in the matrix and used for watershed area closure checks.

[0073] Through intertemporal consistency correction and source state transition analysis, the system outputs a state transition detail table, a state transition summary table, and a state transition quality inspection table.

[0074] S4: Evaluation unit construction and watershed-scale sediment source dynamic index calculation.

[0075] Evaluation units are selected from watersheds, sub-watersheds, gully units, or regular grids. Preferably, flow direction, runoff accumulation, and gully network are extracted based on a digital elevation model, and evaluation units are generated by combining them with existing watershed boundaries.

[0076] Each evaluation unit is assigned a unique number, and its area, average elevation, average slope, topographic relief, distance to active fault, distance to gully, peak ground acceleration, dominant lithology, dominant aspect, and other environmental factors are calculated. Specifically: the evaluation unit area is recalculated using a unified projected coordinate system; the average elevation and average slope are the arithmetic mean of all raster pixels within the unit; the topographic relief is calculated using a fixed analysis window to determine the difference between the maximum and minimum elevations within the unit; the distance to active fault and peak ground acceleration are the weighted average of the source patch areas within the unit; the distance to the gully is the arithmetic mean of the effective pixels in the raster at the nearest gully within the evaluation unit; the dominant lithology, dominant aspect, and other environmental factors are calculated. The lithology and slope aspect categories that account for more than 50% of the area within a unit are prioritized for the guiding slope direction. If no category accounts for more than 50% of the area, the category with the largest area proportion is selected as the dominant category. If multiple categories have the same largest area proportion, the category is determined according to the preset category priority, and the corresponding area proportion is recorded. After determining the dominant category, it is converted to a score in the 0-1 range according to the preset level. For units without active fault data, the distance from the fault is assigned the maximum value of the entire region. For units with no valid pixels in the raster at the distance from the channel, the distance from the channel is assigned the maximum value of the entire region. For evaluation units with too small an area, too many missing periods, or that may lead to abnormal density indicators, the system marks them as units that need to be reviewed in the quality inspection table.

[0077] The multi-period source patches, changing patches, and state transition matrices are spatially intersected with the evaluation units, and the source dynamic indicators for each evaluation unit are statistically analyzed within each time interval. These source dynamic indicators include the area of ​​various sources, the newly added area of ​​various sources, the area of ​​various sources regressing, and the net cumulative area. The metrics include total change area, total change density, change frequency density, persistent area ratio, category transformation ratio, disaster chain direction transformation ratio, state transition entropy, dominant source type, and dominant transfer path. The area of ​​each source type is the sum of the areas grouped by source type after the intersection of a single-period source patch and the evaluation unit. The newly added area of ​​each source type is obtained by superimposing patches from two adjacent periods, extracting areas where there was no source in the previous period but the corresponding type of source appeared in the next period. The area of ​​each source type that has receded is obtained by superimposing patches from two adjacent periods, extracting areas where the corresponding source existed in the previous period but disappeared in the next period. The above newly added and receding areas are determined and statistically analyzed based on the 0 or non-0 source transformation rules of the 9×9 source state transition matrix.

[0078] Net cumulative area Calculation: ; Total area of ​​change Calculation: ;

[0079] Total change density Calculation: ;

[0080] Frequency density of change Calculation: ;

[0081] in, This represents the newly added area of ​​material sources in evaluation unit i within time interval t, and the newly added area of ​​various material sources. By spatially superimposing source patches from two adjacent periods, regions that were in a source-free state (coded as 0) in the previous period and were coded as source types 1-8 in the next period are extracted, and the areas are summed up according to source type. The total newly added source area of ​​the evaluation unit is the sum of the newly added areas of all categories within the unit. The area representing the receding source is obtained by spatially superimposing source patches from two adjacent periods. The area with source status codes 1-8 in the previous period and no source status code 0 in the next period is extracted. The total receding area of ​​the unit is obtained by summing the areas according to the source category. The receding area of ​​each type of source is the summed area of ​​the receding area of ​​the corresponding type. This represents the number of patches that have changed. Changed patches include four categories: newly added, regressed, migrated, and category-transformed patches, excluding patches of the same type that persist across different periods. All changed patches are extracted by overlaying source patches from two periods, then intersecting with the evaluation unit space to count the total number of changed patches within each unit. . This represents the length of adjacent time intervals, preferably in years; when all adjacent time intervals are the same length... Take a uniform time length, and when adjacent time intervals have different lengths, assign values ​​according to the corresponding time lengths.

[0082] Persistent area ratio The calculation formula is: ;

[0083] in, This represents the area of ​​the same type of material source that continuously exists in evaluation unit i within time interval t. This represents the area of ​​the combined source of the previous period and the current period within the time interval for evaluation unit i; in other embodiments, The above calculation formula can also use the total area of ​​the previous period's source as the denominator, but a unified approach should be used in the same evaluation task; when When it is 0, The value is 0, and the area without material source is recorded in the quality inspection form.

[0084] Category conversion rate The calculation formula is: ; In the formula, A set of combinations of category transformation states; and Representing the evaluation units respectively In the Period and the The state categories of the period are consistent with the meaning of the state symbols in step S3 above, but it should be noted that state 0 does not participate in the calculation of the category transformation ratio, so neither of them can take the value of 0 here; the summation symbol represents the accumulation of the overlapping surfaces corresponding to all ordered state combinations in the set.

[0085] Disaster chain direction conversion ratio The calculation formula is: ;

[0086] In the formula, the double summation represents traversing the set. All source categories in the internal source region With sets All types of sediment in the internal channels For all pairwise combinations, sum up the corresponding spatial overlap areas. Among them, Representing the Source region source state category during the period Representing the The category of channel transport or accumulation state during the period; preferably, These correspond to six types of hazard sources in order: glaciers, glacial lakes, glacial till, rock wall deposits, dangerous glacial rock masses, and landslides. These correspond to two types of disaster depositional landforms: channel deposits and depositional fans. For the first Within the evaluation unit, the first Period Category The corresponding spatial distribution area of ​​the material source, For the first Within the evaluation unit, the first Period Category The corresponding spatial distribution area of ​​the accumulation volume, This is a spatial area calculation operator. In this embodiment of the invention, the study area is divided into closed mountain watersheds according to topographic watersheds. The cross-temporal spatial overlap area is used to characterize the spatial correlation degree of the transformation from the source state of the source area to the state of channel deposits or accumulation fans. It is a proxy indicator for the transformation of the disaster chain direction and is not equivalent to the determination of the material source or the actual transport trajectory alone. When interpreting the results, a comprehensive judgment is made in combination with the topographic confluence path, channel connectivity conditions and manual verification results.

[0087] For the first The total area of ​​each evaluation unit is used as the indicator normalization benchmark to characterize the proportion of the disaster chain direction transformation area within the evaluation unit to the total area of ​​the unit. Using evaluation units with uniform scale and homogeneous geomorphic attributes, the evaluation benchmark for all units can be unified by normalizing the total area of ​​the entire region. This avoids the problems of inconsistent unit benchmarks and ineffective regional horizontal comparisons caused by the normalization of effective source area, ensuring the comparability of disaster chain transformation activity in different spatial units. In this embodiment, the time-series monitoring data uses four adjacent time intervals: 2016-2018, 2018-2020, 2020-2022, and 2022-2024, each interval being 2 years. Therefore, this embodiment uniformly uses a time interval length of 2 years. For other study areas or evaluation tasks, adjacent time intervals can be set to equal or unequal lengths depending on data availability. When using unequal time intervals, standardization is performed according to the actual time length of each interval, and a unified time unit and calculation caliber are used within the same evaluation task to ensure indicator comparability. The spatial development scale and activity evaluation index of disaster chain direction transformation is used to characterize the spatial distribution characteristics of chain transformation of mountain disaster source-deposit body. It does not involve disaster risk intensity classification. Therefore, it adopts equal weight calculation for various transformation combinations, which can uniformly characterize the spatial transformation and development characteristics of disaster chains of different types of sources.

[0088] State transition entropy The calculation formula is: =- ;

[0089] in, This represents the area ratio of the k-th effective state transition path of evaluation unit i within time interval t; K is the total number of all effective state transition paths of evaluation unit i within time interval t; i and t have the same meaning as above, representing the evaluation unit number and the adjacent time interval from period t to period t+1, respectively. Step S3 constructs a 9×9 source state transition matrix, with a total of 81 potential transition paths for all (a,b) combinations. This embodiment of the invention defines 9 source states: a, b {0,1,2,…,8}, where a represents the source state in period t, and b represents the source state in period t+1 (the next period), and the meanings of a and b here are consistent with those in step S3; each combination of a transforming into b constitutes a source state transformation path. An effective path is equal to the actual overlapping transformation area of ​​a transforming into b within this evaluation unit and time period. The determination rule is: >0; If the area of ​​intersection of the two states is greater than 0, then this (a→b) is a valid path and participates in the entropy calculation; if the area of ​​intersection is equal to 0, then this transformation has never occurred and is an invalid path, which is excluded from the summation. For example, if a watershed only has two transformations, namely, the glacier remains unchanged and the glacier transforms into glacial moraine, and the area of ​​the remaining 79 paths is all 0, then there are only 2 valid paths.

[0090] All the above dynamic indicators are calculated from the standard source patch library, evaluation unit boundaries, and 9×9 state transition matrix. , , , All variables are based on the area recalculated in the projected coordinate system. When the area of ​​the evaluation unit is too small, the source of a certain period is missing, or the density index is abnormal due to the overlap of patches in the same year, the system retains the original record and marks it as needing to be reviewed, and does not use it to directly modify the evaluation results.

[0091] The evaluation unit construction and dynamic index statistics module is invoked, using the evaluation watershed as the basic evaluation unit, to statistically analyze the dynamic indicators of material sources for each evaluation watershed in each adjacent time interval. These dynamic indicators include net cumulative area, total changed area, total changed density, changed frequency density, persistent area ratio, category transformation ratio, disaster chain direction transformation ratio, state transition entropy, dominant material source type, and dominant transfer path.

[0092] S5: Calculation of dynamic conversion intensity index and regular supply potential index.

[0093] Building upon step S4, the regular index calculation module is invoked. Following the normalization function, negative normalization function, and weight constraints in the regular index calculation steps of this invention, the Dynamic Transformation Intensity Index (DTI), Dynamic Source Accumulation Index (DASCI), Topographic Coupling Index (TTCI), Channel Connectivity Index (CCI), and Regular Source Potential Index (SPI) are calculated respectively. These indices correspond to four control links in the glacial hazard chain source supply process: the DASCI index corresponds to dynamic source accumulation, the TTCI index corresponds to topographic tectonic background, the CCI index corresponds to channel connectivity conditions, and the DTI index corresponds to inter-period transformation intensity. Specifically, the Dynamic Source Accumulation Index characterizes the scale and sustainability of source replenishment, the Topographic Coupling Index characterizes the instability background, the Channel Connectivity Index characterizes the probability of source entering channels, and the Dynamic Transformation Intensity Index characterizes the activity level of multi-period source state transitions. Figure 6 This is a flowchart of the SPI, PML, and ML-SPI fusion evaluation process, showing the fusion of rule-based indicators such as DASCI, TTCI, CCI, and DTI to form SPI, the output of the time-series weakly supervised machine learning model as the PML probability, and the fusion weights. By fusing SPI and PML, ML-SPI is obtained, which then completes the process of source potential ranking, classification, and high potential cell identification.

[0094] The proportion of new material sources and the proportion of category conversion obtained based on step S4. , proportion of disaster chain direction conversion State transition entropy and continuous area ratio The Dynamic Conversion Intensity Index (DTI) is calculated using the following formula: ; Where N represents the normalization function, to These are the weighting coefficients. The proportion of newly added material sources in evaluation unit i within time interval t is represented by the area of ​​newly added material sources. Divide by the area of ​​the evaluation unit to obtain. And: ; In a preferred embodiment, , , , , The initial values ​​are 0.20, 0.25, 0.30, 0.15 and 0.10 respectively, and can be calibrated by rolling time verification. This set of weights is used as the r group in the unified parameter vector Θ and is optimized according to the parameter space constraint and group iterative calibration method described in S8.

[0095] The normalization function N(x) preferably uses quantile-truncation normalization, as shown in the formula: = ( ); in, and They represent the indicators to be normalized in the training samples, respectively. The 5th percentile and the 95th percentile; This means truncating the calculation result to the range of 0 to 1. When 0 is taken, Take 1 in the first case, and take 1 in the other cases. .when If the indicator is missing entirely, the normalized value of the indicator is set to 0, and the corresponding field is recorded as abnormal in the quality control table.

[0096] For negative factors such as distance from the channel and distance from the active fault, where smaller values ​​indicate higher source potential, a negative normalization function is used: ; Based on the dynamic indicators of the source material, topographic structural factors, and channel connectivity factors of the evaluation unit, the dynamic source material accumulation index (DASCI), topographic structural coupling index (TTCI), channel connectivity index (CCI), and dynamic conversion intensity index (DTI) are calculated respectively, and then fused to obtain the regular source material potential index (SPI).

[0097] The formula for calculating the Dynamic Source Accumulation Index (DASCI) is as follows: ; Represents net cumulative density, derived from net cumulative area. Divide by the unit area to obtain; Indicates the total change density. Indicates the frequency density of change. Indicates the proportion of continuous area. This indicates the contribution rate of the dominant sediment source. Within the same evaluation unit and the same time interval, the sediment source with the largest area among the eight types is the dominant sediment source type for that unit in that period. It equals the area of ​​the largest material source within the current period in the evaluation unit, divided by the total area of ​​all material sources in the unit during the current period. to These are the weighting coefficients, and ; , , , , The initial values ​​are 0.25, 0.25, 0.20, 0.15, and 0.15 respectively, and can be automatically calibrated and optimized through rolling time verification. This set of weights is used as group d in the unified parameter vector Θ and optimized according to the parameter space constraints and group iterative calibration method described in S8.

[0098] The formula for calculating the Topographic Coupling Index (TTCI) is as follows: ; in, The average elevation of the evaluation unit is the arithmetic mean of all DEM grids within the unit. This represents the average slope, calculated as the arithmetic mean of the slope grid within a cell. To represent the terrain relief, the difference between the maximum and minimum elevations within a 5×5 grid window is used for calculation. It represents the distance to the active fault, preferably the area-weighted average of the distance from the source patch to the nearest active fault within the evaluation unit; The value represents the peak ground acceleration, preferably the area-weighted average of the peak ground acceleration grid values ​​within the evaluation unit. The slope aspect category is assigned a value, which divides the slope aspect into 9 categories: plane, north, northeast, east, southeast, south, southwest, west, and northwest. The original score is preset according to the degree of weathering instability in the glacier area, and then linear normalization is used to map it to the 0-1 interval. The south and southeast sunny slopes have the highest scores, and the plane has the lowest scores. The lithology category is assigned a value, which is divided into five groups according to the rock mass’s resistance to weathering: hard bedrock, relatively hard rock, interbedded soft and hard rock, weak rock, and loose glacial till or colluvial deposits. Loose deposits have the highest original score, while intact hard rocks have the lowest score. The original scores of all groups are linearly scaled to 0-1. and The scores are converted to category scores of 0 to 1 using a preset level assignment method, and then substituted into the TTCI formula for calculation. They can also be represented using one-hot encoding in the input of machine learning models. to Let be the weighting coefficient, and: In a preferred embodiment, to The initial values ​​are taken as 0.12, 0.18, 0.14, 0.20, 0.14, 0.11, and 0.11 respectively. Automatic calibration can be performed by grid search, random search, or Bayesian optimization combined with rolling time verification. This set of weights is used as group g in the unified parameter vector Θ. It is optimized according to the parameter space constraints and group iterative calibration method described in S8 to select the optimal weight combination that maximizes the comprehensive objective function J.

[0099] The formula for calculating the Channel Connectivity Index (CCI) is as follows: ;

[0100] in, This represents the arithmetic mean of the effective pixels in the raster grid within the evaluation cell, representing the distance from the nearest channel. The preferred option represents the average gradient of the main channel within the evaluation unit; The preferred method represents the cumulative flow at the outlet of the evaluation unit or the maximum cumulative flow in the main channel, and a unified statistical caliber is used within the same study area; This represents the proportion of the source area intersecting with the channel buffer zone to the total source area of ​​the evaluation unit, and is used to convert the source channel connection relationship into a calculable numerical feature. to Let be the weighting coefficient, and: In a preferred embodiment, , , , The initial values ​​are 0.35, 0.20, 0.25, and 0.20 respectively. Automatic calibration is performed using grid search, random search, or Bayesian optimization combined with rolling time verification. This set of weights is group c in the unified parameter vector Θ. It is optimized according to the parameter space constraints and group iterative calibration method described in S8 to select the optimal weight combination that maximizes the comprehensive objective function J.

[0101] The methods for obtaining the above four parameters are as follows: A distance grid of the nearest channel is generated using the channel network as the distance source, and partition statistics are performed at the evaluation cell boundaries. The formula for calculating the channel distance is: ;

[0102] Among them, M i This represents the number of effective pixels in the channel distance raster that participate in the statistics within evaluation unit i. Let represent the distance from the r-th valid pixel to the nearest channel; calculate the arithmetic mean of the distances of all valid pixels within the evaluation unit, and convert the distance units to kilometers to obtain . . To evaluate static environmental factors at the unit scale, no time subscript t is set.

[0103] The channel line intersects with the evaluation unit space and is divided into channel segments. The starting and ending elevations and segment lengths of each segment are extracted using a digital elevation model. The slope of each segment is calculated, and the average slope of the main channel within the evaluation unit is obtained by weighting the segments according to their lengths. The calculation formula is as follows: ;

[0104] in, For the evaluation unit i, the first The slope of the main channel segment is calculated by dividing the absolute value of the elevation difference between the starting and ending points of the segment by the length of the segment. For the corresponding line segment length, The number of main channel segments included in the statistics. The channel with the largest cumulative flow within the evaluation unit is preferred as the main channel; if channel level information is lacking, the longest channel within the evaluation unit is used as the main channel.

[0105] The value is obtained by calculating the cumulative flow grid through the digital elevation model in sequence, including filling depressions, flow direction, and cumulative flow. The maximum cumulative flow of the main channel grid in the evaluation unit is selected. When the outlet point of the evaluation unit is known, the cumulative flow at the outlet point can also be taken.

[0106] Generate a preset distance threshold centered on the channel line. The channel buffer zone, source patches, and evaluation units are spatially intersected. The formula for calculating the channel connectivity ratio is as follows: ;

[0107] in, This represents the union of the source patches of evaluation unit i in the previous and next periods within the time interval t. This represents the channel buffer zone centered on the channel line and with a distance threshold of d0. The distance threshold is uniformly determined based on the spatial resolution, channel width, and mapping error of the remote sensing image. A uniform distance threshold is used within the same study area; when the total source area of ​​the evaluation unit is 0, Set the value to 0 and record the area without material source in the quality inspection table.

[0108] The formula for calculating the Regular Source Potential Index (SPI) is as follows: ;

[0109] Comparing different evaluation units within the same time interval Numerical data are used to divide the region into high, medium, and low sediment supply potential zones, and to identify high-risk watersheds for glacial debris flows and glacial till disasters; multiple periods within the same unit can be compared at different time intervals. Changes can be quantitatively identified to determine areas of increasing or decreasing supply potential, supporting dynamic disaster early warning and temporal evolution analysis.

[0110] in, to Let be the weighting coefficient, and: ; In a preferred embodiment, the internal weights of DTI, DASCI, TTCI, and CCI, as well as the combined weights of SPI, can be weighted equally or use the aforementioned preferred initial values, and calibrated according to the unified parameter space described in S8; when the training samples are insufficient to support automatic calibration, a preset weight combination with non-negative weights in each group and a sum of weights within the group of 1 is used to fuse the weights. The optimal value is 0.6. The aforementioned DASCI, TTCI, CCI, DTI, and SPI are not arbitrary weighted sums but correspond to five calculation stages: dynamic accumulation of sediment sources, topographical background, channel connectivity conditions, intertemporal conversion intensity, and comprehensive source potential. The input fields for each indicator are derived from the standardized database formed by S1 to S4; the weights of r, d, g, c, and w are calibrated only within the constraint that they are non-negative and the sum of the weights within each group is 1, and the fusion weights are... Calibration is performed within the range of 0 to 1, based on the high-activity source identification performance in rolling time validation. When an environmental factor is missing or has insufficient accuracy in the study area, the factor can be set to zero according to the preset missing value rules, replaced by statistical values ​​from similar areas, or the corresponding evaluation unit can be included in the quality control checklist.

[0111] S6: Temporal Weakly Supervised Label Generation and Sample Construction.

[0112] Using the source evolution characteristics, state transition characteristics, regularity index characteristics, and environmental factors of the previous time interval as input, weakly supervised labels are generated based on the actual source change results of the next time interval, and training samples of evaluation units multiplied by the time interval are constructed. Figure 5 This is a schematic diagram of time-series weakly supervised sample construction and leakage prevention training. It shows that the source evolution characteristics, regular index characteristics and environmental factors of the previous time interval are used as inputs, and weakly supervised labels are generated based on the real source change results of the next time interval. It also ensures that the change indicators of the next time interval do not enter the model training method of the corresponding sample input features.

[0113] For any evaluation unit i and time interval t, the input features Only data up to and including time interval t is included; the change index at the next time interval t+1 is used only for label generation. , does not enter the input features of the corresponding sample.

[0114] Criteria for determining highly active power supply units: If unit i satisfies any one of the quantile threshold conditions in time period t+1, then the tag... Indicates high activity, otherwise Indicates low activity; when evaluation unit i meets any of the following conditions in the next time interval t+1, it is marked as a high-activity source unit: The total change density of unit i in the current period is greater than or equal to the qth quantile of the total change density of the whole region; it means that the overall change activity of the source of the unit is at the top (1-q) level of the whole region.

[0115] The frequency density of change in unit i during the current period is greater than or equal to the qth quantile of the frequency density of change in the whole region; this indicates that the frequency of evolution of the source patches in this unit is at a high level in the whole region.

[0116] The dynamic transformation intensity index of unit i is greater than or equal to the qth quantile of the DTI of the entire region; it represents that the complexity of the mutual transformation and evolution of source types is significantly higher than that of most regions.

[0117] The unit's net cumulative density is greater than the qth quantile of the total net cumulative density of the region, and the net cumulative density is strictly greater than 0. The total amount of new material sources exceeds the total amount of material sources that have receded, and the loose material sources are continuously replenished.

[0118] Where i represents the evaluation unit to be judged; j represents the evaluation unit traversal sequence number when calculating the quantile threshold; N represents the total number of evaluation units within the target time interval; Let z represent the q-th quantile of the sample set consisting of all evaluation units' corresponding indicator values ​​within the target time interval t+1, where z represents each indicator value. The value of q ranges from 70% to 90%, preferably 80%. If any of the above conditions are met, then the weak supervision label... =1; otherwise =0. In this way, the present invention does not rely on artificial disaster point labels, but automatically constructs a supply trend prediction sample using the real changes in material sources over multiple periods.

[0119] The quantile threshold q for weakly supervised labels is not set as a fixed constant, but is determined based on the distribution of change intensity of evaluation units within the training time window. In this embodiment, the 80th percentile is preferred to achieve a balance between the number of positive samples and the representativeness of highly active units. When the overall change intensity of the study area is low, the sample size is insufficient, or the positive and negative samples are extremely imbalanced, q can be adjusted within the range of 70% to 90%, and the threshold and sample size used are saved in the model training record.

[0120] The time-series weakly supervised sample construction module is invoked, taking the dynamic indicators, state transition features, regular index features, and environmental factors of the previous adjacent time interval as input features, and using the real source change results of the next adjacent time interval to generate weakly supervised labels.

[0121] Specifically, the system predicts highly active source tags for the adjacent time intervals from 2018 to 2020 based on the characteristics of adjacent time intervals from 2016 to 2018; predicts highly active source tags for the adjacent time intervals from 2020 to 2022 based on the characteristics of adjacent time intervals from 2018 to 2020; and predicts highly active source tags for the adjacent time intervals from 2022 to 2024 based on the characteristics of adjacent time intervals from 2020 to 2022. The change indicators of the latter adjacent time interval are only used to generate tags and are not included in the input features of the corresponding samples, thereby avoiding future information leakage.

[0122] In this embodiment, 3571 clean, weakly supervised samples were constructed, of which 712 were positive samples, representing a positive sample ratio of approximately 19.94%.

[0123] S7: Machine learning for supply trend prediction.

[0124] A machine learning binary classification model is trained based on time-series weakly supervised samples. The input is a standardized feature vector for each evaluation unit, which includes at least DASCI, TTCI, CCI, and DTI (details omitted). The model outputs the predicted probability (PML) of each evaluation unit becoming a highly active source unit in the next stage. The PML value ranges from 0 to 1, with a value closer to 1 indicating a higher risk of future high-activity supply. The machine learning classification model includes random forest, gradient boosting decision tree, extreme gradient boosting tree, histogram gradient boosting tree, or other tree ensemble classification models. The model divides the training and validation sets according to time sequence and obtains the PML based on the probability output mechanism of the selected tree ensemble classification model. For random forest, the PML is output based on the positive class voting ratio or the average positive class probability of each decision tree. Gradient boosting models such as gradient boosting decision tree, extreme gradient boosting tree, and histogram gradient boosting tree can map the original predicted scores to a probability range of 0-1 using a sigmoid function. The model loss function is determined based on the selected model. For probabilistic binary classification models using the binary cross-entropy loss function, this function can be used for training. Gradient boosting tree-like models pass the raw output score of the model through the sigmoid function. = Mapping to the interval 0-1 yields the positive class probability, where... Representing the sigmoid activation function, z is the original predicted score output by the gradient boosting tree, and e is the natural constant. The general formula for the binary cross-entropy is: In the formula: is the weak supervision label for the s-th training sample, taking the value 0 or 1; Output the high-activity prediction probability for the model corresponding to the s-th training sample; s represents the total number of samples; s represents the training sample number.

[0125] The model input features include various source area variation characteristics, state transition characteristics, DASCI, TTCI, CCI, DTI, SPI, elevation, slope, topographic relief, distance to active fault, distance to channel, peak ground acceleration, lithology, aspect, dominant source type, and dominant transfer path. The model output is:

[0126] in, denoted by , i represents the probability that evaluation unit i will become a highly active source unit in the next stage based on the input features of time interval t; f represents the machine learning classification model obtained through training.

[0127] Model training employed either rolling time validation or leave-one-out-of-time validation. Evaluation metrics included the area under the receiver operating characteristic (ROC) curve, mean precision, F1 score, hit rate of the top 20% of high-potential units, recall rate of the top 20% of high-potential units, and recognition rate of highly active units.

[0128] The machine learning supply trend prediction module is invoked to train a tree ensemble classification model based on the aforementioned time-series weakly supervised samples. This embodiment employs a random forest model and uses a feature combination composed of source change characteristics and environmental factors.

[0129] The model outputs the probability that each evaluated watershed will become a highly active water source unit in the next stage, denoted as the Machine Learning Highly Active Water Source Probability (PML). Subsequently, the automatic weight calibration and index fusion module is invoked to fuse the regular water source potential index (SPI) with the PML, resulting in the Machine Learning Calibrated Water Source Potential Index (ML-SPI). The specific meanings, value ranges, and calibration methods of each group of weights to be calibrated and the fused weights are described below.

[0130] S8: Automatic weight calibration and machine learning calibration of the source potential index calculation.

[0131] By fusing the regularized source potential index SPI with the machine learning output probability PML, a machine learning-calibrated source potential index ML-SPI is obtained, calculated as follows: in, The fusion weights represent the machine learning probabilities PML in the ML-SPI calculation, with mathematical values ​​ranging from 0 to 1, preferably within the range of 0.3 to 0.8. For each evaluation task, the optimal fusion weights are determined during the unified parameter calibration process described in S8, based on the rolling time validation results of that task. This method is adopted uniformly across all evaluation units in this evaluation task. This embodiment, after calibration, takes... =0.6; for other study areas or evaluation tasks, the appropriate value can be determined by recalibration. When the sample size is insufficient to support automatic calibration, 0.6 is preferred as the initial value or default value.

[0132] To address the issue of manually setting index weights in traditional source potential assessments, this invention applies rolling time validation results to the internal weights, combined weights, and fusion weights of DTI, DASCI, TTCI, and CCI. Automatic calibration is performed. To ensure consistency between the optimization variables and the calculation formulas of each index, all the weights to be calibrated are grouped into a unified parameter vector, and the comprehensive objective function J is used as the unified calibration basis within the rolling-time verification framework. The unified parameter vector is defined as follows: Among them, group r includes to Group d includes to , Group includes to Group C includes to Group w includes to , The fusion weights are determined by the following formula: Parameter space Ω constraints r, d, The weights of groups c and w are all not less than 0 and the sum of the weights within each group is 1; The mathematical value range is 0 to 1, with an optimal search range of 0.3 to 0.8. The comprehensive objective function is: Where J represents the comprehensive optimization objective, This represents the optimal parameter vector obtained from calibration, where Ω represents the parameter space. The values ​​represent the maximum value, P20 represents the hit rate of the top 20% of high-potential units, AUC represents the area under the receiver operating characteristic (ROC) curve, and AP represents the average precision. AUC, or the area under the ROC curve, is a general discrimination index for binary classification models. It is generated by iterating through the ROC curve with the ML-SPI value corresponding to the candidate parameter combination as the threshold, and the area under the curve is obtained by integration. It can be directly solved by a standard machine learning package, and the value ranges from 0 to 1. The larger the value, the stronger the ability of the fusion index to distinguish between high and low active source units. AP, or average precision, is the area under the integral of the Precision-Recall curve. It is adapted to the sample characteristics of imbalanced positive and negative samples in this invention and quantifies the accuracy of the fusion index in identifying high-active source units. P20, or the hit rate of the top 20% of high-potential units, is calculated by sorting all evaluation units in the validation set according to the ML-SPI calculated by the candidate parameter combination from high to low, selecting the top 20% of units, and counting the proportion of units with weak supervision label Y=1. It is used to measure the screening effect of high-active source units in the fusion result. , , Let represent the weights of the objective function, and: In a preferred embodiment, Take 0.5, Take 0.3, The value is set to 0.2; where α, β, and γ represent only the weights of the three fusion verification metrics P20, AUC, and AP in the comprehensive objective function; for ease of implementation, a grouped iterative calibration method is preferred: with other parameters fixed, r, d, ... The weights within each group (c) are optimized in one round, and then the weights of group (w) are jointly optimized. Repeat the iterations until the overall objective function J no longer increases or reaches the preset number of iterations, then output the result. After each update of the internal weights, the corresponding index and SPI are recalculated. When the corresponding index is used as a model input feature, the model is retrained and the PML is updated using the same time division, and then the ML-SPI and J are calculated. Each optimization can be solved using grid search, random search, or Bayesian optimization, or a joint solution can be performed using random search or Bayesian optimization within the parameter space Ω.

[0133] Furthermore, in the same evaluation task, the optimal parameter vector determined through rolling time validation... The same approach is adopted across all evaluation units without any unit-level adjustments. For evaluation units with missing data, density anomalies, small watershed anomalies, or discrepancies in model rules, only quality control markers are made and they are added to the review checklist; no changes are made. The value of .

[0134] S9: Supply potential classification, quality inspection and verification and output of results.

[0135] The ML-SPI resource potential is classified into levels. The classification employs the equal-interval method, the natural breakpoint method, or the quantile method. Preferably, the non-uniform quantile method based on percentile rank is used. The ML-SPI values ​​of all evaluation units within the study area are sorted in ascending order and divided into four levels according to percentile rank: low, medium, relatively high, and high. Evaluation units with a percentile rank not greater than 50% are considered low resource potential; evaluation units with a percentile rank greater than 50% but not greater than 80% are considered medium resource potential; evaluation units with a percentile rank greater than 80% but not greater than 90% are considered relatively high resource potential; and evaluation units with a percentile rank greater than 90% are considered high resource potential. For evaluation units with the same ML-SPI value, the average rank is used to determine the percentile rank, and evaluation units at the classification boundary are assigned to the lower level corresponding to the boundary. Considering that the source potential classification is mainly used to screen a small number of high potential evaluation units that need to be prioritized for monitoring and verification, the above non-proportional classification is aligned with the screening rules for the top 20% of high potential units, and the top 20% of evaluation units are further divided into two levels: higher and higher, in order to enhance the ability to distinguish the evaluation units at the tail end of the source potential distribution. Figure 7 This is a schematic diagram of the output and quality inspection review, showing the process by which the system outputs a source potential classification map, a list of high-potential source units, an abnormal unit review table, model interpretation results, and a regional evaluation report based on the ML-SPI evaluation results, dominant source types, dominant transfer paths, model prediction probabilities, and quality inspection marks. Figure 8 This is a schematic diagram of the evaluation results and model verification of the embodiment, showing the number of ML-SPI source potential levels, PML probability distribution, rolling verification performance of the model based on PML, the top 10 feature importance items, and the ML-SPI source potential spatial classification results in the embodiment, which is used to illustrate the application effect of the present invention in a specific study area.

[0136] The system further overlays the dominant source type, dominant transfer path, model prediction probability, and feature contribution interpretation results, and outputs a list of high-potential supply units, a supply potential classification map, a machine learning probability map, a dominant source type map, a dominant transfer path map, a dynamic conversion intensity index spatial distribution map, a model interpretation map, and an abnormal unit quality inspection map.

[0137] The system performs quality control and verification marking on the evaluation results. Quality control includes checking for missing periods, low-confidence identified patches, small-area abnormal evaluation units, extreme density variations, abnormal state transitions, and significant discrepancies between model predicted probabilities and rule indices. For evaluation units that trigger quality control rules, the system generates a list of abnormal units, recording the cause of the anomaly, the corresponding time interval, the type of material involved, and the recommended verification method.

[0138] The system retrieves the previously generated study area boundary, multi-period source identification results, state transition statistics, model prediction probabilities, source potential classification indicators, and full quality inspection records. Based on a preset report template, it automatically fills in text, statistical values, and spatial charts to form a regional evaluation report. The system further generates a regional evaluation report based on the study area boundary, multi-period source identification results, state transition results, model prediction results, source potential classification results, and quality inspection records. This report includes a study area overview, data source description, source identification and standardization results, inter-period source evolution characteristics, main state transition paths, distribution of high-potential source units, dominant control factors, units requiring verification, and disaster prevention and monitoring recommendations.

[0139] The system utilizes modules for source potential classification, quality inspection and verification, and results output to sort and classify the ML-SPI of 1192 evaluated watersheds, generating four levels of source potential results: low, medium, relatively high, and high. The system further identifies the top 20% of high-potential watersheds according to the ML-SPI ranking and, combined with data quality, model prediction results, rule index results, and anomaly markers, generates a reliable list of high-potential watersheds.

[0140] The system synchronously outputs source potential classification results, a list of high-potential source units, a spatial result attribute table, a comprehensive quality inspection checklist, and a regional evaluation report. The regional evaluation report includes a study area and data overview, multi-period source standardization results, inter-period state transition and dynamic evolution results, time-series weakly supervised prediction results, ML-SPI source potential evaluation results, high-potential source unit identification results, quality inspection review results, and disaster prevention monitoring recommendations. The system further generates... Figure 8 The illustrated example shows the evaluation results and model verification diagram.

[0141] The following is an explanation of the quality inspection and reliability control of the results after the dynamic source evaluation results of the present invention are generated: The quality inspection and review objects of the method of this invention include patch-level objects, evaluation unit-level objects, and outcome-level objects. Patch-level objects include low-confidence patches, patches with abnormal boundary confidence, small-area patches, overlapping patches of different categories in the same year, and patches requiring review in a given period. Evaluation unit-level objects include evaluation units with excessively small areas, evaluation units with abnormal change density, evaluation units with abnormal state transition area closure, and evaluation units with abnormal state transition paths. Outcome-level objects include evaluation watersheds that are among the top 20% of high-potential lists but have data quality issues, as well as evaluation watersheds where the regular source potential index (SPI) and the machine learning high-activity source probability (PML) show significant discrepancies.

[0142] The input data in this embodiment includes a standard source patch library, a 9×9 source state transition matrix, a dynamic index table, PML prediction results, SPI results, ML-SPI results, a list of high-potential source units, and a spatial result attribute table.

[0143] In this embodiment, the small-area patch verification threshold is preferably 0.001 square kilometers; the evaluation unit area anomaly threshold is determined based on the watershed area distribution of the study area; and the state transition area closure error threshold is set based on the evaluation unit area and the spatial overlay calculation accuracy.

[0144] When the transfer area ratio RATIO is greater than 1 and the difference between RATIO and 1 is less than 1× When the system determines that the error is a floating-point calculation tail difference, it is marked as a low-level quality inspection item; when the difference between RATIO and 1 is greater than or equal to 1× When the system determines that the area ratio that needs to be reviewed is abnormal, it determines the medium- or high-level quality inspection item based on the abnormal area and whether it is included in the top 20% high-potential list.

[0145] The significant discrepancy can be determined based on the normalized difference between PML and SPI, the difference in quantile levels, or the difference in source potential levels. Preferably, when the quantile levels of the two differ by two or more levels, or the normalized difference exceeds a preset threshold, the system determines that there is a discrepancy in the model rules and a review is required.

[0146] The quality inspection level is divided into three levels: low, medium, and high, determined comprehensively based on the anomaly type, degree, scope of impact, whether it is included in the top 20% high-potential list, and whether it may affect the judgment of the supply potential. The same quality inspection type can be divided into different levels under different influence conditions. Low-level quality inspection items mainly correspond to floating-point tail differences, small-area problems that are not included in the top 20% high-potential list and do not affect the overall interpretation, or other local anomalies that do not affect the overall interpretation; medium-level quality inspection items mainly correspond to local overlap, local density anomalies, or boundary problems that require manual verification; high-level quality inspection items mainly correspond to large-area overlap that may affect the high-potential judgment, extreme density, state transition anomalies, significant discrepancies between rule indices and machine learning probabilities, and small-area evaluation units that are included in the top 20% high-potential list or may affect the high-potential judgment due to the area denominator effect. The quality inspection review mainly includes six steps: patch-level quality inspection review, state transition matrix quality inspection review, dynamic indicator and evaluation unit quality inspection review, model and rule result discrepancy check, comprehensive quality inspection list generation, result reliability marking, and report output. Patch-level quality inspection and review: The system performs integrity checks and spatial consistency checks on each patch in the standard source patch library. Integrity checks include checking for missing years, abnormal category codes, zero or negative areas, identification confidence levels below preset thresholds, and boundary confidence levels below preset thresholds. Spatial consistency checks include checking for overlap between different categories of patches from the same year, frequent category jumps between adjacent periods, and significant category changes within the same region that do not conform to inter-period evolution logic. For patches that trigger quality inspection rules, the system records the patch number, year, source type, quality inspection type, quality inspection level, and review recommendations.

[0147] State transition matrix quality inspection and review: The system performs area closure checks, state coding checks, and transition path checks on the 9×9 state transition matrix for each evaluation watershed and each adjacent time interval. The area closure check determines whether the sum of all transition areas in the 9×9 matrix within each evaluation watershed matches the area of ​​that watershed. The state coding check determines whether there are cases where the previous or next state is empty, the state code exceeds the range of 0 to 8, or the category name is unrecognizable. The transition path check identifies abnormal category transition paths, frequent category jump paths, and abnormal paths in the disaster chain direction. For issues such as negative area records, area closure errors, RATIO greater than 1, and overlap of different categories of patches in the same year, the system assigns low, medium, or high-level quality inspection marks and writes them into the state transition quality inspection checklist.

[0148] Dynamic Indicators and Evaluation Unit Quality Inspection and Review: The system checks the attributes of dynamic indicators and evaluation units. The inspection includes checking for outliers in total change density, change frequency density, category transformation ratio, disaster chain direction transformation ratio, state transition entropy, whether the evaluated watershed area is too small, and whether there are missing areas, negative increases or decreases in area, abnormal numbers of change patches, or total change area exceeding a reasonable range. For small-area evaluation units, high change density evaluation units, extreme frequency density evaluation units, and evaluation units with abnormal indicator calculations, the system marks them as evaluation units requiring review and records the reasons for the anomalies.

[0149] Model and rule-based result discrepancy check: The system performs consistency analysis on the results of PML, SPI, and ML-SPI. If a given watershed has a high PML but a low SPI, it indicates that the model has identified a potential high-activity trend from historical evolution samples, but the regular geoscientific index is insufficient to support it. If a given watershed has a high SPI but a low PML, it indicates that the given watershed has a strong topographic structure or channel connectivity background, but a high-activity trend has not yet been shown in historical evolution samples. These situations are not directly judged as errors, but rather marked as units requiring interpretation and review. For given watersheds that simultaneously enter the top 20% high-potential list and have significant model-rule discrepancies, the system prioritizes them for review.

[0150] Comprehensive Quality Inspection Checklist Generation: The system integrates patch-level quality inspection, state transition quality inspection, dynamic indicator quality inspection, model rule divergence quality inspection, and high-potential list review results to form a unified comprehensive quality inspection checklist. The comprehensive quality inspection checklist includes at least the following fields: quality inspection number, quality inspection source, quality inspection type, time interval, year, evaluation watershed number, relevant sediment source type, relevant transfer path, abnormal area, quality inspection level, and review recommendation. The system sorts the quality inspection objects according to their quality inspection level and whether they are included in the top 20% high-potential list. Evaluation watersheds that are included in the top 20% high-potential list and trigger high-level quality inspection rules are prioritized at the top of the review list; evaluation watersheds with only floating-point tail differences or small-area issues are considered low-priority review objects. The comprehensive quality inspection checklist uses quality inspection rule trigger records as the statistical object; when the same evaluation unit triggers different modules or different quality inspection rules simultaneously, separate quality inspection records can be generated, therefore the number of quality inspection records is not equal to the number of deduplicated evaluation units.

[0151] Results Reliability Marking and Report Output: The system integrates comprehensive quality inspection results into the source potential classification map, PML probability map, ML-SPI spatial classification map, high-potential source unit list, and regional assessment report. For high-potential assessment watersheds without obvious quality inspection problems, the system marks them as reliable high-potential units; for high-potential assessment watersheds with high-level quality inspection problems, the system marks them as priority review units; for assessment watersheds with medium- or low-level problems, the system marks them as general review units or suggestive review units. The regional assessment report simultaneously outputs the study area and data overview, multi-period source standardization results, inter-period state transition and dynamic evolution results, time-series weak supervision prediction results, ML-SPI source potential assessment results, high-potential source unit identification results, quality inspection review results, and disaster prevention monitoring recommendations.

[0152] The quality inspection review results of a specific embodiment of the present invention show that the comprehensive quality inspection checklist records a total of 596 quality inspection items, including 97 high-level quality inspection items, 2 medium-level quality inspection items, and 497 low-level quality inspection items.

[0153] In terms of quality inspection types, the comprehensive quality inspection checklist includes 553 small-area evaluation unit reviews, 13 floating-point tail difference records with a transfer area ratio RATIO greater than 1, 8 negative area anomalies, 5 annual statistical records of overlapping patches of different categories in the same year, 4 state transfer area closure errors, 4 watershed area closure maximum errors, 4 records of change density greater than 1, 2 high potential reliability reviews of the top 20%, 1 missing start state code or end state code, 1 unidentified category, and 1 time interval with abnormally high category transformation area, totaling 596 items.

[0154] The small-area evaluation unit review records include 473 low-level and 80 high-level records; the records of overlapping patches of different categories in the same year include 2 medium-level and 3 high-level records; 8 negative area anomalies, 4 records with change density greater than 1, and 2 records of reliability review for the top 20% high potential are all high-level, and the remaining quality inspection types are low-level, thus summarizing 97 high-level, 2 medium-level, and 497 low-level records. The above 80 high-level small-area review records correspond to 40 small-area evaluation units that simultaneously trigger the ML-SPI fusion evaluation result quality inspection rules and the top 20% high potential list review rules, with each of the two quality inspection sources forming 40 records. In addition, overlapping conflicts of different categories of patches in the same year are recorded as issues requiring review, and the area of ​​the original overlapping conflict records after annual summarization is approximately 724.464 square kilometers; negative area records are treated as quality inspection items and are not included as the main evaluation results. In this embodiment, there are 239 top 20% high potential watersheds, of which 237 are reliable high potential watersheds. For evaluation watersheds that are included in the top 20% high potential list and trigger high-level quality inspection rules, the system will list them as priority review targets; for high potential watersheds that have not triggered high-level quality inspection rules, the system will mark them as reliable high potential units.

[0155] The system integrates the comprehensive quality inspection results into the source potential classification map, PML probability map, ML-SPI spatial classification map, high-potential source unit list, and regional evaluation report. Therefore, the quality inspection and review mechanism of this invention can identify issues requiring review at different levels, including source patches, state transition matrices, dynamic indicators, model fusion results, and output results. It can also incorporate the quality inspection type, quality inspection level, and review recommendations into the output list and regional evaluation report, thereby improving the traceability of dynamic source evaluation results and the reliability of engineering applications.

[0156] In one specific embodiment of the present invention, a 9×9 state transition matrix containing states from 0 to 0 is generated for each of the four adjacent time intervals. The state transition statistics are as follows: Between 2016 and 2018, the area of ​​newly added material sources was 1296.545 square kilometers, the area of ​​material source regression was 1050.672 square kilometers, the area of ​​similar continuity was 2863.989 square kilometers, and the area of ​​category transformation was 480.996 square kilometers; between 2018 and 2020, the area of ​​newly added material sources was 1055.798 square kilometers, the area of ​​material source regression was 1894.748 square kilometers, the area of ​​similar continuity was 2267.831 square kilometers, and the area of ​​category transformation was 449.9 square kilometers. 87 square kilometers; During the adjacent time interval from 2020 to 2022, the newly added source area was 2068.717 square kilometers, the source decline area was 865.438 square kilometers, the area of ​​similar continuity was 2508.094 square kilometers, and the category transformation area was 534.582 square kilometers; During the adjacent time interval from 2022 to 2024, the newly added source area was 915.144 square kilometers, the source decline area was 1107.234 square kilometers, the area of ​​similar continuity was 3643.001 square kilometers, and the category transformation area was 542.164 square kilometers.

[0157] In the state transition matrix statistics, the largest category transition path in terms of area is the glacial till to glacier path within adjacent time intervals from 2016 to 2018, covering an area of ​​approximately 242.584 square kilometers. This result, as a statistical output of the state transition matrix, is not directly equivalent to a single real geomorphic process. Interpretation of the results requires consideration of patch boundaries, classification differences, overlapping conflicts between different categories of patches in the same year, and manual verification results.

[0158] The rolling-time validation results based on PML prediction probabilities show that, as Figure 8 In (d), the ROC-AUC (area under the receiver operating characteristic curve) of the random forest model is 0.901, AP (mean precision) is 0.677, and the accuracy of the top 20% of high-potential units sorted by PML is 0.660. Figure 8 As shown in (b), based on 1192 evaluation watersheds with adjacent time intervals from 2022 to 2024, the maximum value of PML was calculated to be approximately 0.981 and the mean value was approximately 0.199; the maximum value of ML-SPI was approximately 0.789 and the mean value was approximately 0.223.

[0159] like Figure 8As shown in (a), the ML-SPI classification results indicate that there are 596 watersheds with low supply potential, 357 with medium supply potential, 119 with relatively high supply potential, and 120 with high supply potential. The top 20% of high-potential watersheds identified by ML-SPI ranking are 239, of which 237 are reliable high-potential watersheds. Reliable high-potential watersheds refer to those included in the top 20% high-potential list that have not triggered high-level quality control rules.

[0160] The quality inspection review results show that this embodiment generated a comprehensive quality inspection checklist of 596 items, including 97 high-level quality inspection items, 2 medium-level quality inspection items, and 497 low-level quality inspection items. The complete composition, level classification, and recording scope of multiple rule triggers for the same evaluation unit are consistent with the aforementioned quality inspection review results. Among them, overlapping conflicts of different types of patches in the same year are recorded as issues requiring review. The area of ​​the original overlapping conflict records after annual summarization is approximately 724.464 square kilometers. Negative area records are treated as quality inspection items and are not included as the main evaluation results. The above quality inspection items indicate that this embodiment retains data issues in the actual spatial overlay process, rather than treating source patches and evaluation units as idealized examples. For evaluation watersheds that are included in the top 20% high potential list and also have high-level quality inspection items, the system does not directly give a definitive conclusion, but instead treats them as priority review objects and requires verification in conjunction with remote sensing imagery, patch boundaries, and necessary field survey data.

[0161] As can be seen from the specific embodiments of this invention, this invention can unify the standardization of multi-period source patches, construction of a 9×9 state transition matrix, calculation of dynamic indicators, time-series weakly supervised prediction, ML-SPI fusion evaluation, source potential classification, quality inspection and verification, and regional evaluation report generation into a single technical process, thereby obtaining source potential classification results at the watershed scale, a list of high-potential source units, anomaly verification prompts, and a regional evaluation report; simultaneously, as Figure 8 As shown in (c), the system outputs the top 10 most important model features, which are used to help explain the main contributing factors of the PML prediction results.

[0162] This invention also discloses a system for identifying and dynamically evaluating the source of materials in glacier disaster chains. The system includes a data access and preprocessing module, a source identification result access and standardization module, a cross-period consistency correction and state transition analysis module, an evaluation unit construction and dynamic index statistics module, a rule index calculation module, a time-series weakly supervised sample construction module, a machine learning source trend prediction module, an automatic weight calibration and index fusion module, a quality inspection and verification module, and a result output module. Figure 2The system module structure diagram shows the composition and relationship of the following modules: data access and preprocessing module, material source identification result access and standardization module, inter-period consistency correction and state transition analysis module, evaluation unit construction and dynamic indicator statistics module, rule index calculation module, time series weakly supervised sample construction module, machine learning supply trend prediction module, weight automatic calibration and index fusion module, quality inspection and verification module, and result output module.

[0163] The data access and preprocessing module is the system execution module. It outputs a basic spatial dataset after unified projection, cropping, and repair, which serves as the basic input for all downstream modules. The data access and preprocessing module is used to acquire the study area boundary, multi-period remote sensing images, multi-period manually interpreted source patches, existing source labeling data, semantic segmentation model output results, and auxiliary geospatial data. It also completes unified coordinate projection, spatial cropping, geometric repair, topology checking, field standardization, time period matching, and quality checking.

[0164] The provenance identification result access and standardization module reads the preprocessed multi-source provenance data, generates a standardized provenance patch library with eight unified codes, and outputs it to the intertemporal state transition analysis module. The provenance identification result access and standardization module is used to receive manual interpretation results, existing provenance labeled data, or semantic segmentation model output results, and uniformly map them into eight provenance types: glaciers, glacial lakes, glacial till, rock wall deposits, glacial hazardous rock masses, landslides, gully deposits, and depositional fans. It also records the category, area, identification confidence, boundary confidence, and data source.

[0165] The intertemporal consistency correction and state transition analysis module constructs a 9×9 state transition matrix based on a standardized source patch library, outputs source temporal evolution data to the evaluation unit statistics module, and performs spatial overlay and category comparison of source patches in adjacent periods to identify newly added, declining, persistent, migrating and category-transformed patches, and constructs a 9×9 source state transition matrix at the evaluation unit level, outputting the dominant transition path and abnormal transformation markers.

[0166] The evaluation unit construction and dynamic indicator statistics module overlays watershed unit and state transition data to calculate a full range of dynamic indicators of material sources. The indicator library is simultaneously supplied to the two major modules of regular index calculation and weakly supervised sample construction. The evaluation unit construction and dynamic indicator statistics module is used to generate watershed, sub-watershed, valley unit or regular grid evaluation units, calculate the environmental factors of the evaluation units, and statistically analyze dynamic indicators such as net cumulative area, total changed area, total changed density, change frequency density, persistent area ratio, category transformation ratio, disaster chain direction transformation ratio, state transition entropy, dominant material source type and dominant transfer path.

[0167] The regular index calculation module reads dynamic indicators and topographic environmental factors, calculates the DTI, DASCI, TTCI, CCI, and SPI regular indices, and stores the index results in the shared index library. The regular index calculation module is used to calculate the dynamic conversion intensity index DTI, dynamic material accumulation index DASCI, topographic structure coupling index TTCI, channel connectivity index CCI, and regular source potential index SPI, and performs missing value processing, quantile truncation normalization, negative factor transformation, and outlier field marking on the indicators.

[0168] The time-series weakly supervised sample construction module retrieves historical time-series dynamic indicators, rule indices, and environmental factors to construct time-series training samples, which are then input into the machine learning prediction module. The time-series weakly supervised sample construction module is used to generate weakly supervised labels by taking the features of the previous time interval as input and the real source change results of the next time interval, constructing training samples of evaluation units multiplied by the time interval, and avoiding the input features of the change indicators of the next time interval into the corresponding sample input features.

[0169] The machine learning source trend prediction module trains the model based on weakly supervised samples and outputs the high-activity source probability PML of each evaluation unit. The PML and SPI are sent to the weight calibration fusion module. The machine learning source trend prediction module is used to train random forest, gradient boosting decision tree, extreme gradient boosting tree, histogram gradient boosting tree or other tree ensemble classification model, outputs the probability PML of the evaluation unit becoming a high-activity source unit in the next stage, and evaluates the model performance through rolling time verification or time window leave-one-out verification.

[0170] The automatic weight calibration and index fusion module automatically optimizes the internal weights of DTI, DASCI, TTCI, and CCI, the combined weights of SPI, and the fused weights based on rolling time validation, and calculates the ML-SPI machine learning calibration source potential index. This module is used to adjust the weights of DTI, DASCI, TTCI, and CCI, the combined weights of SPI, and the fused weights based on the rolling time validation results. Automatic calibration is performed, and the regular source potential index SPI is fused with the machine learning probability PML to obtain the machine learning calibrated source potential index ML-SPI.

[0171] The quality inspection and review module is a full-process post-verification module. It reads all intermediate results, including standardized patches, state transition matrix, dynamic indicators, SPI, PML, and ML-SPI, identifies various anomalies according to preset thresholds, classifies and marks them, and generates a quality inspection checklist. The quality inspection and review module is used to identify issues such as missing periods, low-confidence patches, small-area anomaly evaluation units, extreme change density, state transition anomalies, and excessive differences between model prediction probabilities and regular indices. It also generates a list of anomaly units, the causes of anomalies, and suggested review methods.

[0172] The output module uniformly reads ML-SPI grading results, PML probabilities, source types, transfer paths, DTI spatial indicators, and quality inspection checklists, and generates various spatial maps, high-potential unit lists, and regional evaluation reports in batches. The first eight modules constitute the main sequential time-series calculation process and must be executed in order. The quality inspection and verification module and the output module are post-processing summary modules, reading the intermediate calculation results from all the aforementioned modules to complete verification and result export. All modules rely on a unified standardized spatial database and dynamic indicator library for data interaction. The output module is used to output source potential grading maps, machine learning probability maps, dominant source type maps, dominant transfer path maps, dynamic conversion intensity index spatial distribution maps, model interpretation maps, abnormal unit quality inspection maps, high-potential source unit lists, and regional evaluation report results.

[0173] In another embodiment of the present invention, the evaluation unit can be replaced by a watershed with a sub-watershed, a valley unit, or a regular grid; the source identification result can be provided by manual interpretation results, existing labeled data, or the output result of a semantic segmentation model; the machine learning model can be replaced by a random forest with a gradient boosting decision tree, an extreme gradient boosting tree, or a histogram gradient boosting tree; the source potential classification method can adopt the quantile method, the natural breakpoint method, or the equal interval method.

[0174] The above embodiments are merely illustrative examples for clear explanation and are not intended to limit the implementation. Those skilled in the art will recognize that other variations or modifications can be made based on the above description. It is neither necessary nor possible to exhaustively list all possible implementations. However, obvious variations or modifications derived therefrom are still within the scope of protection of this invention.

Claims

1. A method for identifying the source material of glacier disaster chains and evaluating their dynamic supply, characterized in that, Includes the following steps: S1. Acquire multi-source data of remote sensing, source patches, and topographic structure of the study area from multiple periods, complete preprocessing, and generate a basic spatial dataset; S2. Read the basic spatial dataset, complete the unified category mapping, and generate a standardized source patch library; S3. Using the standardized source patch library as the data source, construct the source state transition matrix by spatially superimposing adjacent patches, and output the time-series evolution statistics for dynamic index calculation. S4. Overlay the source state transition matrix with the evaluation unit layer and statistically analyze the dynamic indicators of multiple sources. S5. Using the dynamic indicators from step S4 combined with topographic factors, solve for four indicators: Dynamic Transformation Intensity Index (DTI), Dynamic Source Accumulation Index (DASCI), Topographic Structure Coupling Index (TTCI), and Channel Connectivity Index (CCI). The weighted sum of these four indicators yields the Regular Source Potential Index (SPI). S6. Extract the source dynamic indicators, DTI, DASCI, TTCI, CCI and SPI time series quantization data from steps S4 and S5 as model input, and use the source change results of the next period to generate unlabeled weakly supervised samples, which are then sent to machine learning training. S7. Based on the unannotated weakly supervised samples, train the tree ensemble model and output the high-activity source probability (PML) of each evaluation unit. S8. Substitute the rule-based source potential index (SPI) and the high-activity source probability (PML) of each evaluation unit into the fusion formula to obtain the ML-SPI. S9. Multi-level power supply potential classification based on ML-SPI.

2. The method for identifying and dynamically evaluating the source of glacier disaster chains according to claim 1, characterized in that: The S3 step constructs a 9×9 source state transition matrix, with rows and columns corresponding to the source states of two adjacent time periods. State code 0 represents a region without source, and codes 1 to 8 correspond to eight types of standardized glacial hazard chain source in sequence. The value of each element of the matrix is ​​the proportion of the area of ​​the corresponding state transition type to the total area of ​​the evaluation unit.

3. The method for identifying and dynamically evaluating the source of glacier disaster chains according to claim 1, characterized in that, The Dynamic Conversion Intensity Index (DTI) is calculated using the following method: ； in, The index represents the dynamic conversion intensity of evaluation unit i within time interval t; N represents the normalization function. to These are the weighting coefficients; This indicates the proportion of newly added material sources in evaluation unit i within time interval t; Indicates the category conversion ratio; Indicates the proportion of disaster chain direction transformation; Represents the state transition entropy; Indicates the percentage of continuous area; The method for calculating the Dynamic Source Accumulation Index (DASCI) is as follows: ； in: This represents the dynamic material source accumulation index of evaluation unit i within time interval t; Indicates net cumulative density; Indicates the total density of change; Indicates the frequency density of change; Indicates the contribution rate of the dominant material source; to These are the weighting coefficients; The method for calculating the terrain structure coupling index (TTCI) is as follows: ； in: This represents the terrain structure coupling index of evaluation unit i; Indicates the average elevation of the evaluation unit; Indicates the average slope; Indicates the degree of terrain relief; Indicates the distance from the active fault; Indicates peak ground acceleration; This indicates the assignment of slope aspect category; Indicates the lithology category assignment; to These are the weighting coefficients. This represents the negative normalization function. ; The method for calculating the Channel Connectivity Index (CCI) is as follows: ； in, This represents the channel connectivity index of evaluation unit i within time interval t; This represents the arithmetic mean of the effective pixels in the raster grid within the evaluation cell, representing the distance from the nearest channel. This represents the average gradient of the main channel within the evaluation unit; This indicates the cumulative flow at the outlet of the evaluation unit or the maximum cumulative flow in the main channel; This indicates the proportion of the source area intersecting with the channel buffer zone to the total source area of ​​the evaluation unit; to These are the weighting coefficients.

4. The method for identifying and dynamically evaluating the source of glacier disaster chains according to claim 1, characterized in that, The method for calculating the rule-based source potential index (SPI) is as follows: ； in: This represents the rule-based supply potential index of evaluation unit i within time interval t; This represents the dynamic material source accumulation index of evaluation unit i within time interval t; This represents the terrain structure coupling index of evaluation unit i; This represents the channel connectivity index of evaluation unit i within time interval t; This represents the dynamic conversion intensity index of evaluation unit i within time interval t; to These are the weighting coefficients.

5. The method for identifying and dynamically evaluating the source of glacier disaster chains according to claim 1, characterized in that, The method for calculating the high-activity source probability PML of each evaluation unit is as follows: based on the unlabeled weakly supervised samples constructed in step S6, a tree ensemble classification model is selected as the prediction model, and the source dynamics index, topographic environmental factors, DTI, DASCI, TTCI, CCI, and SPI time-series quantitative data corresponding to each evaluation unit are used as input features. The model divides the training set and validation set according to the time sequence, and obtains the Probability Output Model (PML) based on the probability output mechanism of the selected tree ensemble classification model; among them, the random forest outputs the PML based on the positive class voting ratio or the average positive class probability of each decision tree.

6. The method for identifying and dynamically evaluating the source of glacier disaster chains according to claim 1, characterized in that, The ML-SPI calculation method is as follows: ； in, This represents the ML-SPI value of evaluation unit i within time interval t; This represents the probabilistic fusion weights for machine learning, with values ​​ranging from 0 to 1; This represents the rule-based supply potential index of evaluation unit i within time interval t; This represents the probability of a highly active source in evaluation unit i within time interval t.

7. A system for identifying and dynamically evaluating the source of materials in glacier disaster chains, characterized in that: It includes a data access and preprocessing module, outputs a unified spatial dataset, and a downstream module for accessing and standardizing object source identification results; The source identification results are accessed and standardized by the module, which reads the spatial dataset to generate a standard source patch library and outputs it to the intertemporal consistency correction and state transition analysis module. The intertemporal consistency correction and state transition analysis module constructs a state transition matrix based on the patch library, and the results are input into the evaluation unit construction and dynamic index statistics module. The evaluation unit construction and dynamic index statistics module calculates a complete set of material source dynamic indicators and transmits them to the rule index calculation module. The rule index calculation module generates the SPI index and simultaneously supplies it to the time-series weakly supervised sample construction module and the automatic weight calibration and index fusion module. The time-series weakly supervised sample construction module generates training samples based on time-series indicators and outputs them to the machine learning supply trend prediction module. The machine learning supply trend prediction module trains the model to output PML, which is then passed to the weight automatic calibration and exponential fusion module. The automatic weight calibration and index fusion module receives SPI and PML, calibrates the internal weights, SPI combined weights, and fusion weights of DTI, DASCI, TTCI, and CCI based on rolling time verification, calculates ML-SPI, and pushes them to the quality inspection and review module and the result output module, respectively. The quality inspection and review module is used to read intermediate data from all the aforementioned modules to identify various anomalies; The output module is used to integrate all data such as ML-SPI and quality inspection marks to output various results.

8. The glacier disaster chain source identification and dynamic supply evaluation system according to claim 7, characterized in that: The data access and preprocessing module outputs a unified spatial dataset to the material source identification result access and standardization module, which then transmits the state matrix, dynamic indicators, and rule indices layer by layer to the machine learning material source trend prediction module. The SPI and PML synchronously input weight automatic calibration and index fusion module generates ML-SPI.

9. The glacier disaster chain source identification and dynamic supply evaluation system according to claim 7, characterized in that: The quality inspection and review module reads source patches, state matrix, SPI, and PML data, identifies missing periods, low-confidence patches, extreme values ​​of indicators, and index grading divergences, classifies them into three levels of quality inspection, and outputs a standardized list of anomalies.

10. The glacier disaster chain source identification and dynamic source evaluation system according to claim 7, characterized in that: The output module integrates ML-SPI classification, PML probability, dominant source material, transfer path, and quality inspection marks to generate multiple spatial classification maps, a list of high-potential units, and a regional assessment report containing disaster prevention and monitoring recommendations.

Images