Ensemble forecasting method for target events considering initial disturbance and model uncertainty

Abstract

The present invention relates to a target event ensemble prediction method that considers initial perturbations and model uncertainties, comprising the following steps: S1, constructing an orthogonal joint conditional nonlinear optimal perturbation model based on initial perturbation variables and model parameter perturbation variables corresponding to the target event; S2, generating initial perturbation samples and obtaining an eigenvector matrix after preprocessing; S3, solving the orthogonal joint conditional nonlinear optimal perturbation model to obtain multiple combinations of initial perturbations and model parameter perturbations; S4, superimposing the model parameter perturbations obtained in step S3 and the dimensionally upgraded initial perturbations on the model parameters and initial fields of an earth system model, and obtaining multiple sets of perturbation prediction results through numerical integration; S5, averaging the control prediction results and the multiple sets of perturbation prediction results to obtain an ensemble prediction result for the target event. Compared with the prior art, the present invention can further improve the accuracy and reliability of forecasting target events.

Description

Technical Field

[0001] The present invention belongs to the fields of computer science and meteorological science and technology, and in particular relates to a target event ensemble forecasting method taking into account initial disturbances and pattern uncertainties. Background Art

[0002] Earth System Models (ESMs) are the primary tools for numerical simulation and prediction of weather and climate. However, the current predictive capabilities of ESMs are primarily affected by two types of uncertainty: the first is initial uncertainty, which stems from inaccuracies in observational data and assimilation systems, resulting in initial conditions that only approximately represent the true atmospheric state and inevitably contain initial errors; the second is model uncertainty, which is related to the inherent limitations of ESMs, including but not limited to subgrid-scale parameterization, artificially set empirical parameters, and truncation errors of numerical methods. Due to the chaotic nature of the atmospheric system, these uncertainties can be rapidly amplified, leading to a complete loss of forecast skill some time in the future. Therefore, forecasting methods that rely on a single determinism become unreliable.

[0003] To address these uncertainties, researchers have begun shifting from single deterministic forecasting methods to probabilistic forecasting methods, namely ensemble forecasting techniques. Ensemble forecasting estimates the probability density distribution of future atmospheric states by constructing a set of perturbation members. By averaging the prediction results of multiple ensemble members, unpredictable factors can be filtered out while retaining the common components between members. In addition, the ensemble discreteness provided by this method can describe the uncertainty of the actual atmospheric state, thereby enhancing the reliability of the forecast results. A successful ensemble forecast aims to achieve an accurate estimate of the atmospheric probability distribution by generating perturbations that can accurately represent the forecast uncertainty. For example, Chinese invention patent CN116520457B discloses a regional ensemble forecast method for major meteorological activities under complex terrain. This method considers the initial value perturbation under the uncertain influence of complex terrain and geomorphological characteristics in different regions on the numerical weather forecast model, and combines the random perturbation technology of physical processes related to time and space scales to create a corresponding rapid, high-resolution regional ensemble forecast system for extreme weather and major meteorological activities, effectively improving the work efficiency and forecast accuracy of the regional ensemble forecast process. However, existing ensemble forecasting methods often only consider initial uncertainty or model uncertainty when constructing perturbations, or simply superimpose these two types of uncertainty. These methods ignore the potentially complex nonlinear relationships between initial conditions and model parameters. This fails to fully capture forecast uncertainty, particularly when dealing with extreme weather events or complex climate patterns. Therefore, it is necessary to design an ensemble forecasting method to further improve the accuracy and reliability of forecasts of target events. Summary of the Invention

[0004] The purpose of the present invention is to overcome the defects of the above-mentioned prior art and provide a target event ensemble forecasting method taking into account initial disturbances and pattern uncertainties, so as to further improve the accuracy and reliability of forecasting target events.

[0005] The purpose of the present invention can be achieved by the following technical solutions:

[0006] The present invention provides a target event ensemble prediction method considering initial disturbance and model uncertainty, which is used to predict target events in any mode of an earth system model, comprising the following steps:

[0007] S1. Based on the initial disturbance variables and pattern parameter disturbance variables corresponding to the target event, as well as their corresponding range constraints, an orthogonal joint conditional nonlinear optimal disturbance model is constructed;

[0008] S2, collecting simulated field data of the Earth system model, generating initial disturbance samples, and obtaining eigenvector matrices after preprocessing;

