Rainfall prediction method based on SHAP meteorological characteristic analysis

Abstract

The invention relates to the technical field of rainfall prediction, and discloses a rainfall prediction method based on SHAP meteorological feature analysis, which performs feature screening by using an SHAP framework to improve model efficiency, and deeply integrates interpretability analysis into a whole process of model construction and business decision. According to the method, physical logic of model decision making is revealed by quantifying feature contributions and comparing action directions of the feature contributions in different samples; and through a dynamic fusion strategy based on interpretability analysis, the transparency of a prediction result and the visualization of a decision basis are realized. Finally, the technical effect of remarkably improving the credibility, operability and self-adaptive decision-making ability of a short-time heavy rainfall prediction model in a complex weather situation is achieved, and artificial intelligence prediction is promoted to be converted from pure pursuit of accuracy to a normal form of interpretability, trustworthiness and optimizability.

Description

Technical Field

[0001] This invention relates to the field of precipitation forecasting technology, and in particular to a precipitation forecasting method based on SHAP meteorological feature analysis. Background Technology

[0002] Rainfall is a crucial link in the water cycle, and timely and accurate rainfall forecasting is essential for understanding climate patterns, preventing floods, and supporting water resource management and sustainable development. Short-duration heavy rainfall (more than 20 mm per hour) is a typical small-to-medium-scale severe convective weather phenomenon characterized by its suddenness, small spatial scale, short duration, and rapid evolution. Despite its brief duration, it easily triggers secondary disasters such as urban flooding, flash floods, and the rapid rise of small and medium-sized rivers, posing a serious threat to people's lives, property, and infrastructure. Therefore, it has become a key focus and challenge in meteorological forecasting and disaster prevention and mitigation.

[0003] Currently, short-term heavy precipitation forecasting in operational practice mainly relies on two traditional methods: one is numerical weather prediction (NWP) models, which are based on solving physical equations and can provide large-scale circulation backgrounds. However, due to limitations in grid resolution and initial field errors, their ability to capture the fine details of small- and medium-scale convective systems is insufficient. The other is the judgment based on meteorological experts' experience, which involves subjective analysis based on numerical forecasts and multi-source observations. Its effectiveness is greatly influenced by the forecaster's experience, and its reproducibility and objectivity are relatively weak. Therefore, traditional methods still have significant limitations in forecasting the location, intensity, and timeliness of short-term heavy precipitation.

[0004] In recent years, with the improvement of the spatiotemporal resolution of atmospheric reanalysis data, a series of physical fields describing thermal, dynamic, and water vapor conditions (such as temperature, humidity, vertical velocity, convective available potential energy (CAPE), and water vapor flux) have been widely used in precipitation mechanism research and forecasting modeling. Studies have shown that these physical factors are closely related to the occurrence and development of short-duration heavy precipitation. Meanwhile, machine learning methods, due to their powerful nonlinear fitting and feature learning capabilities, have shown significant advantages in the field of precipitation prediction. For example, the existing technology, published under CN114896894A, discloses a machine learning-based precipitation prediction method and system that uses Support Vector Machines (SVM) and its optimized variants for precipitation prediction. By introducing atmospheric circulation factors (such as the Niño1+2 exponent) as input features, it constructs GA-SVM or GA-LSSVM models, improving the accuracy of monthly precipitation prediction. In addition, the existing technology disclosed in CN116306215A is a precipitation prediction method and system based on ensemble learning. It also studies the use of ensemble learning strategies to integrate multiple machine learning models such as random forest, XGBoost and extremely random tree, and introduces multi-source satellite and reanalysis data to further improve the accuracy of quantitative precipitation estimation through optimal interpolation fusion.

[0005] The aforementioned and existing related technologies often suffer from the following shortcomings: Although machine learning methods have made positive progress in precipitation forecasting, most current models still exhibit "black box" characteristics, meaning that the model prediction process lacks physical interpretability, making it difficult to clarify the specific contributions of each input feature to the prediction results and their interrelationships. For forecast objects like short-duration heavy precipitation, whose physical mechanisms are clearly defined, model interpretability not only affects the reliability of forecast results but also restricts the depth of feature selection, model optimization, and operational applications. Therefore, how to enhance model interpretability and clarify the role of key meteorological features while maintaining the predictive performance of machine learning has become a critical issue that urgently needs to be addressed in current research on intelligent forecasting of short-duration heavy precipitation. Summary of the Invention

