Model Uncertainty Representation for a Convection-Allowing Ensemble Prediction System Based on CNOP-P

Lu WANG, Xueshun SHEN, Juanjuan LIU, and Bin WANG,5

1Asset Operation Centre, China Meteorological Administration, Beijing 100081, China

2State Key Laboratory of Severe Weather, Chinese Academy of Meteorological Sciences, Beijing 100081, China

3Numerical Weather Prediction Center/National Meteorological Center,China Meteorological Administration, Beijing 100081, China

4Institute of Atmosphere Physics, Chinese Academy of Sciences, Beijing 100029, China

5Center for Earth System Science, Tsinghua University, Beijing 100084, China

ABSTRACT Formulating model uncertainties for a convection-allowing ensemble prediction system (CAEPS) is a much more challenging problem compared to well-utilized approaches in synoptic weather forecasting. A new approach is proposed and tested through assuming that the model uncertainty should reasonably describe the fast nonlinear error growth of the convection-allowing model, due to the fast developing character and strong nonlinearity of convective events. The Conditional Nonlinear Optimal Perturbation related to Parameters (CNOP-P) is applied in this study. Also, an ensemble approach is adopted to solve the CNOP-P problem. By using five locally developed strong convective events that occurred in pre-rainy season of South China, the most sensitive parameters were detected based on CNOP-P, which resulted in the maximum variations in precipitation. A formulation of model uncertainty is designed by adding stochastic perturbations into these sensitive parameters. Through comparison ensemble experiments by using all the 13 heavy rainfall cases that occurred in the flood season of South China in 2017, the advantages of the CNOP-P-based method are examined and verified by comparing with the well-utilized stochastically perturbed physics tendencies (SPPT) scheme. The results indicate that the CNOP-P-based method has potential in improving the under-dispersive problem of the current CAEPS.

Key words:CNOP-P,convective scale,model uncertainty,ensemble forecastforecast

1.Introduction

Convection-allowing models (CAMs) have been receiving increasing attention in recent years owing to their potential in describing the detailed characteristics of mesoscale systems (Baldauf et al., 2011; Seity et al., 2011; Clark et al.,2016). Usually, the grid size of a CAM is 1-4 km, but this range of grid spacing is regarded as the “physics grey zone”(Gerard, 2007; Honnert et al., 2016). Within this grid size range, the physical processes are partly resolved, and some remain parameterized. It is still not fine enough to resolve convective cells, subgrid-scale orographic gravity wave drag, as well as turbulence diffusion. These bring uncertainties to the model, which exist in association with physical parameter-izations and parameters in physical schemes. In addition,interactions between dynamics and physics or among different physical schemes in CAMs may become even more complicated owing to the increase in number of degrees of freedom in modeling the atmospheric phenomena at convection-allowing scales. Hohenegger and Schar (2007) pointed out that there is a much stronger nonlinear interaction between dynamics and physical processes owing to the rapidly developing character of local strong events. Thus,CAMs have their limitations in predicting the initiation and development of strong events only by using a single deterministic forecast, including their location, intensity and duration.

Due to its effective probabilistic forecast information,the ensemble approach is believed necessary for the purpose of fully exploiting the potential of CAMs. Indeed, convection-allowing ensemble prediction systems (CAEPSs)play an important role in local high-resolution numerical weather prediction. Several leading forecast centers have started to run a CAEPS operationally in recent years, including COSMO-DE-EPS (Peralta et al., 2012), MOGREPS-UK(Hagelin et al., 2017), and AROME-EPS (Bouttier et al.,2012). CAEPSs have also been developed as part of NOAA’s Hazardous Weather Spring Experiment, such as SSEF (Clark et al., 2012) and NCAR’s EnKF-based EPS(Schwartz et al., 2015). See Table A1 for further details.

However, the above CAEPSs have a common problem,known as under-dispersion. Typically, as the forecasts proceed, the ensemble spread among the ensemble members becomes too small to ensure that the observed weather pattern could be contained within the ensemble (Mascaro et al.,2010; Weyn and Durran, 2018).