[0009] S3, solving the orthogonal joint conditional nonlinear optimal perturbation model in a low-dimensional orthogonal subspace to obtain a plurality of combinations of initial perturbations and pattern parameter perturbations, and mapping the initial perturbations in each combination back to the high-dimensional space through the eigenvector matrix obtained in step S2;

[0010] S4, superimposing the model parameter perturbation obtained in step S3 and the initial perturbation after dimensionality upgrade on the model parameters and initial field of the Earth system model respectively, and obtaining multiple sets of perturbation prediction results through numerical integration;

[0011] S5. Taking an average of the control forecast result and the multiple groups of disturbance forecast results obtained in step S4 to obtain an ensemble forecast result of the target event, wherein the control forecast result is a result of a natural integration of the earth system model without adding disturbance.

[0012] Furthermore, in step S1, the objective function of the orthogonal joint conditional nonlinear optimal perturbation model is specifically expressed as follows:

[0013]

[0014] Among them, J is the objective function of the model, is the optimal joint perturbation, s 0i is the initial perturbation, p i is the mode parameter perturbation, S0 and P represent the original initial field and the default mode parameter value respectively, M t (P)(S0) represents the state of the nonlinear system, i.e., the Earth system model M, integrated to time t when no disturbance is added. t (P+p i)(S0+s 0i ) indicates that when adding the joint perturbation (s 0i ,p i ), the nonlinear system, that is, the earth system model M, is integrated to the state at time t, Ω i is a constraint condition.

[0015] Furthermore, the constraint condition Ω i The specific expression is as follows:

[0016]

[0017] Among them, Ω i is the constraint condition for the i-th optimization, δ s is the constraint range of the initial perturbation, δ p is the constraint range of the mode parameter perturbation.

[0018] Furthermore, in step S2, the preprocessing includes a normalization operation and a dimensionality reduction operation.

[0019] Furthermore, the Z-Score method is used for normalization.

[0020] Furthermore, principal component analysis is used to perform dimensionality reduction.

[0021] Furthermore, in step S3, the orthogonal joint conditional nonlinear optimal perturbation model is solved in a low-dimensional orthogonal subspace by a particle swarm optimization algorithm, a simulated annealing algorithm, a genetic algorithm or a differential evolution algorithm.

[0022] Furthermore, in the process of solving the orthogonal joint conditional nonlinear optimal perturbation model using the particle swarm optimization algorithm, the particle update method is as follows:

[0023]

[0024] in, and are the position and velocity of the jth particle at the kth iteration, respectively. The final particle position is the disturbance, ω k is the coefficient of speed at the kth iteration, which is initially set to 1 and gradually decreases in each iteration. c1 and c2 are cognitive and social learning factors, respectively. r1 and r2 are random numbers drawn from a uniform distribution in the interval [0,1]. X pb is the optimal position of the particle so far, X gb is the global optimal position of the entire particle swarm at the kth iteration.

[0025] Furthermore, the coefficient of speed at the kth iteration is k The specific expression is as follows:

[0026]

[0027] Among them, ω max and ω min are the maximum inertia weight and the minimum inertia weight, respectively, iter max is the maximum number of iterations, and k is the current number of iterations.

[0028] Furthermore, in step S4, the specific process of obtaining multiple sets of disturbance prediction results through numerical integration is as follows:

[0029] S401, selecting the first n combinations of initial perturbations and model parameter perturbations obtained by solving step S3, superimposing the model parameter perturbation and the initial perturbation after dimension increase in each combination on the model parameters and initial field of the Earth system model, and obtaining n perturbation forecasts by numerical integration;

[0030] S402. Superimpose the disturbances that are opposite to the model parameter disturbance and the initial disturbance after dimensional upgrading in each combination on the model parameters and initial field of the Earth system model respectively, obtain n disturbance forecasts through numerical integration, and finally obtain 2n disturbance forecast results.

[0031] Compared with the prior art, the present invention has the following beneficial effects:

[0032] 1. The present invention constructs an orthogonal joint conditional nonlinear optimal perturbation model based on the initial perturbation variables and model parameter perturbation variables corresponding to the target event, as well as their corresponding range constraints. It can focus on the dynamic coordinated growth of the initial perturbation and the model parameter perturbation, solve the perturbation combination with the maximum development in different orthogonal subspaces, and optimize the initial perturbation and the model parameter perturbation at the same time. It can more accurately and comprehensively capture the uncertainty in the ensemble forecast, thereby improving the accuracy and reliability of the forecast of the target event.