[0006] The technical problem to be solved by this invention is that existing technologies have the disadvantage of not being able to enhance the interpretability of the model and clarify the role of key meteorological features while maintaining the predictive performance of machine learning. To this end, we propose a precipitation prediction method based on SHAP meteorological feature analysis.

[0007] To achieve the above objectives, this application adopts the following technical solution: a precipitation prediction method based on SHAP meteorological characteristic analysis, comprising the following steps:

[0008] S1. Based on meteorological reanalysis data and precipitation data from ground observation stations, a short-term heavy precipitation dataset containing various meteorological and physical characteristics is constructed. The dataset includes positive and negative samples.

[0009] S2. Based on the first machine learning model combined with the SHAP framework, calculate the Shapley value of each feature in the dataset, and perform fine-tuning based on the Shapley value to obtain a subset of key features; the fine-tuning includes calculating and comparing the SHAP contribution direction of each feature under positive and negative samples;

[0010] S3. Using the aforementioned key feature subset, train at least two different types of second machine learning models;

[0011] S4. Analyze the trained second machine learning model based on the SHAP framework to quantify the contribution of each feature in the key feature subset to the model's prediction decision.

[0012] S5. Combine the trained second machine learning model with the numerical weather prediction model, and based on the interpretability analysis results of step S4, use a dynamic decision algorithm to perform multi-model fusion to generate a short-term heavy precipitation prediction model for real-time forecasting.

[0013] Furthermore, S1 includes the following steps:

[0014] S11. Based on hourly precipitation observation data from ground automatic stations, determine positive samples for short-term heavy precipitation;

[0015] S12. Based on meteorological reanalysis grid data, calculate and extract various physical diagnostic quantities covering thermal, dynamic, water vapor and stratification conditions to form an initial feature pool.

[0016] S13. The physical diagnostic quantities in the initial feature pool are spatially interpolated to the location of the ground automatic station, and spatiotemporally matched with the precipitation observation labels to construct a preliminary dataset.

[0017] S14. Using a process-oriented sampling strategy, for each regional short-term heavy precipitation event, extract the corresponding positive samples and corresponding negative samples to construct a balanced dataset for model training and testing.

[0018] Furthermore, step S2 includes the following steps:

[0019] S21. Based on a dataset containing all candidate physical features, train an initial machine learning model and apply standard SHAP analysis to calculate the SHAP value of each feature, obtaining an initial global feature importance ranking based on the average absolute SHAP value.

[0020] S22. Introduce the time-delay embedding SHAP analysis algorithm to construct a time-delay feature vector for each physical feature, which includes the state of the current time and multiple consecutive historical time moments before it.

[0021] S23. Based on the SHAP contribution values ​​of each feature at different time delays obtained in step S22, form its time-series contribution sequence, and perform weighted integration on the time-series contribution sequence to calculate a comprehensive importance index that reflects the continuous influence of the feature during the key evolution period.

[0022] S24. Calculate the average SHAP value of each feature in the positive samples of short-term heavy precipitation and the negative samples of non-precipitation; combine the comprehensive importance index with the difference between the average SHAP values ​​of the positive and negative samples to make a judgment and select the core feature subset with the strongest discrimination ability.

[0023] Furthermore, the time-delay embedding SHAP analysis algorithm in step S22 is specifically as follows:

[0024] For each physical feature Construct time-delay eigenvectors ,in For time step, The time delay window length;

[0025] A custom weighted Gaussian kernel function is used when calculating the SHAP value. To handle the similarity between time delay vectors, where For time decay weight, The attenuation coefficient is... For kernel width parameter, Represents the sample Features The time-delay sequence, Represents the sample Features The time-delay sequence.

[0026] Furthermore, step S23 specifically includes:

[0027] For features Its time-series contribution sequence is Its comprehensive importance The calculation method is as follows: ,in To correspond to different time delays The integration weight.

[0028] Furthermore, step S3 specifically includes:

