Sensitivity Analysis of Ensemble Simulations on a Torrential Rainfall Case over South China Using Multiple PBL and SL Parameterizations

The exchange of moisture, heat and momentum through turbulent mixing that occurs between the earth's surface layer (SL) and the free troposphere is vital for the initiation and development of convective clouds and associated precipitation ^[1-9]. Synoptic-scale processes can now be explicitly anticipated reasonably well by numerical weather prediction (NWP) models, while such vital turbulent processes within the planetary boundary layer (PBL) are still primarily implicitly represented by parameterization. The low representation of lower-troposphere turbulence is one of the major sources of forecast inaccuracy in NWP models ^{[3, 8, 10-13]}. Therefore, a skillful parameterization of PBL turbulence is of great importance for the prediction of clouds and convective activities.

PBL parameterization in both weather and climate models has been extensively explored for decades. The forecast skill of NWP models has been demonstrated to be sensitive to the formulation of vertical mixing. In addition, to obtain a good hurricane simulation, PBL schemes have been found to be of the same importance as the schemes of cloud microphysics ^[14-16]. Many PBL parameterizations have been developed to represent the turbulence process and its impacts on global or regional NWP models ^[17-22]. For example, in the latest version of the well-known Weather Research and Forecasting Model (WRF ^[23]), there are 13 PBL parameterization scheme options. Each of these options is formulated by different prognostic variables, diagnostic variables and cloud mixing variables^{[16, 24-32]}. These options provide multiple mechanisms for the model to characterize the effects of subgrid turbulent processes in the PBL.

SL parameterization, on the other hand, is closely linked to PBL parameterization. It characterizes the surface properties and inputs them into the PBL. In the WRF model, there are multiple parameterizations for different suites of PBL and SL schemes. They are of great importance for the initiation and development of clouds, especially in a convective environment. Most previous studies mainly focus on the impacts of PBL schemes on numerical simulations while seldom discuss the influence of SL parameterization. In the present study, the sensitivity of the model results to SL parameterization is also investigated to assess the relative contributions of different SL parameterizations to the simulation of a torrential rainfall case.

The structure and height of PBL have strong seasonal and diurnal variations in China ^[33]. The PBL height is generally higher in spring and summer than that in autumn and winter. This seasonal variation is found to be negatively correlated with the seasonal variations in surface pressure and lower-troposphere stability and positively correlated with the seasonal variations in near-surface wind speed and temperature. The height of the PBL tends to peak in the early afternoon, and the diurnal amplitude of the PBL height is higher in the northern and western subregions of China than in other subregions. In addition, the PBL height is also largely influenced by terrain features and major land-water boundaries ^[33]. These findings suggest that the description of the turbulence processes within the PBL should differ seasonally and regionally. The selection of PBL and SL parameterizations in NWP models may vary according to the regional underlying surface and thermodynamic and kinetic conditions.

The impacts of PBL and SL parameterizations on regional weather simulations have been extensively investigated, including over the southern Great Plains in the United States ^[34], the tropical ocean ^[35], northern China ^[36], western China ^[37], and the Yangtze and Huai River Basins ^[38]. In south China, which is characterized by a complex underlying surface, the annual precipitation amount reaches ~2000 mm ^[39], and the torrential rainfalls during the annually first rainy season (April to June) are usually associated with the PBL process. However, the impacts of PBL and SL parameterizations on regional torrential rainfall over south China during the annually first rainy season have rarely been investigated, which was critical for the accurate numerical prediction for the torrential rainfall to prevent loss of life and property caused by it.

This work aims to investigate the influence of the different suites of the PBL and SL parameterizations available in the WRF model on a torrential rainfall case occurred in south China and to identify the sources of bias and errors that affect the initiation and development of convection. The rest of this paper is organized as follows. Section 2 describes the data, NWP model and the evaluation methods used in this study. An overview of the torrential rainfall case is introduced in section 3. The results of the evaluation and diagnosis are presented in section 4. Finally, the conclusion is given in section 5.