[0033] 2. The present invention adopts the principal component analysis method to perform dimensionality reduction operation on the initial disturbance sample to obtain the eigenvector matrix, and then solves the orthogonal joint conditional nonlinear optimal disturbance model in the low-dimensional space to obtain multiple combinations of initial disturbances and pattern parameter disturbances, and maps the initial disturbances in each combination back to the high-dimensional space through the eigenvector matrix. The complex orthogonal joint conditional nonlinear optimal disturbance optimization problem can be converted into a low-dimensional optimization problem, which significantly speeds up the solution of the optimization problem.

[0034] 3. The present invention solves the orthogonal joint conditional nonlinear optimal perturbation model through the particle swarm optimization algorithm, which can effectively address the challenge of directly solving the optimization problem in the absence of gradient information, and cleverly circumvents the computational obstacles that traditional optimization algorithms may encounter in complex earth system models without adjoint modes.

[0035] 4. The present invention superimposes the model parameter disturbance and the initial disturbance after dimensionality upgrade in each combination on the model parameters and the initial field of the earth system model model, and obtains n disturbance forecasts through numerical integration. At the same time, the disturbance that is the opposite of the model parameter disturbance and the initial disturbance after dimensionality upgrade in each combination is superimposed on the model parameters and the initial field of the earth system model model, and obtains n disturbance forecasts through numerical integration. Finally, 2n disturbance forecast results are obtained, and then the average value is taken with the control forecast result to obtain the ensemble forecast result, which can ensure that all disturbance biases are 0, further improving the accuracy and reliability of the ensemble forecast results of the target event. BRIEF DESCRIPTION OF THE DRAWINGS

[0036] Figure 1 is a flow chart of the method of the present invention;

[0037] Figure 2 This is a schematic diagram of the geographical location of the initial disturbance area in the embodiment;

[0038] Figure 3 Schematic diagram of the results of the comparative experiment in the embodiment, wherein the black legend Control is a control experiment, which is a single deterministic forecast without adding disturbances, the dark blue legend O-CNOP-Is-EF is an ensemble forecast only for the initial uncertainty, the green legend O-CNOP-Ps-EF is an ensemble forecast only for the model parameter uncertainty, the light blue legend O-CNOP-(I / P)s-EF is an ensemble forecast that considers two uncertainties through different members, the orange legend O-CNOP-(I+P)s-EF is an ensemble forecast that simultaneously considers two uncertainties in a simple superposition manner, and the red legend OC-CNOPs-EF is the method proposed in this application.

[0039] (3a) is the RMSE evolution diagram of the comparative experiment in the embodiment, the vertical axis RMSE represents the root mean square error, and the horizontal axis Day represents time.

[0040] (3b) is a comparison chart showing the improvement of the ensemble forecast compared to the control forecast in the comparative experiment in the embodiment.

[0041] (3c) is the evolution diagram of the ratio of RMSE to Spread (set dispersion) of the comparative experiment in the embodiment, Rartio represents the ratio,

[0042] (3d) is the Talagrand plot of the comparative experiment in the embodiment. DETAILED DESCRIPTION

[0043] The present invention is described in detail below with reference to the accompanying drawings and specific embodiments. This embodiment is implemented based on the technical solution of the present invention, and provides a detailed implementation method and specific operation process, but the protection scope of the present invention is not limited to the following embodiments.

[0044] Example:

[0045] This embodiment provides a target event ensemble prediction method that considers initial perturbations and model uncertainties, also known as the orthogonal combined conditional nonlinear optimal perturbations (OC-CNOPs), for predicting target events in any mode of Earth System Models (ESMs). The method includes the following steps:

[0046] S1. Based on the initial disturbance variables and pattern parameter disturbance variables corresponding to the target event, as well as their corresponding range constraints, an orthogonal joint conditional nonlinear optimal disturbance model is constructed.