[0029] S31. Reconstruct the training set and test set using the key feature subset selected in step S2;

[0030] S32. Select three gradient boosting tree algorithms, namely CatBoost, XGBoost and LightGBM, as the second machine learning model respectively.

[0031] S33. For each second machine learning model, within its preset hyperparameter space, the grid search cross-validation method is used to optimize the parameters to obtain the optimized model.

[0032] Furthermore, step S4 specifically includes:

[0033] S41. Input the test set data into each of the trained second machine learning models, and use the SHAP framework to calculate the SHAP value of each feature in the key feature subset in each model.

[0034] S42. Based on the calculated SHAP value, perform a global interpretability analysis to obtain the overall contribution ranking and distribution of each feature in the model;

[0035] S43. For both positive and negative samples, calculate the average SHAP value of the key features, and use the interaction effect calculation function of SHAP to decompose the total contribution of the features into main effects and interaction effects, so as to analyze the independent and synergistic effects of the features.

[0036] Furthermore, step S5 includes the following steps:

[0037] S51. Construct a static model fusion benchmark based on historical data;

[0038] S52. Calculate the gridded forecast field for the future time period output by the numerical weather prediction model into the feature values ​​corresponding to the key feature subset, and input them into each second machine learning model to obtain their respective probability forecast results;

[0039] S53. Execute the dynamic decision algorithm to determine the final forecast output based on the degree of divergence between the forecast results of each model and the real-time feature confidence based on SHAP.

[0040] Furthermore, the dynamic decision-making algorithm in step S53 specifically includes:

[0041] S531. Calculate the divergence between the predicted probabilities of each second machine learning model. ,in Let be the predicted probability of the m-th model. The mean of the predicted probabilities of all models. Number of models;

[0042] S532. For each second machine learning model, calculate its confidence score for the current forecast time. ,in Score the historical performance of the model. The current feature matching degree based on the SHAP value. These are the weighting coefficients;

[0043] S533, Based on the aforementioned degree of divergence With preset threshold , Based on the comparison results, choose to implement one of the following decisions:

[0044] when At that time, the static fusion result or the model average result shall be adopted;

[0045] when When, adopt confidence scores The highest-performing model is used as the dominant model, with its forecast results and SHAP interpretation as the final output;

[0046] when At that time, adopt the following The weighted average of the predicted probabilities of each model is given.

[0047] Furthermore, the current feature matching degree based on the SHAP value in step S532 The calculation method is as follows:

[0048] ,

[0049] in, This represents the set of the top K most important features in the historical SHAP analysis of model m. The feature vector of the current forecast time The SHAP value of the j-th feature in the model. Let be the average SHAP value of feature j on historical positive samples.

[0050] The technical effects and advantages of this invention are as follows:

[0051] This invention systematically solves the "black box" problem prevalent in existing machine learning precipitation prediction models by constructing a dataset, selecting key features based on SHAP, training multiple models, conducting deep interpretable analysis, and integrating with business processes. This method not only utilizes the SHAP framework for feature selection to improve model efficiency, but more importantly, it deeply integrates interpretability analysis into the entire process of model construction and business decision-making. By quantifying feature contributions and comparing their roles in different samples, it reveals the physical logic of model decisions; and through a dynamic fusion strategy based on interpretability analysis, it achieves transparency in prediction results and visualization of decision-making basis. Ultimately, it significantly improves the credibility, operability, and adaptive decision-making capabilities of short-term heavy precipitation prediction models under complex weather conditions, promoting a paradigm shift in artificial intelligence forecasting from simply pursuing accuracy to "interpretable, trustworthy, and optimizable" models. Attached Figure Description

[0052] The disclosure of this invention is illustrated with reference to the accompanying drawings. It should be understood that the drawings are for illustrative purposes only and are not intended to limit the scope of protection of this invention. In the drawings, the same reference numerals are used to refer to the same parts:

[0053] Figure 1 This is an overall flowchart of the present invention;

[0054] Figure 2 The ranking of the importance of Shapley values ​​and the bee heatmap of 72 features in 5 categories in short-term heavy precipitation prediction according to the present invention;

[0055] Figure 3 This is a map showing the research area and site distribution of this invention. Detailed Implementation