Previous studies have pointed out that CAEPSs could benefit through improving the representation of model uncertainties (Johnson et al., 2014; Romine et al., 2014). To date,several approaches have been proposed in representing model uncertainties, including multi-model, multi-physics,multi-parameter, and stochastic physics. The multi-model approach considers the uncertainty of the dynamical core and parameterization by using different models (Iversen et al., 2011; Melhauser et al., 2017). Multi-physics and multiparameter methods are used to describe the uncertainties embedded in the physical parameterizations by changing schemes (Charron et al., 2010; Hacker et al., 2011) or their empirical parameters (Gebhardt et al., 2011). Stochastic physics is more theoretical compared with the above approaches,in which the model uncertainties are described by introducing random perturbations into the physical processes, such as the well-utilized Stochastic Perturbed Parameterization Tendencies (SPPT) scheme (Buizza et al., 1999;Christensen et al., 2017), the Stochastic Total Tendency Perturbation scheme (Hou et al., 2006, 2008), and the Random Parameters/Stochastically Perturbed Parameterizations scheme (Bowler et al., 2008; McCabe et al., 2016; Ollinaho et al., 2017). In addition, considering the kinetic energy dissipation at the truncation scale, the Stochastic Kinetic Energy Backscatter scheme (Shutts, 2005; Berner et al., 2009) is another stochastic method that expresses the model structural uncertainty. Recently, a new approach based on the nonlinear forcing singular vector has been proposed to represent the combined effect of different sources of model uncertainties (Tao and Duan, 2019; Qin et al., 2020).

The above-mentioned model perturbation approaches have been widely utilized in global ensemble forecasts for a long time, and recently applied in CAEPSs for forecasting severe convective weather, tropical cyclone intensity, and so on. However, it is found that the above model perturbation approaches have limited effects in dealing with the under-dispersion problem of CAEPSs (Fresnay et al., 2012;Baker et al., 2014). Such limited effects may come from deficiencies in formulating the model uncertainties at the convection-allowing scales. Compared with global ensemble prediction systems, the representation of model uncertainty in CAEPSs is thought be lacking in systematic research and a theoretical basis, therefore making it an important issue worthy of further study.

As known, the main target of a CAEPS is to capture the probability of occurrence of rapidly evolving small- to mesoscale strong events. From the ensemble forecast viewpoint,the rapidly nonlinear growing perturbations should be described in the CAEPS’s formulations both for initial errors and model uncertainties. Mu et al. (2003) proposed the concept of conditional nonlinear optimal perturbation(CNOP). The CNOP represents the maximum nonlinear evolution of initial perturbations satisfying constraint conditions at the prediction time, which is thought to be greatly important for the predictability (Wang and Tan, 2010).Since Mu et al. (2003), CNOP has been applied in studying the various predictability issues, such as ENSO predictability, target observations, and so on. Duan and Huo (2016)also proposed the approach of orthogonal CNOPs, which extended the CNOP to ensemble forecasting. Furthermore,considering the predictability caused by the model error,Mu et al. (2010) extended the CNOP approach to search for the optimal model parameter perturbations, which resulted in the largest departure from a given reference state at a prediction time with physical constraints. This approach is referred to as CNOP-P, which has been extensively employed in understanding the climate model sensitivity to model parameters and investigating the impact of parameter errors on ENSO predictability, amongst other problems (Duan and Zhang, 2010; Yu et al., 2012; Yin et al.,2015).

As discussed above, a good CAEPS should reasonably describe the rapidly growing nonlinear errors or maximum growing perturbations at a short prediction time either due to initial errors or model uncertainties. In this sense, the CNOP-P approach may have potential to detect and describe the maximum perturbation growth caused by the model parameter uncertainties for CAEPSs. Compared with current formulations of model uncertainties in CAEPSs, the CNOP-P approach is expected to be more objective.However, the original computation method of CNOP-P requires integration of an adjoint model to calculate the cost function gradient. This limits the practical use of CNOP-P in operational applications. Wang and Tan (2010) proposed a pragmatic method to determine the CNOP by using an ensemble approach, which is easier to implement and more time-saving than the adjoint-based method. Thus, our research seeks to exploit a novel method to formulate the model uncertainties for CAEPSs by using the CNOP-P concept and the approach to determining the CNOP proposed by Wang and Tan (2010).