[0047] Target events are those that have a significant impact on weather or climate forecasts, such as the El Niño-Southern Oscillation (ELNO) In this example, the NAO event is selected as the target event for ensemble forecasting, and the Community Earth System Model (CESM) is used to simulate and forecast NAO events. This model includes multiple modules such as the atmosphere, land, ocean, sea ice, and land ice, and is a fully coupled complex model.

[0048] Six initial disturbance variables with significant impact on NAO were selected, including zonal wind, meridional wind, temperature, specific humidity, surface pressure and surface geopotential height. The disturbance region was defined as the North Atlantic region, which is between 90°W and 40°E and 20°N and 80°N. Figure 2 To ensure the rationality of atmospheric circulation simulation, the constraint range of the initial perturbation variable δ s It is set to 10% of the ground state atmospheric dry energy, as shown in the following formula:

[0049]

[0050] Where D is the disturbance region, i.e. the North Atlantic region, U0, V0, T0 and SP0 represent the initial conditions of zonal wind, meridional wind, temperature and surface pressure respectively, and C p is the specific heat at constant pressure, which is 1005.7 J·kg -1 K-1 , T r is the reference temperature, which is 270K, R a is the ideal gas constant, which is 287.05 J·kg -1 K -1 ,π r is the reference static pressure, which is 1000hPa.

[0051] Eleven sensitive variables closely related to NAO events were selected as model parameter disturbance variables. The detailed description, default value, and range of each parameter are shown in Table 1.

[0052] Table 1 Detailed description of parameter disturbance variables, their default values ​​and parameter ranges

[0053]

[0054] The objective function of the OC-CNOPs method is constructed to ensure that the key variables affecting the NAO forecast can be optimized. In the actual forecast process, the initial error and the model parameter error affect each other. Therefore, when constructing the objective function, one cannot optimize one type of perturbation alone. Instead, one should focus on the dynamic coordinated growth of the initial perturbation and the model parameter perturbation and optimize both types of perturbations simultaneously. Assuming that the combined perturbation of the initial and model parameters is (s 0i , p i ), the OC-CNOPs optimization problem is to solve a set of perturbation combinations with maximum development in different orthogonal subspaces. To solve OC-CNOPs, first calculate the global (i.e., first) C-CNOP through optimization methods; then, calculate the second C-CNOP in a subspace orthogonal to the first C-CNOP; then, calculate the third C-CNOP in a subspace orthogonal to the first two C-CNOPs; and so on, calculate the i-th C-CNOP in a subspace orthogonal to the first i-1 C-CNOPs. Specifically, the i-th C-CNOP is a combination of perturbations that satisfies the following optimization problem:

[0055]

[0056]

[0057] Among them, J is the objective function, s 0i is the initial perturbation, p i is the mode parameter perturbation, S0 and P represent the original initial field and the default mode parameter value respectively, M t (P)(S0) represents the state of the nonlinear system, i.e., the Earth system model M, integrated to time t when no disturbance is added. t (P+p i )(S0+s 0i) indicates that when adding the joint perturbation (s 0i ,p i ), the nonlinear system, that is, the earth system model M, is integrated to the state at time t, Ω i is the constraint condition, δ s is the constraint range of the initial perturbation, δ p is the constraint range of the mode parameter perturbation. OC-CNOPs solves a set of optimal (i.e., with the largest development) initial and parameter joint perturbations in different orthogonal subspaces, namely Therefore, there is a greater chance of capturing information about extreme disturbances, which is very beneficial for forecasting extreme weather or climate events.

[0058] This embodiment constructs the objective function of the OC-CNOPs method based on the key characteristics of NAO events to ensure that it can optimize the key factors affecting NAO forecasts. The objective function of OC-CNOPs is designed to quantify the difference in the North Atlantic Oscillation Index (NAO index, NAOI) between the disturbance state and the reference state. Considering that NAO events are divided into positive and negative phases, this embodiment develops two different objective functions for OC-CNOPs that develop into NAO+ and NAO- events, respectively:

[0059]

[0060] Among them, z is the initial and parameter combined perturbation, NAOI refer Indicates the NAOI value of the reference state without any disturbance. + ) and NAOI(z - ) respectively represent the superposition of z + and z - The NAOI value of the perturbation state after the perturbation is maximized by The optimal combination of perturbations that leads to a NAO+ event can be determined. Similarly, by minimizing the objective function The optimal combination of perturbations that develop into NAO-events can be identified.