[0056] It is readily understood that, based on the technical solution of this invention, those skilled in the art can propose various interchangeable structural methods and implementations without altering the essential spirit of the invention. Therefore, the following detailed embodiments and accompanying drawings are merely illustrative examples of the technical solution of this invention and should not be considered as the entirety of the invention or as limitations or restrictions on the technical solution of this invention.

[0057] This embodiment uses Anhui Province, China as an example, referring to... Figure 1 As shown, the present invention provides a technical solution: a precipitation prediction method based on SHAP meteorological feature analysis, comprising the following steps:

[0058] Step S1: Construct a multi-source spatiotemporal matching dataset of short-term heavy precipitation, as follows:

[0059] Target data (label) acquisition and quality control: Data was collected from 974 automatic weather stations under the Anhui Provincial Meteorological Bureau (including 16 municipal stations, 64 county stations, and 894 township stations; station distribution is shown in the appendix). Figure 3 Hourly precipitation observation data during the flood season (usually May to September) from 2021 to 2024. Short-duration heavy precipitation events are defined as: hourly precipitation (R) ≥ 20.0 mm at a single station. Samples meeting this condition are labeled as positive samples (label y=1), otherwise as negative samples (label y=0).

[0060] To ensure sample quality, quality control was performed on precipitation observation data: hourly automatic weather station precipitation data were spatiotemporally matched with the corresponding weather radar combined reflectivity factor product. For each station reporting R≥20mm, radar reflectivity values ​​were extracted for 10 time intervals (6-minute intervals) within one hour before the start of precipitation. Only when the average of these 10 reflectivity values ​​was greater than 20dBZ was the sample finally confirmed as a valid positive sample of short-duration heavy precipitation, thus eliminating some false alarms.

[0061] Feature data (forecasting factors) calculation: The ERA5 reanalysis dataset (spatial-temporal resolution: 0.25° × 0.25°, hourly) was obtained from the European Centre for Medium-Range Weather Forecasts (ECMWF) and was contemporaneous with precipitation observations in the same region. Based on the three-dimensional gridded data (temperature, air pressure, humidity, wind field, geopotential height, etc.) provided by ERA5, and according to meteorological principles, a total of 72 physical diagnostic quantities covering 5 major categories were calculated as the initial feature pool for the model input. A detailed list of these 72 physical quantities across 5 major categories is provided in Table 1.

[0062] Table 1 shows 72 physical characteristic quantities in 5 major categories based on ERA5 data input.

[0063]

[0064]

[0065] Table 1 explains: This table systematically lists all the meteorological and physical features used in constructing the prediction model. These features comprehensively cover the thermodynamic, dynamic, water vapor, and stratification conditions closely related to the development of severe convection, providing a complete multi-dimensional input information foundation for the model to learn the complex formation mechanism of short-duration heavy precipitation.

[0066] Spatiotemporal matching of data and dataset construction:

[0067] The calculated two-dimensional grid data (0.25°×0.25°) of 72 physical quantities were interpolated to the geographical locations of 974 observation stations using bilinear interpolation, forming a 72-dimensional feature vector for each station on an hourly basis.

[0068] The processed feature data and the quality-controlled precipitation label data are precisely matched by station and time to form the initial dataset.

[0069] To address the issue of imbalanced samples (positive samples being far fewer than negative samples), a process-based sampling strategy is employed to construct a balanced training set: A regional short-duration heavy rainfall event is defined as the occurrence of short-duration heavy rainfall at any two or more stations in any of the 80 districts and counties of Anhui Province. For each event, all stations experiencing heavy rainfall are designated as positive samples. Simultaneously, a sample size 1.2 times the number of positive samples is randomly selected from stations that did not experience heavy rainfall during the event (including the 80 district / county stations and city-level stations) and from stations randomly selected outside the event period. This strategy ensures the spatiotemporal representativeness of the negative samples and avoids the bias that might be introduced by random undersampling.

[0070] All samples are sorted chronologically, and the data from the first 80% of the time period is selected as the training set, while the last 20% is used as an independent test set to evaluate the generalization performance of the model.

[0071] Step S2: Select key physical features based on the CatBoost and SHAP framework, as follows:

[0072] This step utilizes machine learning models and interpretability analysis to select the core factor set that contributes most to the prediction of short-term heavy rainfall from high-dimensional features.

[0073] Training the initial CatBoost model: Using the training set (containing all 72 features) constructed in step S1, train a CatBoost classifier with default hyperparameters and evaluate its basic performance (such as AUC score) using five-fold cross-validation to ensure that the model has initial learning capabilities.

[0074] Standard SHAP Global Contribution Analysis: Using the TreeExplainer interpreter within the SHAP framework, the mean absolute SHAP value of each feature across all samples on the training set is calculated. This value reflects the overall magnitude of the feature's influence on the model output. Based on this value, all 72 features are sorted in descending order to obtain an initial ranking of global feature importance.

[0075] To capture the evolution signals of meteorological elements before the occurrence of severe convection, a time-delay embedding SHAP analysis algorithm is introduced, specifically including:

[0076] Time-delay feature construction: For each physical feature Construct a time-delay eigenvector In this embodiment, a time step is set. hours, time delay window length That is, to examine the characteristic state of the current moment and the previous 5 hours.

[0077] Improved kernel function calculation: During the SHAP value calculation, a custom weighted Gaussian kernel function is used to handle the time delay vector when evaluating the similarity between two samples.

[0078]

[0079] in, It is a time decay weight. It is the attenuation coefficient (set to 0.5 in this embodiment), used to control the attenuation rate of the influence of historical time features on the current similarity assessment; This is the kernel width parameter, used to adjust the sensitivity of similarity to differences in eigenvalues. This function assigns greater weight to recent eigenvalues, consistent with the fact that nowcasting relies more on prior meteorological knowledge of recent conditions.

[0080] Time-varying contribution and overall importance: The algorithm described above is used to obtain the contribution of each feature. Different time delays SHAP contribution sequence on .feature The overall importance Defined as the weighted sum of the absolute values ​​of each time delay contribution: ,in It uses integrated weights (in this embodiment, it is set to uniform weights). This algorithm can uniquely identify those precursory signal features that begin to contribute significantly several hours before precipitation occurs, such as the accumulation of water vapor or slow changes in stability.