The national radar mosaics produced by the Chinese Meteorological Administration (CMA) were used to analyse the evolution of the precipitation system. The hourly rainfall estimates, sea level pressure, 2-min averaged surface wind speed and direction, and the surface temperature were obtained from the quality-controlled automated weather stations (AWSs) provided by the CMA. There were 2910 AWSs in south China (20° - 26° N, 109° - 118° E), with an average distance between them of approximately 10 km, and the observations are available every 1 h ^[40]. The NCEP FNL (Final) Operational Global Analysis data (gridded on 1°× 1°) provide atmospheric variables four times a day (0000, 0600, 1200 and 1800 UTC).

The Advanced Research WRF Model (version 4.0) ^[23] was used in this study. Simulations were performed over two-way nested domains with the horizontal grid spacing of 9 km in the innermost domain (Fig. 1). There were 33 terrain-following hydrostatic-pressure vertical levels in all domains, with a top level of 50 hPa. The model was initialized at 0800 BJT on 17 April 2011 with the initial and lateral boundary conditions provided by the NCEP FNL analysis. Both domains used the same physical parameterization schemes, including the WRF double-moment 6-class microphysics scheme ^[41], the rapid radiative transfer model with GCM applications (RRTMG) longwave radiation scheme ^[42], and the Grell-Devenyi ensemble cumulus parameterization scheme ^[43]. To investigate the sensitivity of different PBL and SL schemes in this case, the available suites of PBL and SL parameterization schemes were chosen to form an ensemble simulation (Table 1).

Figure 1.  The two domains for the WRF simulations. Terrain heights (shaded; metres above the mean sea level) are also shown.

The performance of the ensemble was examined using objective evaluation methods, including the threat score (TS), equitable threat score (ETS), probability of detection (POD), miss ratio (MISS), false alarm ratio (FAR), and root-mean square error (RMSE). The TS measures the fraction of observed and / or forecasted events that are correctly predicted. Its value depends on the climatological frequency of events (poorer scores for rarer events) since some hits can occur purely due to random chance. As an improved TS, the ETS measures not only the fraction of observed and / or forecasted events that are correctly predicted but also adjusts for hits associated with random chance. The ETS is unbiased because it penalizes both misses and false alarms in the same way. The POD is the ratio of correct forecasts to the total number of observations. The MISS measures the fraction of missed forecasts to all forecasts of the event, while the FAR answers the question of what fraction of the predicted hit events do not actually occur. The higher the values of the TS, ETS, and POD and the lower values of the MISS and FAR indicate a better forecast. The RMSE measures the difference between the observation and the forecast using the same equation as Eq (1) in Cohen et al. ^[8], which is formulated as Model forecasts and observations are indicated by ${Y_t^s}$ and ${Y_t^a}$, respectively, while T represents the domain grid size. Specifically, the RMSE can theoretically be as low as zero, representing a perfect forecast, with progressively higher values indicating poorer forecast quality. In the present study, the simulation results were interpolated onto the AWSs using the Cressman interpolation method^[44] to calculate the domain-wide (20°- 26° N, 109° - 118° E) TS, ETS, POD, MISS, FAR and RMSE of the rainfall accumulation. The simulation that deviates the most from the observation will be analysed to investigate the error sources of the PBL and SFL schemes while forecasting torrential rainfall in this case, both from a dynamical and thermodynamic perspective.

The case took place from 0900 to 2300 BJT(Beijing Time; BJT = UTC + 8 h) on 17 April 2011. It was a typical torrential rainfall event that occurred over Guangdong Province (Fig. 2). The focused mesoscale convective system (MCS) originated from the northeast of Guangxi Province in the early morning at 17 BJT on April 2011 (not shown). It moved to the northwest border of Guangdong Province at 0900 BJT on 17 April 2011 (Fig. 2a) and approached the Pearl River Delta (PRD) region (21.5° - 23.5° N, 112° - 115° E), with an evident enhancement (in terms of composite radar reflectivity) in the following hours (Fig. 2a - d). It maintained quasi-stationary, affecting the PRD region for another 8 h (Fig. 2d and 2f), before it moved toward the shore and dissipated (Fig. 2f).