[0061] S2. Collect the simulated field data of the Earth system model, generate the initial disturbance samples, and obtain the eigenvector matrix after preprocessing.

[0062] Since the dimensions of the initial variables in ESMs are usually large, it is necessary to reduce the dimensions of the initial perturbation variables in step A.2. An appropriate dimensionality reduction method should be selected to transform the high-dimensional space problem into a low-dimensional space. Considering that the dimensions of the initial variables in ESMs are usually large, for example, in order to reduce the dimensionality of the solution space, the principal component analysis (PCA) method is used to reduce the dimensionality of the sample space. The formula is:

[0063] X′X′ T V′=V′Σ′

[0064] Where X′ is the initial perturbation sample, V′ is the obtained eigenvector, and Σ′ is the corresponding eigenvalue. Σ′ is sorted from largest to smallest, and the first n eigenvectors are selected to form the eigenvector matrix, which serves as the low-dimensional feature space for subsequent optimization.

[0065] Initial perturbation samples are obtained by performing equal-interval interpolation on the EMS simulation fields. This example samples the 20-year free-running CESM simulation results at daily intervals, extracting the winter (December to February) simulation fields, resulting in 90 × 20 = 1800 simulation fields. These winter simulation fields are then interpolated at 15-day intervals to generate (90 - 15) × 20 = 1500 initial perturbation samples. The Z-Score normalization method is used to map variables of different scales and ranges to a unified range to eliminate dimensional differences between variables. The calculation formula is:

[0066]

[0067] Among them, A min With A max Represent the minimum and maximum values ​​of each initial disturbance variable, A norm Represents the preprocessed sample. S3. Solve the orthogonal joint conditional nonlinear optimal perturbation model in the low-dimensional space using an intelligent algorithm to obtain multiple combinations of initial perturbations and pattern parameter perturbations, and map the initial perturbations in each combination back to the high-dimensional space using the eigenvector matrix obtained in step S2.

[0068] Typically, OC-CNOPs can be solved by descending the gradient of the objective function with respect to the perturbation. However, most complex ESMs lack corresponding adjoint patterns to obtain gradient information. To address this issue, intelligent algorithms can be used to solve OC-CNOPs, such as simulated annealing (SA), particle swarm optimization (PSO), genetic algorithms (GA), and differential evolution (DE).

[0069] This embodiment uses the particle swarm optimization algorithm PSO to solve the orthogonal joint conditional nonlinear optimal perturbation model until the iteration termination condition is met (the maximum number of iterations is reached). The particle update method is as follows:

[0070]

[0071] in, and are the position and velocity of the jth particle at the kth iteration, respectively. The final particle position is the disturbance. c1 and c2 are cognitive and social learning factors, respectively, and are usually set to 2. r1 and r2 are random numbers drawn from a uniform distribution in the interval [0,1]. X pb is the optimal position of the particle so far, X gb is the global optimal position of the entire particle swarm at the kth iteration, ω k is the coefficient of the speed at the kth iteration, which is initially set to 1 and gradually decreases in each iteration, ω max and ω min They are the maximum inertia weight and the minimum inertia weight, which are set to 1 and 0.5 respectively in this embodiment. max is the maximum number of iterations.

[0072] S4. Superimpose the model parameter perturbation obtained in step S3 and the initial perturbation after dimensionality increase on the model parameters and initial field of the Earth system model, and obtain multiple sets of perturbation prediction results through numerical integration. The specific process is as follows:

[0073] S401, selecting the first eight combinations of initial perturbations and model parameter perturbations obtained by solving step S3, superimposing the model parameter perturbations and the initial perturbations after dimension upgrade in each combination on the model parameters and initial fields of the Earth system model, and obtaining eight perturbation forecasts by numerical integration;

[0074] S402. In order to ensure that all perturbation biases are 0, it is also necessary to superimpose the perturbations that are opposite to the model parameter perturbations in each combination and the initial perturbations after dimensionality upgrade on the model parameters and initial fields of the Earth system model respectively. By numerical integration, 8 perturbation forecasts are also obtained, and finally 16 perturbation forecast results are obtained.

[0075] S5. Taking an average of the control forecast result and the multiple groups of disturbance forecast results obtained in step S4 to obtain an ensemble forecast result of the target event, wherein the control forecast result is the result of the natural integration of the earth system model without adding disturbance.