[0081] Refined screening based on SHAP comparison of positive and negative samples: Calculate the average SHAP value of each feature in the positive and negative sample subsets. An ideal discriminative feature should exhibit a large positive contribution in positive samples and a large negative contribution (or close to zero) in negative samples. Combined with global importance ( By comparing the contributions of positive and negative samples, the top 20 features with the most significant contributions and the clearest physical meaning were selected, forming a subset of key features. The feature importance ranking and examples of positive and negative sample contributions obtained through this step are shown in Table 2.

[0082] Table 2 shows the top 20 most important features obtained using the CatBoost-Shapley method.

[0083]

[0084]

[0085]

[0086] Table 2 illustrates the results of feature optimization. SHAP_ABS_MEAN measures the overall magnitude of the feature's influence; SHAP_MEAN indicates the directional contribution of increasing feature values ​​to the prediction results (positive values ​​induce precipitation, negative values ​​inhibit precipitation). By comparing the SHAP_MEAN values ​​of positive and negative samples, the bidirectional effect of features on the forecast results can be clearly seen. For example, the 500 hPa vertical velocity contributes positively when precipitation is predicted and negatively when no precipitation is predicted, demonstrating clear physical meaning and making it a typical strong discriminant factor.

[0087] Step S3: Train and optimize multiple machine learning models based on key feature subsets, as detailed below:

[0088] Data reconstruction: Using the 20 key features selected in step S2, reconstruct the feature portions of the training and test sets.

[0089] Multi-model training and hyperparameter optimization: Classification models were constructed using three gradient boosting tree algorithms: CatBoost, XGBoost, and LightGBM. For each algorithm, the GridSearchCV method was used to optimize the system within its key hyperparameter space (see Table 3). The AUC (area under the curve) was chosen as the evaluation metric to find the parameter combination that performs optimally on the training set and has strong generalization ability.

[0090] Table 3 Hyperparameter ranges for three types of machine learning models

[0091]

[0092] Table 3 explains: This table lists the core hyperparameters used to construct the three base models and their tuning ranges. By performing systematic grid search cross-validation on these parameters, the optimal configuration for the short-duration heavy precipitation prediction task can be found for each model, ensuring model performance and preventing overfitting.

[0093] Model performance benchmark evaluation: The three optimized models were evaluated on independent test sets, and their precision, recall, and F1 score were recorded as benchmarks for model performance and references for subsequent fusion. These metrics were calculated based on the confusion matrix (see Table 4).

[0094] Table 4 is a schematic table of confusion matrix.

[0095]

[0096] Table 4 explains that the confusion matrix is ​​a fundamental tool for evaluating the performance of classification models. This table shows four possible combinations of predicted results and the true situation in a binary classification problem: TP, FN, FP, and TN. Accuracy can be calculated based on this matrix. Recall rate and F1 score Key indicators such as...

[0097] Step S4: In-depth interpretability analysis of the model decision-making process, as detailed below:

[0098] Global interpretability analysis: For the three optimized models, their SHAP values ​​are calculated on the test set. A global feature importance bar chart is plotted to show the ranking of the average absolute SHAP values ​​of each feature. A SHAP beehive plot is also plotted to visually demonstrate the relationship between the magnitude (color) of each feature value and the SHAP value (horizontal position), and how the overall distribution of features affects the predicted output.

[0099] Sample category specificity analysis and contribution decomposition: Calculate the average SHAP value of key features in positive and negative samples of the test set respectively, and compare them, as described in step S2, to deepen the understanding of the direction of feature effect.

[0100] Furthermore, this invention employs the interaction effect calculation function of SHAP to decompose the total contribution of a feature to the model prediction into main effects (the contribution generated independently by the feature) and interaction effects (the contribution generated by the synergistic effect of the feature with other features). By performing contribution decomposition on positive and negative samples respectively, the analysis results shown in Table 5 can be obtained.

[0101] Table 5. Decomposition of the total effect, main effect, and interaction effect of the top 5 features in positive and negative samples.

[0102]

[0103] Table 5 illustrates the decomposition of feature contributions. For example, for the A index in the positive sample, its total contribution (0.327) can be decomposed into a large main effect (0.288) and a small interaction effect (0.038), indicating that its predictive power mainly stems from its independent effect. However, for the A index in the negative sample, its interaction effect (0.093) is significant, suggesting that when predicting no rain, its effect relies more on the synergistic inhibition relationship with other features. This decomposition provides forecasters with more detailed decision-making insights, revealing the differences in dominant physical processes under different weather conditions.

[0104] Step S5: Generate a business prediction model that incorporates interpretable information, as follows:

[0105] Constructing a historical static fusion benchmark: The predicted probability outputs of the three base learners on the historical test set are used as new features to train a random forest model as a static fusion processor. This model provides a performance fusion benchmark on historical data.

[0106] Deploying a real-time business forecasting process:

[0107] Feature forecast acquisition: In the operational environment, access the gridded forecast field for the next 72 hours with a time resolution of 3 hours from ECMWF (or other operational numerical models). Based on the definition formulas of the 20 key features determined in step S2, calculate the feature values ​​for each forecast time in real time.

[0108] Base model prediction: The calculated real-time feature values ​​are input into three pre-trained base learners to obtain three independent short-term heavy precipitation probability forecasts: .

[0109] Perform multi-model dynamic decision-making:

[0110] a. Calculate the forecast divergence: Calculate the divergence between the forecast probabilities of the three base models. Here, D is the forecast divergence index among models, which is used to quantify the degree of consistency between the prediction results of the three base models; These represent the predicted probabilities output by the CatBoost, XGBoost, and LightGBM models, respectively. The arithmetic mean of these three probabilities. The larger the value, the higher the uncertainty between models.

[0111] b. Calculate the real-time feature confidence based on SHAP: Using the SHAP interpreter trained for each base model in step S4, calculate the feature vector for the current forecast time. The SHAP value under this model. Define the model. Confidence score for the current situation The calculation is as follows:

[0112] .

[0113] in, It is an adjustment weight used to balance historical performance weights and current context confidence weights; It is a model The AUC score (normalized) on the historical test set represents its historical performance; This is the current feature matching degree.

[0114] This is the current feature matching degree, calculated using the following formula:

[0115] .

[0116] here, It is a model Ranked among the most important in historical analysis The feature set, Current features The SHAP value of this model, It is a feature The average SHAP value on historical positive samples (obtainable from the analysis results in step S4). The smaller the value, the more similar the current combination of environmental features is to the heavy precipitation pattern in the model's memory, and the higher the model's confidence in the current forecast.

[0117] c. Dynamic decision-making rules:

[0118] Low divergence scenario ( If there is a high degree of agreement among the models, then the result of the static fusion unit or the average probability of the three models is directly adopted. , along with a consistent SHAP explanation.

[0119] High divergence scenarios If the discrepancies between the models are significant, then a confidence-based decision is initiated. A confidence score is selected. The model with the highest prediction probability is used as the dominant model. As the final forecast result, the SHAP interpretation plot of the dominant model is highlighted.

[0120] Medium divergence scenario ( The probabilities of three models can be used. The weighted average of the weights is used as the final output.

[0121] Operational Product Output: The final forecast product not only includes gridded occurrence probabilities and landing area forecasts, but also a mandatory SHAP interpretation diagram of the dominant or integrated model selected by the dynamic decision-making algorithm. This diagram clearly identifies the top contributing physical features, their contribution values, and directions (positive / negative), providing conclusions in natural language or visualization, such as the main basis for this forecast being unusually abundant lower-level water vapor and strong upward motion in the middle level. This significantly enhances the decision-making support value of AI products in forecasters' work.

[0122] The technical scope of this invention is not limited to the content described above. Those skilled in the art can make various modifications and variations to the above embodiments without departing from the technical concept of this invention, and all such modifications and variations should fall within the protection scope of this invention.

Claims

1. A precipitation prediction method based on SHAP meteorological feature analysis, characterized in that, Includes the following steps: S1. Based on meteorological reanalysis data and precipitation data from ground observation stations, a short-term heavy precipitation dataset containing various meteorological and physical characteristics is constructed. The dataset includes positive and negative samples. S2. Based on the first machine learning model combined with the SHAP framework, calculate the Shapley value of each feature in the dataset, and perform fine-tuning based on the Shapley value to obtain a subset of key features; the fine-tuning includes calculating and comparing the SHAP contribution direction of each feature under positive and negative samples; S3. Using the aforementioned key feature subset, train at least two different types of second machine learning models; S4. Analyze the trained second machine learning model based on the SHAP framework to quantify the contribution of each feature in the key feature subset to the model's prediction decision. S5. Combine the trained second machine learning model with the numerical weather prediction model, and based on the interpretability analysis results of step S4, use a dynamic decision algorithm to perform multi-model fusion to generate a short-term heavy precipitation prediction model for real-time forecasting.

2. The precipitation prediction method based on SHAP meteorological feature analysis according to claim 1, characterized in that: S1 includes the following steps: S11. Based on hourly precipitation observation data from ground automatic stations, determine positive samples for short-term heavy precipitation; S12. Based on meteorological reanalysis grid data, calculate and extract various physical diagnostic quantities covering thermal, dynamic, water vapor and stratification conditions to form an initial feature pool. S13. The physical diagnostic quantities in the initial feature pool are spatially interpolated to the location of the ground automatic station, and spatiotemporally matched with the precipitation observation labels to construct a preliminary dataset. S14. Using a process-oriented sampling strategy, for each regional short-term heavy precipitation event, extract the corresponding positive samples and corresponding negative samples to construct a balanced dataset for model training and testing.

3. The precipitation prediction method based on SHAP meteorological feature analysis according to claim 1, characterized in that: Step S2 includes the following steps: S21. Based on a dataset containing all candidate physical features, train an initial machine learning model and apply standard SHAP analysis to calculate the SHAP value of each feature, obtaining an initial global feature importance ranking based on the average absolute SHAP value. S22. Introduce the time-delay embedding SHAP analysis algorithm to construct a time-delay feature vector for each physical feature, which includes the state of the current time and multiple consecutive historical time moments before it. S23. Based on the SHAP contribution values ​​of each feature at different time delays obtained in step S22, form its time-series contribution sequence, and perform weighted integration on the time-series contribution sequence to calculate a comprehensive importance index that reflects the continuous influence of the feature during the key evolution period. S24. Calculate the average SHAP value of each feature in the positive samples of short-term heavy precipitation and the negative samples of non-precipitation; combine the comprehensive importance index with the difference between the average SHAP values ​​of the positive and negative samples to make a judgment and select the core feature subset with the strongest discrimination ability.

4. The precipitation prediction method based on SHAP meteorological feature analysis according to claim 3, characterized in that: The time-delay embedding SHAP analysis algorithm in step S22 is specifically as follows: For each physical feature Construct time-delay eigenvectors ,in For time step, The time delay window length; A custom weighted Gaussian kernel function is used when calculating the SHAP value. To handle the similarity between time delay vectors, where For time decay weight, The attenuation coefficient is... For kernel width parameter, Represents the sample Features The time-delay sequence, Represents the sample Features The time-delay sequence.

5. A precipitation prediction method based on SHAP meteorological feature analysis according to claim 3 or 4, characterized in that: Step S23 specifically involves: For features Its time-series contribution sequence is Its comprehensive importance The calculation method is as follows: ,in To correspond to different time delays The integration weight.

6. The precipitation prediction method based on SHAP meteorological feature analysis according to claim 1, characterized in that: Step S3 specifically includes: S31. Reconstruct the training set and test set using the key feature subset selected in step S2; S32. Select three gradient boosting tree algorithms, namely CatBoost, XGBoost and LightGBM, as the second machine learning model respectively. S33. For each second machine learning model, within its preset hyperparameter space, the grid search cross-validation method is used to optimize the parameters to obtain the optimized model.

7. The precipitation prediction method based on SHAP meteorological feature analysis according to claim 1, characterized in that: Step S4 specifically includes: S41. Input the test set data into each of the trained second machine learning models, and use the SHAP framework to calculate the SHAP value of each feature in the key feature subset in each model. S42. Based on the calculated SHAP value, perform a global interpretability analysis to obtain the overall contribution ranking and distribution of each feature in the model; S43. For both positive and negative samples, calculate the average SHAP value of the key features, and use the interaction effect calculation function of SHAP to decompose the total contribution of the features into main effects and interaction effects, so as to analyze the independent and synergistic effects of the features.

8. The precipitation prediction method based on SHAP meteorological feature analysis according to claim 1, characterized in that: Step S5 includes the following steps: S51. Construct a static model fusion benchmark based on historical data; S52. Calculate the gridded forecast field for the future time period output by the numerical weather prediction model into the feature values ​​corresponding to the key feature subset, and input them into each second machine learning model to obtain their respective probability forecast results; S53. Execute the dynamic decision algorithm to determine the final forecast output based on the degree of divergence between the forecast results of each model and the real-time feature confidence based on SHAP.

9. A precipitation prediction method based on SHAP meteorological feature analysis according to claim 8, characterized in that: The dynamic decision-making algorithm in step S53 specifically includes: S531. Calculate the divergence between the predicted probabilities of each second machine learning model. ,in Let be the predicted probability of the m-th model. The mean of the predicted probabilities of all models. Number of models; S532. For each second machine learning model, calculate its confidence score for the current forecast time. ,in Score the historical performance of the model. The current feature matching degree based on the SHAP value. These are the weighting coefficients; S533, Based on the aforementioned degree of divergence With preset threshold , Based on the comparison results, choose to implement one of the following decisions: when At that time, the static fusion result or the model average result shall be adopted; when When, adopt confidence scores The highest-performing model is used as the dominant model, with its forecast results and SHAP interpretation as the final output; when At that time, adopt the following The weighted average of the predicted probabilities of each model is given.

10. A precipitation prediction method based on SHAP meteorological feature analysis according to claim 9, characterized in that: The current feature matching degree based on SHAP value in step S532 The calculation method is as follows: ， in, This represents the set of the top K most important features in the historical SHAP analysis of model m. The feature vector of the current forecast time The SHAP value of the j-th feature in the model. Let be the average SHAP value of feature j on historical positive samples.

Images