Figure 2.  Snapshots of the observed composite radar reflectivity (shaded; dBZ) at (a) 0900, (b) 1100, (c) 1300, (d) 1500, (e) 1700, and (f) 2300 BJT on 17 April 2011. The box in each panel indicates the Pearl River Delta region. The red ellipse in (a) shows the focused mesoscale convective system in this study. The red ellipses in (d) show the rain band over coastal Guangdong, and the stratiform precipitation over northern Guangdong.

The southeastward-moving MCS caused severe rainfall along its path, leaving a visible northwest - southeast band of rainfall accumulation (Fig. 3a). The extreme rainfall tended to be concentrated in the PRD region (Fig. 3a and 3b), with the hourly rainfall amount exceeding 100 mm at three in situ surface weather stations (Fig. 3b). The peak hourly rainfall amount was 154.8 mm at station 713524. The temporal evolution of the hourly rainfall amount (averaged within the PRD region) showed two evident peaks (Fig. 3c). The major peak occurred at 1400 BJT and the secondary peak occurred at 2000 BJT. The maximum hourly rainfall at station 713524 reached 77.3 mm, exceeding the 95% percentile of hourly rainfall from 1971 to 2016 over the contiguous Guangdong Province ^[45]. The hourly rainfall peak at the stations with rainfall amounts exceeding 100 mm presented a 1- or 2-h delay compared to the domain-averaged peak depending on their distance to the leading edge of the MCS. By comparing the temporal evolutions of the composite radar reflectivity and the domain-averaged hourly rainfall, we found that the first peak in the hourly rainfall was caused by the linearly developing MCS, which moved towards coastal south China from 1300 to 1700 BJT (Fig. 2c-e), and the second peak was associated with the dissipating MCS, which moved away from land (Fig. 2f).

Figure 3.  (a) Rainfall accumulation from 0900 to 2300 BJT on 17 April 2011 obtained from rain gauges. The rectangle denotes the Pearl River Delta (PRD) region. (b) The same as that in (a) but from 0900 to 2300 BJT on 17 April 2011 over the PRD region. The symbols in the bottom right corner (asterisks, squares and triangles) represent the top three hourly rainfall amounts. (c) Time series of hourly rainfall (units: mm) at three stations, 713524 (asterisk), 713522 (square), and 713528 (triangle), and the PRD averaged hourly rainfall (grey) from 0900 to 2300 BJT on 17 April 2011.

The MCS that contributed to the torrential rainfall was associated with a lower-level mesoscale vortex and the associated shearline. At 0800 BJT on 17 April, a cyclonic circulation was observed from 925 to 850 hPa (labelled by "C" in Fig. 4a and 4d). Below the mesoscale vortex, the surface was characterized by a low sea level pressure center (Fig. 5a). The nearly zonal shearline (denoted by the dashed line in Fig. 4a) extended from the centre of Guangxi Province to the east of Guangdong Province, providing a favorable synoptic-scale lift for the development of the MCS. The PRD region was dominated by the surface low and southwesterlies (Fig. 5a) featuring high equivalent potential temperature (θ_e), PBL warm temperature advection, precipitable water (Fig. 4a and 4d). Such multiscale environmental conditions favored the development of the existing MCS. At 1400 BJT, the shearline oriented nearly parallel to the coastline of Guangdong Province. It directly influenced the PRD region during this period with the supply of moist and warm air (Figs. 4b, 4e, and 5b), contributing to the maintenance of the MCS (Fig. 2d)