[0076] The evaluation indicators are formulated by fully considering the key characteristics of the target event to ensure that the evaluation indicators can comprehensively and objectively reflect the performance and effect of the ensemble forecasting method. In this embodiment, the evaluation indicators used include RMSE (root mean square error), the ratio of RMSE to ensemble spread (Spread), and the Talagrand plot.

[0077] Specifically, RMSE can measure the standard of the difference between the predicted value and the observed value, and is a key indicator of accuracy. The smaller the RMSE value, the higher the accuracy of the forecast result.

[0078] The ensemble dispersion is the root mean square error of each ensemble member relative to the ensemble average, and its calculation formula is:

[0079]

[0080] in, is the ensemble average forecast result, P n is the forecast result of the nth member. The ratio of RMSE to Spread of an ideal ensemble forecasting system should be close to 1.

[0081] The Talagrand graph is formed by sorting the N set prediction members in ascending order to form N+1 intervals, and then calculating the probability of the observation falling into the kth interval:

[0082]

[0083] Among them, S k represents the cumulative number of observations that fall into the kth interval, and M represents the total number of valid samples. For a reliable ensemble forecast system, observations should fall into each interval with equal frequency. Therefore, the flatter the Talagrand plot, the more reliable the ensemble forecast.

[0084] To verify the effectiveness of the above method, this example uses RMSE, the ratio of RMSE to spread, and the Talagrand plot as evaluation indicators, conducts statistical analysis on 10 NAO cases, and compares the prediction results of the control experiment (Control) and the five ensemble experiments. The comparison results are as follows: Figure 3 shown.

[0085] The control experiment (Control) is a single deterministic forecast without adding disturbances. The five ensemble forecast experiments are: ensemble forecast only for initial uncertainty (O-CNOP-Is-EF), ensemble forecast only for model parameter uncertainty (O-CNOP-Ps-EF), ensemble forecast that considers two uncertainties through different members (O-CNOP-(I / P)s-EF), ensemble forecast that considers two uncertainties at the same time in a simple superposition manner (O-CNOP-(I+P)s-EF), and the OC-CNOPs ensemble forecast (OC-CNOPs-EF) proposed in this application. Figure 3 As shown in (3a), OC-CNOPs-EF (i.e., the method proposed in this application) presents the lowest RMSE value; Figure 3 As shown in (3b), OC-CNOPs-EF showed the most significant improvement compared with the control experiment; Figure 3As shown in (3c), compared with several other ensemble forecasting methods, the RMSE and Spread ratio of OC-CNOPs-EF is closer to 1; Figure 3 As shown in (3d), the Talagrand plot of OC-CNOPs-EF shows a flatter distribution. The above results show that the proposed method OC-CNOPs-EF can more accurately reflect the uncertainty in the forecast, and its ensemble forecast results are more reliable.

[0086] If the above method is implemented in the form of a software functional unit and sold or used as an independent product, it can be stored in a computer-readable storage medium. Based on this understanding, the technical solution of the present invention, or the part that contributes to the prior art, or the part of the technical solution, can be embodied in the form of a software product. The computer software product is stored in a storage medium and includes several instructions for enabling a computer device (which can be a personal computer, server, or network device, etc.) to execute all or part of the steps of the method described in each embodiment of the present invention. The aforementioned storage medium includes: U disk, mobile hard disk, read-only memory (ROM, Read-Only Memory), random access memory (RAM, Random Access Memory), disk or optical disk, and other media that can store program code.

[0087] The above description of the embodiments is intended to facilitate understanding and use of the invention by those skilled in the art. It will be apparent that those skilled in the art can readily make various modifications to these embodiments and apply the general principles described herein to other embodiments without requiring inventive effort. Therefore, the present invention is not limited to the above-described embodiments. Improvements and modifications made by those skilled in the art based on the disclosure of the present invention, without departing from the scope of the present invention, should be within the scope of protection of the present invention.

Claims