The rest of the paper is organized as follows. Section 2 briefly describes the CNOP-P approach utilized in this paper and the experimental design. The results of sensitivity tests and the impact of the CNOP-P approach on the ensemble forecasts are reported in section 3. Finally, conclusions and future work are summarized in section 4.

2.Methods

2.1.CNOP-P and its computation

According to Mu et al. (2010), in a nonlinear prediction modelM,can denote a prediction from the initial time 0 toτwith the parameter vector P and the initial value U0. The prediction increment y' at the prediction timeτis defined as the departure from the reference state as a result of parameter perturbations p':

The optimal parameter perturbations p' should correspond to maximum the prediction increments y' with a constraint of || p'||≤σ(σis a small positive number, given here as 1). The p' can be obtained through solving the following equation:

In order to reduce the dependency ofJ(p') on initial values, Yin et al. (2015) introduced the ensemble-averaged type of cost function, which is constructed by usingNshortterm time series predictions with independent initial values U0,k(k=1,2,…,N) as follows:

The above cost function can be rewritten into a minimization problem for convenience of solution:

The key to this minimization problem is the calculation of its gradient with respect to the parameter perturbation,which can be expressed as:

As derived by Wang and Tan (2010) and Yin et al.(2015), the solution of Eq. (5) can be obtained by the following ensemble approach.

Firstly, for each initial value U0,k,nparameter perturbation samples p'i(i=1,2,…,n) are used to generatencorresponding prediction increment samples y'i,k(i=1,2,…,n;k=1,2,…,N). The parameter perturbations and prediction increment matrices are denoted as P' and Y'k, respectively:

The second step is to obtain a series of orthogonal basis functions of P' by using singular value decomposition:

where, Apis the left singular vector of P', which is taken as the orthogonal basis function.

The third step is to estimate the gradient of the cost function based on ensemble projection. With Apas the projection matrix, each parameter perturbation can be projected to the ensemble space:

whereαis the parameter perturbation after projection.

According to the tangent linear relationship between p'and y'k, each prediction increment can also be projected to the same ensemble space:

whereM'is the tangent linear model ofM.

Then, a statistical model between p' and y'kcan be estabished:

The last step is to compute CNOP-P and detect the most sensitive parameters. The Spectral Projected Gradient 2 (SPG2) optimization algorithm of (Birgin et al., 2001) is used to obtain a result iteratively. To guarantee nonlinearity of the optimal solution, GRAPES-Meso (Global/Regional Assimilation and Prediction System, mesoscale version),with a horizontal grid size of 3 km, is used to calculate the cost function value. Meanwhile, the gradient of the cost function and constraint condition of parameter perturbations are required during iteration. The constraint condition is expressed as:

There aremparameters involved in the calculation of CNOP-P. Perturbations are standardized by dividing their corresponding default values. In this way, the greater the perturbation weight of a parameter in the CNOP-P, the greater the relative sensitivity of this parameter.

2.2.Experimental design

2.2.1.Model configuration

The model utilized here is the operational GRAPESMeso, version 4.1, with a convection-allowing resolution setup. This model solves the fully compressible nonhydrostatic equations using a semi-implicit, semi-Lagrangian integration scheme on a spherical Arakawa-C grid horizontally,and Charney-Phillips staggered grids vertically. It has 50 vertical levels up to 10 hPa. For the 3 km resolution, the 30-s time step is adopted.

Physical parameterization schemes used here include the Rapid Radiation Transfer Model longwave radiation scheme (Mlawer et al., 1997), the Dudhia shortwave radiation scheme (Dudhia, 1989), the Monin-Obukhov surface layer scheme, the Noah land surface scheme (Chen and Dudhia, 2001), the MRF boundary layer scheme (Hong and Pan,1996), and the WSM6 cloud microphysics scheme (Hong and Lim, 2006). The cumulus convection scheme is not applied.

The initial conditions are generated by the GRAPES Cloud Analysis System (GCAS) (Qu et al., 2010) with a 3-km resolution, which has been tested in GRAPESMeso_3km from 2015 to the present (Zhu et al., 2017). By using GCAS, less model spin-up time in the precipitation forecast was observed. NCEP Final Operational Global Analysis data at a spatial resolution of 1° are used as the background. Based on this background, GCAS integrates data sources from the Doppler weather radar 3D mosaic reflectivity data at a spatial resolution of 0.01°, blackbody temperature and total cloud amount data from the FY-2G geostationary satellite. The lateral boundary conditions are also provided by the NCEP Final Operational Global Analysis with a 1° horizontal resolution.

Because CAEPSs are mainly designed for obtaining the forecast probability of convection, a region with frequent occurrence of locally developed convection is preferred as the modeling area. As is known, most of the precipitation during the pre-summer rainy season in South China is of a convective nature with mesoscale organizational characteristics.Therefore, the model domain covers South China(18.5°-27.5°N, 105.0°-118.5°E) in this study. In order to reduce the influences from the boundaries on the CNOP-P computation, a smaller area (20°-26°N, 106.5°-117°E) is selected for calculating the CNOP-P, as shown in Fig. 1.

2.2.2.Experimental design for detecting the sensitive physical parameters based on CNOP-P

For the purpose of better isolating the physical processes influencing the initiation and development of convection, five locally developed rainstorm cases in 2017 that occurred in South China are chosen for constructing the ensemble to compute the CNOP-P, i.e.,N=5 in Eq. (3).These five cases occurred under weak synoptic forcing.Here, the 12-h cumulative precipitation is used to calculate the cost function for detecting the most sensitive parameters, which means that CNOP-P could lead to the maximum variations in precipitation during the first 12 h of model integration. Details of the five cases are given as Cases 1-5 in Table 1. Then, a supplementary experiment is carried out to confirm the robustness of the above parameter detection.Five heavy rainfall cases under strong synoptic forcing are selected to compute the CNOP-P, which are given as Cases 6-10 in Table 1. It is worth mentioning that all the forecasts in Table 1 are initialized about 6 h before the onset of rainfall events. These cases are classified as well as documented by the forecasters in the Central Meteorological Bureau (Zhang et al., 2018).

Fifteen parameters as shown in Table 2 are selected from the boundary layer and microphysics, which are thought to be the most important two processes in influencing convective events. Theoretically, all the physical parameters should be selected to perturb for calculating the CNOP-P. However, due to the high computing costs, we only selected 15 parameters in this study. According to previous research on CAEPSs (Baker et al., 2014; Gebhardt et al., 2011; McCabe et al., 2016), it was believed that this selection would not lose the generality. When we generate the parameter perturbation samples, only one parameter is changed each time, and the perturbed value is ±50% of its default

Fig. 1. Model forecast domain (D1: 18.5°-27.5°N,105.0°-118.5°E) and computation domain for CNOP-P (D2:20°-26°N, 106.5°-117°E).

Table 1. Twenty-three typical heavy rainfall cases that occurred in South China in 2017.

value.

2.2.3.Ensemble experiment design

By applying the CNOP-P, a group of the most sensitive parameters can be detected from the 15 parameters,which correspond to the maximum forecast error growth caused by the parameter uncertainties. Based on these most sensitive parameters, a novel formulation of model uncertainty is designed, and the ensemble experiment design is given below.

As is known, currently existing systems typically use ensemble sizes of around 10 or more members, like MOGREPS-UK, AROME-EPS, NCAR’s EnKF-based EPS,etc. Due to the available computer resources, a 10-member ensemble run was designed in this study, which included one control member and nine perturbed members. For the purpose of investigating the advantages of this new approach in describing model uncertainty, we established the ensemble prediction system by only taking into account the model uncertainty. The well-utilized SPPT scheme is adopted as the benchmark, which represents the typical model uncertainty formulation in current operational global ensemble and regional ensemble systems (Bouttier et al., 2012; Leutbecher et al., 2017). Two sets of ensemble forecasts are carried out with different model perturbation methods. The first set of ensemble forecasts uses SPPT scheme, which is called EXP1 for short. The second uses the CNOP-P-based approach, which is called EXP2 for short. Next, the details of the two model uncertainty formulations are described.

In the SPPT scheme, a stochastic perturbation is added to the net tendency term of physical parameterizations:

Table 2. The chosen physical parameters for computing CNOP-P.

where X0is the net tendency term, Xpis the perturbed net tendency term, andφis a 3D random field. The scheme firstly introduces a random field generator with certain space structures and time-correlated characteristics. This generator is an expansion of spherical harmonics along the horizontal directions based on first-order autoregressive random processes. Then, the random number is multiplied into the net tendency of all the physical processes. Aimed at maintaining the numerical stability, vertical profiles are introduced to reduce the perturbations near the surface and the model top. The following parameter configurations are based on Yuan et al. (2016): the random field follows a Gaussian distribution with a mean value of 0, a standard deviation of 0.27, a perturbation range of [0.2, 1.8], and temporal decorrelation scale of 6 h. Because of time limitations,this paper does not consider the tuning of the SPPT scheme,or its reformulation, e.g., by using the iSPPT formulation(Christensen et al., 2017). It should be noted that the SPPT scheme used here has been well tuned for the GRAPES mesoscale ensemble forecast system with a 15-km resolution,and has been operationally implemented since 2014 (Yuan et al., 2016).

For the CNOP-P-based approach, based on the detected sensitive parameters as described above, a random parameter algorithm (McCabe et al., 2016) is introduced to construct the CNOP-P-based formulation of the model uncertainty.

By giving a physically reasonable range [Pmin,Pmax] for varying to each parameter, the detected sensitive parameters are first mapped onto the range [?1, 1]:

whereηis the mapped value, andPminandPmaxare the minimum and maximum values.

By using the first-order auto-regression model, the mapped parameters are stochastically perturbed and updated at time stept. The parameters are perturbed stochastically and independently at regular intervals throughout the relations:

whereηtis the value of mapped parameters at timet,ris the autocorrelation coefficient, andεtis the stochastic shock term.εtis formulated in terms ofrand stochastic numberAt.ris expressed by the update time interval Δtand the characteristic time scaleτp. According to McCabe et al. (2016), the update time interval is set to 5 min and the characteristic time scale is set to 7 d. Then,ηtwill be mapped back to original range using Eq. (14).

2.2.4.Verification and metrics

For the purpose of systematically comparing the forecast skills using the two different approaches of model uncertainty formulation, 13 heavy rainfall cases that occurred during the flood season of South China in 2017 were selected to carry out reforecast experiments. These cases are given as Cases 11-23 in Table 1, including heavy rainfall events with strong and weak synoptic forcing.

The near-surface variables, including 2-m specific humidity (Q2m), 10-m wind speed (V10m) and 2-m temperature(T2m), as well as hourly precipitation (R), are mainly evaluated in this study. Hourly observations from automatic weather stations are used in the verification.

The verification metrics include the ratio of spread to root-mean-square error (RMSE), outlier frequency, and continuous ranked probability score (CRPS). In an ideal ensemble, the ratio of spread to RMSE (hereafter, S/R)should be as close to 1 as possible; that is, the spread can represent the error evolution characteristics of the system. In reality, spread is often less than the RMSE owing to the observation error. Here, observation error is not taken into account in the verification in our study:

where (ˉ) denotes the mean value of all grid points,fiis the forecast of ensemble memberi,nis the number of ensemble members,fmeanis the ensemble mean,f0is the control forecast, andois the observation value.

Outlier frequency is used to evaluate the reliability of the CAEPS. Suppose there areNvalid samples andnensemble members. After sorting the forecast of ensemble members from small to large,n+1 intervals are generated.Sidenotes the occurrence frequency of observation in theith interval. The Talagrand distribution map can be drawn from the probability distribution:

In an ideal Talagrand distribution, there is not much difference between the probabilities of each interval, but in reality, the two ends are larger. Outlier frequency is a measure of the value at both ends of the distribution, which denotes the average probability that the observations fall outside the range of the forecasts. The ideal value is 1/(n+1).

CRPS is a probabilistic skill score to measure the differences in cumulative distribution functions between observations and forecasts. It can be expressed as the integral of the Brier score for all possible thresholds of the meteorological variables:

whereF(x) andFo(x) represent the cumulative distribution functions of the ensemble forecast and observation, respectively. The value for an ideal forecast of CRPS is 0. A lower CRPS indicates a more skillful forecast.

3.Results

3.1.Sensitive parameter detection

As described in section 2.1, iteration is necessary in solving Eq. (5) for detecting the sensitive parameters. In this study, the iteration is stopped when the cost function becomes almost unchanged (< 1%), and the spectral projection gradient decreases by five to six orders of magnitude.Table 3 shows the results of the iterative process. It is found that the above criteria start being met from step 5 of iteration, and then the parameters remain stably ranked. Here,the top eight parameters are regarded as the most sensitive ones, denoted as P1, P13, P7, P5, P3, P10, P4 and P2 in Table 2. We attempt to provide a physical explanation of these eight sensitive parameters.

In the MRF scheme, the PBL and free atmosphere vertical diffusion is calculated based on the local K approach.The asymptotic mixing length (P1) is a physical quantity describing the size of the large energy-containing eddies in turbulence mixing, and is thus crucial in determining the vertical diffusion intensity. The mixing length scale is detected as the number one sensitive parameter, illustrating the rather important role of boundary layer processes in affecting locally developed convections. P2, the critical Richardson number, is a key parameter for calculating boundary layer height. P3 is a profile shape exponent for calculating the momentum diffusion coefficient, which determines the magnitude of the vertical diffusion as well as the height at which the diffusion reaches a maximum. P4, the Prandtl number coefficient, is involved in computing eddy diffusivity for temperature and moisture within the boundary layer. In the WSM6 microphysics scheme, the constants for the terminal velocity of rain and snow (P5, P7) play important roles in the vertical distribution of hydrometeors and auto-conversion in clouds. P10, the density of graupel, varies greatly,and can lead to obviously different characteristics of cloud microphysics. Setting P10 as a constant may underestimate or overestimate the precipitation efficiency of clouds. The maximum value of the cloud ice diameter, named P13 here,is related to the maximum mass of ice crystals, which is critical in auto-conversion from ice to snow, and further influences the rainfall intensity for the strong convection cases.

Table 3. The iteration of CNOP-P computation.

These eight sensitive parameters are applied to address the model uncertainty using Eqs. (14)-(17). Table 4 gives the range of perturbations for these parameters. It should be noted that this range is a little different from the perturbed range in the calculation of CNOP. To ensure the representativeness and independence of samples, we perturbed ±50%of the defaults to generate the parameter perturbation samples in computing CNOP-P. Here, however, the perturbations need to be closer to the actual values in the whole variation ranges, which should be large enough to add reliable variability into the forecasts, and small enough to make the forecasts reasonable. Although these exact thresholds are not clear, we define the values by referring to previous research (McCabe et al., 2016; Ollinaho et al., 2017).

The following are the results of the supplementary experiment, which, using five heavy rainfall cases under strong synoptic forcing, was in sharp contrast to the above in which we selected five cases with weak synoptic forcing. The sensitive parameters detected by this experiment are P1, P13, P5,P7, P2, P3, P14 and P12, as shown in Table 5. Only two parameters changed compared with the experiments shown inTable 3. The result indicates that the detection and sensitivity ranking of parameters is relatively stable and robust.

Table 4. Sensitive parameters and their ranges.

3.2.Evaluation of the ensemble forecasts

The purpose of this paper is to introduce the CNOP-P-based model uncertainty formulation and understand its advantages for CAEPSs. For a relatively fair comparison,we chose a total of 13 heavy rainfall cases to understand the systematic effects, and calculated the averaged differences of metrics over a forecasting time of 6-24 h between EXP1 and EXP2, i.e., EXP2 minus EXP1, to investigate the added value of the CNOP-P-based approach. According to the definition of these metrics, EXP2 performs better when the difference of S/R is greater than 0, or the differences of outlier frequency and CRPS are less than 0, and vice versa. Figure 2 shows the 6-24-h averaged differences of S/R on four variables for the 13 cases, withT2m,V10m,Q2mandRshown in Figs. 2a-d, respectively. From all the positive values of S/R differences in Fig. 2, it is obvious that the CNOP-P-based method, i.e., EXP2, exhibits systematic improvement with respect to the S/R relationship for all the near-surface variables and precipitation. This result indicates that the CNOPP-based method has great potential in improving the underdispersive problem of current CAEPSs.

Similar to Fig. 2, Fig. 3 gives the 6-24-h averaged differences of outlier frequency for each variable. As mentioned above, outlier frequency denotes the average probability that the observations fall outside the range of the forecasts.As shown in Fig. 3, the time-averaged differences of outlier frequency forT2m,V10mandQ2mexhibit negative values for all 13 cases. Only two exceptions are found in the precipitation outlier frequency difference (Fig. 3d). It can be concluded that, on average, the CNOP-P-based method can reduce the probability that the observations fall outside the range of the forecasts, thus enhancing the reliability of the CAEPS.

The CRPS is frequently used to assess the respective accuracy of two probabilistic forecasting models. The 6-24-h time-averaged differences of CRPS for the two probabilistic forecasts are shown in Fig. 4 for each variable. Accord-ing to the definition of CRPS, a negative time-averaged difference means a better probabilistic forecast. From Fig. 4, it can be found that, for all four variables, about 11 out of 13 cases show negative differences, indicating a better forecast corresponding to the cumulative distribution function of the observations.

Table 5. The iteration of CNOP-P computation using five cases under strong synoptic forcing as initial conditions.

Fig. 2. The 6-24-h time-averaged differences of S/R (spread/RMSE) for 13 cases. Black lines represent EXP2’s S/R on four variables, and columns represents the difference between EXP2 and EXP1 (i.e. EXP2 minus EXP1). (a) 2-m temperature; (b) 10-m wind speed; (c) 2-m specific humidity; (d) hourly precipitation.

Fig. 3. As in Fig. 2 but for outlier frequency.

Fig. 4. As in Fig. 2 but for CRPS.

As shown in Figs. 2-4, the above verification results indicate that the CNOP-P-based method can bring more spread, and has relatively reliable and skillful probabilistic forecasts compared with the SPPT scheme. That is, perturbing the most sensitive parameters or physics may have more potential to improve CAEPSs than perturbing the net physics tendency.

Next, we check the different effects of the two approaches in the vertical direction. Figure 5 shows the vertical profiles of ensemble spread for the 24-h forecast of zonal wind, temperature and specific humidity, the profiles of which are the averages of the 13 cases. For the zonal wind field in Fig. 5a, EXP2 can obtain larger spread within the boundary layer, while EXP1 seems to be more dispersive in the middle atmosphere. This may partly reflect the characteristics of the two approaches in formulating the model uncertainties, but on the whole the difference is not significant. For temperature and specific humidity in Fig. 5b and Fig. 5c, EXP2 has consistently and statistically significant larger spread over the troposphere compared with EXP1. This again confirms the merits of the CNOP-P-based approach,because through detecting the most sensitive parameters from CNOP-P, the process-specific and event-oriented model uncertainties can be realized. More specifically, focusing on the CAEPS of convective events in this study, the most sensitive parameters embedded in the PBL and microphysics are detected by using CNOP-P. Thus, the formulation of model uncertainty of these two key processes for convection, which has solid theoretical basis, can bring more spread for temperature and humidity over the troposphere.

Finally, it should be noted that some ensembles might not necessarily be under-dispersive in terms of all parameters, so that our study may be more about creating “a right kind of spread ” rather than increasing the spread in absolute terms. A sensitivity test was conducted to explain why we chose to perturb the first eight sensitive parameters in the ranking of Table 3. Five comparison experiments were implemented, each of which perturbed the first 4, 6, 8, 10 and 12 parameters in the sensitivity ranking of Table 3. We chose Case 1 in Table 1 for ensemble prediction, and calculated the differences of three metrics over the 6-24-h forecasting time (Figs. 6-8). From the results of the three verifica-tion metrics, we found that the perturbation scheme with eight parameters is more effective than that with four or six parameters. However, increasing to more than eight parameters does not obviously improve the forecasts. It can be seen that the uncertainty of parameters is mainly due to the important sensitive parameters.

Fig. 5. Vertical profiles of ensemble spread for the 24-h forecasts of 13 cases: (a) zonal wind (units: m s?1); (b)temperature (units: °C); (c) specific humidity (units: g kg?1).

4.Summary and remarks

The model perturbation methods in CAEPSs lack consideration of the rapidly nonlinear growing errors at the convective scale, which is partly responsible for the under-dispersion problem of CAEPSs. This paper seeks to address this problem by designing a more reasonable formulation of model uncertainty based on the CNOP-P approach. The CNOP-P can search for the optimal model parameter perturbations, which result in the largest departure from a given reference state at a prediction time with physical constraints.Focusing on locally developed strong convective events that occurred in the pre-rainy season of South China, the most sensitive parameters were detected using CNOP-P, which resulted in the largest departure of heavy rainfall forecasts for five heavy rainfall cases. Then, a formulation of model uncertainty was designed by adding stochastic perturbations into these sensitive parameters. Comparison ensemble experiments between the CNOP-P-based method and the well-utilized SPPT scheme were conducted, in which the ensemblecost funtin.system only considered the model uncertainties. These experiments were carried out using all the 13 heavy rainfall cases that occurred in the flood season of South China in 2017.Our conclusions and remarks are as follows:

The CNOP-P approach is very robust in detecting the sensitive parameters when applied to the convection-allowing-scale model.

The CNOP-P-based method has potential in improving the under-dispersive problem of current CAEPSs, and has more reliable forecast skill forT2m,V10m,Q2mand precipitation compared with the SPPT scheme.

The CNOP-P-based method brings more spread for humidity and temperature over the troposphere owing to its process-specific and event-oriented formulation of model uncertainties.

Although the cases we selected in this study are specific to South China and may therefore limit the conclusions, the experiments with these 13 heavy rainfall cases in the 2017 flood season strongly suggest a very bright future in terms of the application of the CNOP-P-based approach in CAEPSs.

There is, however, still a practical and critical issue to resolve in terms of how to select the parameters if two of the eight parameters are different when different cases are used. The CNOP-P approach has been proven theoretically in previous studies and in some climate model applications,and it is believed in practical NWP applications that a stable group comprising the most sensitive parameters could be obtained if the cases are chosen to share the similar physical features or processes. The supplementary experiment in our study illustrates this issue, i.e., the parameters based on the “strong forcing ” cases were similar to those based on the “weak forcing” cases. Much research is needed in which the CNOP-P approach is further applied in operational CAEPSs-such as, its application in heavy rainfall events in the mid-to-high latitudes, the impact of parameter perturbation constraints, the use of spatiotemporal random perturbation approach when perturbing the parameters, and the choice of cost functions.

Fig. 6. The 6-24-h time-averaged S/R (spread/RMSE) for the case initialized at 1200 UTC 6 May.EXP_4 (red), EXP_6 (blue), EXP_8 (black) EXP_10 (grey), and EXP_12 (green) are five ensemble forecasts that perturbed the first 4, 6, 8, 10 and 12 parameters, respectively: (a) 2-m temperature; (b) 10-m wind speed; (c) 2-m specific humidity; (d) hourly precipitation.

Fig. 7. As in Fig. 6 but for outlier frequency.

Fig. 8. As in Fig. 6 but for CRPS.

Acknowledgements.We sincerely appreciate the constructive comments and suggestions of the anonymous reviewers. This work was supported by the National Key Research and Development Program of China (Grant No. 2017YFC150 1904).

APPENDIX

Table A1. Summary of operational convection-allowing ensemble prediction systems.