Figure 4.  NCEP FNL analysis on 17 April 2011. (a-c) Warm temperature advection (shaded; K s^-1) with horizontal winds at 925 hPa. The dashed black lines indicate the location of shearlines. (d-f) Equivalent potential temperature (contoured at intervals of 5 K with the contour of 340 K in purple and 350 K in red) and precipitable water (shaded; mm) with horizontal winds at 850 hPa. Grey shadings denote the areas with terrain heights above 925/850 hPa. The letter C in (a) and (d) indicates the low-level mesoscale vortex. The rectangle denotes the Pearl River Delta (PRD) region.

Figure 5.  Surface analysis of temperature (shaded; ℃, wind (vector; m s^-1), and sea level pressure (blue contours; hPa) on 17 April. The rectangle denotes the Pearl River Delta (PRD) region.

The above analysis suggests that this torrential rainfall event was closely associated with the thermodynamic and dynamical conditions in the PBL and SL. Therefore, the performances of NWP models in simulating this event are highly associated with the PBL parameterization and its corresponding SL parameterization. To assess how the PBL and SL schemes affect the simulation of the torrential rainfall event, an ensemble consisting of the available suites of PBL and SL schemes in the WRF was performed and presented in section 4.

The performance of the ensemble was first evaluated to examine the predictability of this torrential rainfall case. Good and bad members were further distinguished to investigate the possible influence that different suites of the PBL and SL parameterizations available in WRF model have on simulating the torrential rain over south China. To identify the sources of bias and errors that affect the initiation and development of the convection, we distinguished the good and bad members according to the following criteria. (a) The simulation captured the convective system over coastal area of Guangdong and the stratiform precipitation system over northern Guangdong during the mature stage of the MCS (Fig. 2d). (b) The simulated accumulative rainfall captured the peak rainfall regions at the PRD region as those observed (Fig. 3a).

Figure 6 shows the simulated rainfall accumulation for the 12 WRF ensemble members, as described in Table 1. Compared to the observations, 11 out of the 12 members captured the major spatial patterns of the rainfall accumulation, which was primarily concentrated over the PRD. The observed feature of the northwest-southeastward rain band was also visible in some members, which was indicated by a secondary rainfall centre to the northwest of the PRD (i.e., Figs. 6b, d, f, g, h, i, and l). All of the members whose rainfall accumulations were concentrated over the PRD showed an underestimated rainfall amount (Fig. 6a-k; termed as a good member). The member that overestimated the rainfall amount also largely overestimated the spatial coverage of rainfall around the coastal area of Guangdong Province (i.e., member s10p10, as shown in Fig. 6l; termed as a bad member). Overall, most of the simulated rainfall accumulations were reasonably comparable with the observations, suggesting a relatively high predictability for this torrential rainfall event.

Figure 6.  Rainfall accumulation (shaded; mm) from 0900 to 2300 BJT on 17 April 2011 corresponding to the 12 WRF ensemble members, as shown in Table 1.

Further analysis shows that the evolution of the MCS that was responsible for the torrential rainfall was well captured by the simulations. The 11 good members not only reproduced the MCS (denoted by the ellipse in Fig. 2a) to the northwest of Guangdong Province but also generally well simulated the well-organized MCS parallel to the coastline and the stratiform precipitation system over the north of Guangdong Province at its mature stage (Figs. 7a-k and 2d). Although the bad member s10p10 also captured the orientation of the MCS, it modelled a more intense and compact linear precipitation system. In addition to the simulation failure of the scattered convection, the bad member largely overestimated convection to the southwest of the PRD (Fig. 7l).

Figure 7.  Simulated composite radar reflectivity (shaded; dBZ) at 1500 BJT on 17 April 2011 for the 12 WRF ensemble members, as shown in Table 1.

The relative impacts of PBL and SL schemes on simulating the torrential rain case was examined in this study. This issue was addressed by the comparison of quantitative differences between paired simulations defined by the same PBL or SL scheme. The paired simulations that have the same PBL scheme but different SL schemes were categorized into the"SL-diff"group (e.g., pair s1p5 and s5p5, pair s2p5 and s5p5, as shown in Table 2), while the paired simulations with the same SL scheme but different PBL schemes were categorized into the"PBL-diff"group (e.g., pair s1p1 and s1p5, pair s1p1 and s1p7, as shown in Table 2). There is a total of 5 pairs of simulations in the SL-diff group and 14 pairs in the PBL-diff group. The"pair-RMSE"was calculated to quantitively measure the difference between the simulations in one pair, following the same equation introduced in section 2.2. Generally, the pair-RMSE in the PBL-diff group varied from 9.27 to 14.05 mm, with an average of 11.31 mm, while the pair-RMSE in the SL-diff group varied from 7.57 to 11.09 mm, with an average of 8.69 mm. These results suggest that the simulations with the same PBL schemes are more similar to each other than those with different PBL schemes, indicating a more significant impact of the boundary layer process on torrential rainfall.

Several objective evaluation methods were applied to provide a thorough assessment of the performance of the ensemble. The RMSE, which was used to indicate the deviation between the observations and the simulation, generally showed an approximate 0.5 mm difference for the whole domain for most of the simulations (Fig. 8a). Most of the simulations presented a POD value of ~0.8 at the threshold of 0.1 mm, a ~0.6 POD value at the threshold of 5 mm, a ~0.5 POD value at the threshold of 10 mm, and a ~0.4 POD value at the threshold of 25 mm (Fig. 8b). These results suggested that the ensemble reproduced the actual rain area with a reasonable rainfall amount. The TSs were also comparable among most of the simulations (Fig. 8c). With lower MISS and FAR values (Fig. 8e and 8f), the TS at the threshold of 0.1 mm could reach 0.6. At higher thresholds (i. e., 25 mm), due to the higher MISS and FAR values, the TSs of all simulations were smaller than 0.2 (Fig. 8c, 8e, and 8f). In contrast, the non-bias of the ETS allows the scores to be more fairly compared among different regimes by penalizing both miss and false alarms in the same way. The ETS values were significantly smaller than the TS values due to the stricter definition of the ETS and its double penalization (Fig. 8c and 8d). Most of the ETS values at the threshold of 0.1 mm were significantly smaller than those at a higher threshold (Fig. 8d) because it is easier to correctly forecast the rain occurrence in wet versus dry cases. Thus, random chance comprises a larger proportion of the hits. When the forecasted hits remove random chance, the actual hits would possibly lead to the reduction in the ETS. The ETSs at the threshold of 10 mm and 25 mm can reach 0.2, suggesting a good performance of the simulations. Overall, the simulations show relatively high skill scores, indicating a high predictability of this rainfall event.

Figure 8.  Objective evaluation of rainfall accumulation from 0900 to 2300 BJT on 17 April 2011. (a) Root mean square error (RMSE), (b) probability of detection (POD), (c) threat score (TS), (d) equitable treat score (ETS), (e) missing ratio (MISS), and (f) false alarm ratio (FAR). The grey, red, light blue, blue and green bars indicate the thresholds of total rainfall, as well as 0.1, 5, 10 and 25-mm rainfall, respectively.

By separately comparing the performances of each simulation, the member s10p10 was found to be significantly different from others. Consistent with the analysis in section 4.1.2, the quantitative differences among the 11 good members are comparable in terms of the RMSE, POD, TS, ETS, MISS, and FAR (Fig. 8). The RMSE of member s10p10, which uses the TEMF scheme for both the PBL and SL parameterizations, was evidently larger than that of other members. Such a significant RMSE difference is mainly induced by the overestimation of both spatial coverage and amount of rainfall. With a larger area of intense rainfall compared with the observation, the member s10p10 would always achieve a higher POD and a lower MISS because the POD and MISS are only sensitive to hits when ignoring false alarms (Fig. 8b and 8e). The TS for the member s10p10 at the threshold of 0.1 mm was significantly higher than that for the other members, which shows a similar feature to the POD (Fig. 8b and 8c). The possible reason for this similarity is that the TS only concerns the forecasts that count; that is, the correct negative values have been removed from consideration. The ETS of the member s10p10, on the contrary, was significantly smaller than that of other members at various thresholds (Fig. 8d). Such a smaller ETS is a result of the severe overestimation of rainfall, which is also indicated by its higher FAR (Fig. 8f).

It is well acknowledged that the evaluation result is highly dependent on the evaluation methods. A comprehensive method should provide a more reliable result. In this case, the ETS and RMSE took not only the rainfall coverage into account but also the rainfall amount. The ETS not only equally punishes both the miss and false alarm but also considers random chance. Consequently, the ETS is widely used in evaluating forecast performance. The results based on the ETS and RMSE show that the performance of member s10p10 is significantly different from that of other members given the overestimation of rainfall occurrence and rainfall amount. To investigate the mechanism of the evident simulation bias of member s10p10, some environmental conditions will be further discussed in the next section.

To investigate why the bad member s10p10 (using the TEMF scheme) produced intensive rainfall over a vast spatial coverage, a good member (s1p1) that applied widely used YSU PBL scheme and the revised MM5 SL scheme ^[46] (Table 1) was analysed as a reference. Fig. 9 shows the evolution of the observed and simulated precipitation system in terms of reflectivity. As described in section 3, the observed precipitation system that was responsible for the torrential rainfall arrived at the northwest border of Guangdong Province at 0900 BJT, which was also captured by both members s1p1 and s10p10 (Fig. 9a). At 0930 BJT, intense convection was initiated to the south of the focused precipitation system and rapidly grew both in size and intensity for the member s10p10 in half an hour (Fig. 9b3). This faked convection rapidly developed to become a long squall line in the next 3 h and became the dominant precipitation system in Guangdong Province (Figs. 9d and e3). It almost covered the western Guangdong Province and caused severe rainfall over the study area (Fig. 6l). Meanwhile, the evolution of the simulation using the s1p1 scheme was similar to that via the observations (Fig. 9), except that this simulation was delayed one hour compared to the observation. Overall, the overestimated convection in member s10p10 was a result of the spurious CI, which subsequently and rapidly grew upscale and caused torrential rainfall in Guangdong Province.

Figure 9.  Composite radar reflectivity (shaded; dBZ) at (a) 0900, (b) 1000, (c) 1100, (d) 1200, and (e) 1300 BJT on 17 April 2011 for the (top) observations, (middle) ensemble member s1p1, and (bottom) member s10p10. The black line in (b3) indicates the location of convection initiation and the vertical sections for Fig. 10.

Figure 10 shows the south-north vertical section of equivalent potential temperature (θ_e) across the spurious CI cell in member s10p10. At the initial model time, the θ_e patterns were similar in both the good member s1p1 and bad member s10p10 (Fig. 10a and e). The θ_e and wind fields in member s1p1 evolved subtly in the following hours (Fig. 10a - d), while in member s10p10 at 0900 BJT, southerlies at ~0.8-1.2 km above sea level were visibly enhanced and induced a high - θ_e tongue tilting along the mountain slope. This pattern also converged with the near-surface low - θ_e northerlies, and the high - θ_e southerlies were uplifted to achieve a buoyant updraft and initiate deep moisture convection (Fig. 10g and 10h). Consequently, the spurious CI simulated by the bad member s10p10 was a result of the overestimated high-θ_e elevated air to the south of the CI location.

Figure 10.  Vertical cross sections of equivalent potential temperature (shaded; K) and horizontal winds along the black line in Fig. (8b3) for (a-d) the ensemble member s1p1 and (e-h) ensemble member s10p10 at (a, e) 0800, (b, f) 0830, (c, g) 0900, and (d, h) 0930 BJT on 17 April 2011. The terrain heights are shaded in white at the bottom.

Further analysis was performed to examine the simulated environmental conditions that favoured the spurious CI and the subsequent convective development. Fig. 11 presents the differences in some commonly used thermodynamic variables at the CI time between the good member s1p1 and the bad member s10p10. The thermodynamic environment in member s10p10 was overall moister and more unstable than that of member s1p1. The surface temperature of s10p10 was generally ~ 2 K higher than that of s1p1, especially around the position of the spurious CI and the PRD region (Fig. 11a). The simulated specific humidity in member s10p10 was also found to be ~1.5 g kg^-1 greater than that of member s1p1 near the CI location (Fig. 11b). With warmer and moister PBL conditions in s10p10, the convective available potential energy (CAPE) (or equivalent potential temperature) in s10p10 was significantly greater than that of s1p1, especially around the CI region (Fig. 11c and d). Such a thermodynamically favourable PBL environment was greatly beneficial for the subsequent CI and the rapid upscale growth of triggered convection. It seems that the TEMF scheme tended to produce a warmer and moister PBL environment, which was also suggested by Wang et al. (2014) ^[36], in which the WRF simulations using different PBL schemes were assessed. The vertical profile extracted at the CI location also showed that the major difference in thermodynamic conditions between the members s10p10 and s1p1 existed at low levels (Fig. 12). The ambient temperature below ~700 hPa in the bad member s10p10 was evidently higher than that for the good member s1p1.

Figure 11.  Difference between the ensemble members s1p1 and s10p10 (s10p10 - s1p1) for (a) 2-m temperature (T2; K), (b) 2-m specific humidity (Q2; g kg^-1), (c) surface-based convective available potential energy (CAPE; J kg^-1), and (d) 2-m equivalent potential energy (θ_e; K) at 0930 BJT. The black triangles indicate the locations of spurious convection initiation in the member s10p10.

Figure 12.  Vertical profiles for the ensemble members s1p1 (red) and s10p10 (blue) for (a) temperature (T; K), (b) specific humidity (Qv; kg kg^-1), (c) equivalent potential temperature (θ_e; K), and (d) moist static energy (MSE; J kg^-1) at the locations indicated by the triangles in Fig. 11.

Previous studies have examined the performance of the TEMF scheme from idealized and realistic experiments ^{[30, 47]}. The NWP model using the TEMF parameterization tends to dry the subcloud layer and moisten the lower cloud layer more than the models that purely run for large eddy simulation (LES) (i.e., without PBL parameterization). There was a tendency for the TEMF to move more moisture out of the subcloud layer and into the lower cloud layer, drying the subcloud layer and moistening the lower cloud layer slightly more than the situation in the LES framework ^{[30, 47]}. Such processes lead to a strong moisture contrast in the vertical profile, which may deviate from the realistic profile. On the other hand, the results of the present study also suggest that the TEMF scheme tends to produce greater instability (i. e., higher temperature and larger CAPE) than other schemes for this subtropical case. More importantly, the TEMF scheme used the mass-flux (MF) approach for the parameterization of shallow and deep moist convection ^[48]. The MF closure was able to represent what we refer to as nonlocal transport due to the strong thermals. The main advantage of this method was that it naturally allowed for a scheme for the cloud-topped boundary layer by allowing the moisture in the strong updraft to condense. In such a case, the updraft model is always active and determined independently of whether these updrafts become cloud-core updrafts. Therefore, when the updraft motion was triggered by the terrain, the condensation of moisture released latent heating in the lower to middle levels. The spurious CI in this study was a result of positive feedback among lower-level convergence, upward motion, latent heat release, and surface pressure decrease.

Compared to this new applied scheme of WRF, the traditional PBL scheme, such as the YSU scheme (s1p1 in this study) applied in WRF model since 2004 ^[24] was more widely used. This scheme was an improved vertical diffusion package with a nonlocal turbulent mixing coefficient, compared to what Hong and Pan implemented ^[16], which revealed a consistent improvement in the skill of precipitation forecasts over the continental United States ^[49]. The YSU scheme was proved to produce a more realistic structure of the PBL and its development. Consequently, it did a better job in reproducing the convective inhibition and CAPE, which would reduce the widespread light precipitation and improve some characteristics such as the intensity of convection ^[24]. During the annually first rainy season, the prevailing low-atmosphere monsoonal flows with moist and warm air provided favorable thermodynamic conditions for the CI in Guangdong Province. Additionally, the interaction of onshore monsoons and the mountainous orography over this region was conducive to the mesoscale lift. In such environmental conditions, the improved representation of convective inhibition and CAPE in the YSU scheme would be expectant to successfully prevent the faked CI produced in TEMF experiment. Therefore, the YSU scheme may be more suitable for simulating the torrential rain over South China during the monsoon season.

This study examined the impacts of the planetary boundary layer (PBL) and surface layer (SL) schemes on a simulated torrential rainfall event over south China by the ensemble-based analysis. The ensemble simulations consist of the available suites of the PBL and SL schemes in the WRF model. Twelve ensemble members were assessed against the simulated rainfall pattern, MCS evolution and forecast skill scores (i. e., RMSE, POD, TS, ETS, MISS, and FAR). Further comparison between the bad and good members revealed the possible mechanism of the simulation failure for the bad member.

The predictability of the torrential rainfall case was high based on the spatial patterns of the simulated rainfall accumulation and the evolution of the MCS that contributed to rainfall. Most of the simulations captured the observed maximum rainfall centre over the PRD. Moreover, the objective evaluation also indicated the high predictability, as the ETS at the threshold of 25 mm reached as much as 0.2. It is worth noting that most of the simulations underestimated the rainfall accumulation. To better illustrate the impact of the model resolution in simulating this case, a similar ensemble run with the same configuration but a resolution of 3 km in the inner domain was performed. However, the results showed that the higher-resolution simulations presented some noisy, scattered convection in almost all of the simulations that were not observed (not shown).

The impact of the objective evaluation methods on assessing the quality of the simulation was also discussed in this study. The metric for a better evaluation method should include all properties in the contingency table (hits, misses, false alarms, and correct negatives). The evaluation methods, such as probability of detection (POD), missing ratio (MISS), false alarm ratio (FAR), and threat score (TS), consisted of only two or three properties, which might cause bias in the verification statistics. On the other hand, the non-bias in the equitable threat score (ETS) allowed scores to be compared more fairly across different regimes. Therefore, it is more reliable to use the ETS in the verification of rainfall in NWP models.

In the ensemble simulations, the member using the total energy mass flux (TEMF) scheme was significantly different from the others, which was characterized by more extensive rainfall coverage and a more intensive MCS. The comparison among the TEMF simulation, the reference simulation, and observations showed that the favourable thermodynamic supply and terrain lift triggered a spurious convection initiation that was not observed or simulated by other members. This spurious convection then developed and matured under favourable conditions with high θ_e air masses. According to the discussion on the theory of the TEMF scheme, the spurious CI was due to the mechanism by which the TEMF scheme tended to produce a warmer and moister PBL environment. These results suggested that accurately expressing the thermodynamic environment in the PBL might avoid overestimating the simulated rainfall.

The torrential rainfall events over south China are often associated with the boundary layer and surface layer process. Therefore, understanding the predictability and mechanism of the PBL and SL schemes in simulating torrential rain over south China is not only of great importance in operational forecast but also a useful guide for localizing the PBL and SL parameterizations. Additional research is needed to reveal how to improve the PBL and SL schemes in different regimes. Although the simulations and analyses presented in this study are based on one case, they provide a possible error source for the PBL and SL schemes in simulating torrential rainfall over south China. Future studies will conduct long-term evaluations of the various parametrization schemes and their coupling in the WRF model to identify the sources and biases in current PBL and SL schemes, which would be helpful for improving the PBL and SL schemes to forecast rainfall over south China via the operational forecast system.