1. A target event ensemble prediction method considering initial disturbances and model uncertainty, used to predict target events in any model of the Earth system model, characterized in that: The following steps are involved: S1. Based on the initial disturbance variables and pattern parameter disturbance variables corresponding to the target event, as well as their corresponding range constraints, an orthogonal joint conditional nonlinear optimal disturbance model is constructed; S2, collecting simulated field data of the Earth system model, generating initial disturbance samples, and obtaining eigenvector matrices after preprocessing; S3, solving the orthogonal joint conditional nonlinear optimal perturbation model in a low-dimensional orthogonal subspace to obtain a plurality of combinations of initial perturbations and pattern parameter perturbations, and mapping the initial perturbations in each combination back to the high-dimensional space through the eigenvector matrix obtained in step S2; S4, superimposing the model parameter perturbation obtained in step S3 and the initial perturbation after dimensionality upgrade on the model parameters and initial field of the Earth system model respectively, and obtaining multiple sets of perturbation prediction results through numerical integration; S5. averaging the control forecast result and the multiple groups of disturbance forecast results obtained in step S4 to obtain an ensemble forecast result of the target event, wherein the control forecast result is a result of a natural integration of the Earth system model without adding disturbances; In step S1, the objective function of the orthogonal joint conditional nonlinear optimal perturbation model is specifically expressed as follows: Among them, J is the objective function of the model, is the optimal joint perturbation, s 0i is the initial perturbation, p i is the mode parameter perturbation, S0 and P represent the original initial field and the default mode parameter value respectively, M t (P)(S0) represents the state of the Earth system model M integrated to time t when no disturbance is added. t (P+p i )(S0+s 0i ) indicates that when adding the joint perturbation (s 0i ,p i ), the Earth system model M is integrated to the state at time t, Ω i is a constraint condition; The constraint condition Ω i The specific expression is as follows: Among them, Ω i is the constraint condition for the i-th optimization, δ s is the constraint range of the initial perturbation, δ p is the constraint range of the mode parameter perturbation.

2. A target event ensemble forecasting method considering initial disturbance and pattern uncertainty according to claim 1, characterized in that: In step S2, the preprocessing includes normalization and dimensionality reduction.

3. A target event ensemble forecasting method considering initial disturbance and pattern uncertainty according to claim 2, characterized in that: The Z-Score method was used for normalization.

4. The target event ensemble forecasting method considering initial disturbance and pattern uncertainty according to claim 2, characterized in that: The principal component analysis method is used for dimensionality reduction.

5. The target event ensemble forecasting method considering initial disturbance and pattern uncertainty according to claim 1, characterized in that: In step S3, the orthogonal joint conditional nonlinear optimal perturbation model is solved in a low-dimensional orthogonal subspace by using a particle swarm optimization algorithm, a simulated annealing algorithm, a genetic algorithm or a differential evolution algorithm.

6. A target event ensemble forecasting method considering initial disturbance and pattern uncertainty according to claim 5, characterized in that: In the process of solving the orthogonal joint conditional nonlinear optimal perturbation model using the particle swarm optimization algorithm, the particle update method is as follows: in, and are the position and velocity of the jth particle at the kth iteration, respectively. The final particle position is the disturbance, ω k is the coefficient of speed at the kth iteration, which is initially set to 1 and gradually decreases in each iteration. c1 and c2 are cognitive and social learning factors, respectively. r1 and r2 are random numbers drawn from a uniform distribution in the interval [0,1]. X pb is the optimal position of the particle so far, X gb is the global optimal position of the entire particle swarm at the kth iteration.

7. A target event ensemble forecasting method considering initial disturbance and pattern uncertainty according to claim 6, characterized in that: The speed coefficient ω at the kth iteration k The specific expression is as follows: Among them, ω max and ω min are the maximum inertia weight and the minimum inertia weight, respectively, iter max is the maximum number of iterations, and k is the current number of iterations.

8. The target event ensemble forecasting method considering initial disturbance and pattern uncertainty according to claim 1, characterized in that: In step S4, the specific process of obtaining multiple sets of disturbance prediction results through numerical integration is as follows: S401, selecting the first n combinations of initial perturbations and model parameter perturbations obtained by solving step S3, superimposing the model parameter perturbation and the initial perturbation after dimension increase in each combination on the model parameters and initial field of the Earth system model, and obtaining n perturbation forecasts by numerical integration; S402. Superimpose the disturbances that are opposite to the model parameter disturbance and the initial disturbance after dimensional upgrading in each combination on the model parameters and initial field of the Earth system model respectively, obtain n disturbance forecasts through numerical integration, and finally obtain 2n disturbance forecast results.

Images