Uncertainty in Tropical Cyclone Intensity Predictions due to Uncertainty in Initial Conditions

Over the past few decades, tropical cyclone (TC) track forecasts have been improved substantially (Duan et al., 2012; Qian et al., 2012; Cangialosi and Franklin, 2016), whereas TC intensity forecasts have obtained little improvement (DeMaria et al., 2014). The discrepancy between TC track and intensity forecast accuracy may be due to the fact that TC tracks are largely dependent on the large-scale environment, which is increasingly better resolved by global models, whereas intensity is largely determined by internal dynamics and smaller-scale processes that are often poorly resolved in global models and some regional models. TC intensity forecasts are typically compromised by considerable uncertainty, with errors in initial conditions as one of the main sources of forecast error.

Using a perfect model framework, Emanuel and Zhang (2016) attempted to quantify the intrinsic predictability of TC intensity and showed that TC intensity forecast error is dominated by initial condition errors out to at least a few days. They further showed that the growth of intensity error is at least as sensitive to the specification of inner-core moisture as it is to that of the wind field (Emanuel and Zhang, 2017). Further studies of idealized TCs have shown that when a small random moisture perturbation is added to the initial boundary layer, the result is a non-negligible spread in TC intensity (e.g., van Sang et al., 2008). There is an intrinsic uncertainty in the prediction of TC maximum intensity. In the presence of environmental shear, initial random noise may result in changes in the timing of TC rapid intensification. The larger the vertical wind shear, the larger the uncertainty in TC intensity forecasts (Zhang and Tao, 2013).

Studies of actual TCs show similar results to those from simulations of idealized TCs. A probabilistic analysis of a non-developing tropical disturbance showed that deep moisture and high CAPE in the initial conditions are the two most important factors for TC genesis (Sippel and Zhang, 2008). Another study of the same tropical disturbance showed that as a result of the effects of initial condition errors on moist convection, small and arguably unobservable differences in initial conditions not only result in different routes to TC genesis, but also that these differences can determine whether or not a tropical disturbance will develop significantly (Zhang and Sippel, 2009). A study of Hurricane Karl (2010) also showed that initial condition errors at the meso scale pose large challenges in the prediction of TC genesis (Poterjoy and Zhang, 2014).

A study of the predictability of Hurricane Humberto (2007) showed that mid-level moisture, low-level convective instability, and the strength of a front are the key factors that led to spread in forecasts of the storm’s strength (Sippel and Zhang, 2010). In model simulations of the non-intensifying Tropical Storm Erika (2009), mid-level relative humidity, absolute vorticity, and the distribution of convection all played a role in determining whether or not the storm intensified (Munsell et al., 2013). In simulations of Hurricane Nadine (2012), an initial intensity bias was most likely the cause of final intensity errors (Munsell et al., 2015). Sensitivity experiments of Hurricane Edouard (2014) revealed that some of the variation in the timing of the hurricane’s rapid intensification could be attributed to differences in its initial intensity (Munsell et al., 2017).

Although many studies of the uncertainty and predictability of TC intensity have been conducted, these studies have mainly focused on idealized TCs and Atlantic hurricanes. With regards to their physical origin, TCs in the Northwest Pacific and South China Sea differ from Atlantic TCs. In this work, we use Weather Research and Forecasting (WRF) initial perturbation ensemble forecasts to investigate the intensity forecast uncertainty of the No. 9 TC Rammasun in 2014, which is the strongest TC to have made landfall in China during the past 50 years. Section 2 describes the forecast model, the data used as the ensemble initial conditions, and the data used to verify the forecast results. In section 3, we present an overview of TC Rammasun, and in section 4 we present results from the ensemble forecast. A discussion and summary are given in section 5.

We use the WRF model (v3.4) configured with three static, two-way nested domains. The horizontal grid spacing of the outer, middle, and inner domain is 27, 9, and 3 km, respectively, and the inner domain covers the region (7.3°−26.4°N, 104.2°−129°E). All forecast results shown are for the inner domain. Vertically, the model employs the terrain-following stationary pressure coordinate and there are 28 levels, with the pressure at the top level being 50 hPa.

The model’s three domains use consistent physical process schemes, except for the convection parameterization scheme. It is commonly thought that when the horizontal resolution reaches 4−5 km, models can resolve cumulus convection rather than relying on its parameterization (Tripoli, 1992; Gopalakrishnan et al., 2013). Therefore, no convection parameterization scheme is applied in the inner domain, whereas the Kain−Fritsch convection parameterization scheme is applied in the outer and middle domains. The other physical process schemes used include the Yonsei University boundary parameterization scheme, the WRF single-moment 6-class microphysics scheme, the Dudhia shortwave radiation scheme, the RRTM longwave radiation scheme, and the Noah land surface scheme.

We use a 21-member ensemble including a control forecast member and 20 initial perturbation members. The initial conditions for the 21 members are obtained from the initial fields of the National Centers for Environmental Predication (NCEP) Global Ensemble Forecast System (GEFS) with a resolution of T190L28. This system consists of 21 members, and the method of ensemble transform with rescaling is used to generate the initial perturbations of the 20 perturbation members. The initial fields of NCEP GEFS are initialized to form the initial fields for the 21 ensemble members used in this study. The lateral boundary conditions for the 21 members are also obtained from the initial fields of NCEP GEFS at several different times. The initialized data at these times provide the lateral boundary conditions for the model’s outer domain, and the lateral boundary conditions are updated every 6 h. Corresponding to the horizontal resolution of the three domains, the step length of integration is 90, 30, and 10 s in the outer, middle, and inner domain, respectively, and the integration time is 48 h. The observational data used to verify the forecast results are from the TC best-track dataset compiled by the Shanghai Typhoon Institute of the China Meteorological Administration (Ying et al., 2014).

TC Rammasun formed over the ocean to the east of the Philippines at 0600 UTC 12 July 2014, reaching the South China Sea on 16 July. It reached its strongest intensity at 0600 UTC 18 July, when at its center the maximum wind velocity was 72 m s⁻¹ and the lowest pressure was 888 hPa. Rammasun subsequently made landfall on the coast of China three times. Based on the lowest recorded value of central pressure—which occurred when Rammasun made landfall in Wenchang, Hainan—Rammasun became the strongest TC to have made landfall in China during the past 50 years. Detailed information about TC Rammasun can be found in Wang and Zeng (2018a) and Wang and Zeng (2018b).

We focus on a forecast period from 0000 UTC 16 July to 0000 UTC 18 July, a total of 48 h, when Rammasun moved away from the Philippines over the South China Sea towards China, before making landfall in Wenchang, Hainan. Rammasun underwent a process of intensification from typhoon to strong typhoon and then to super-strong typhoon.

Errors in initial conditions may result in uncertainty in TC forecasts. The root mean of difference total energy (RM-DTE) is used to characterize forecast uncertainty (Sippel and Zhang, 2008). The larger the RM-DTE, the greater the forecast uncertainty. RM-DTE is calculated as

DTE in this equation is a commonly used measure of ensemble spread (Mitchell et al., 2002; Zhang et al., 2006; Meng and Zhang, 2007) and is calculated as

where u', v' and T ' denote the difference between a member and the ensemble mean of the zonal velocity component, the meridional velocity component and the temperature, respectively, and $k = {C_p}/{T_{\rm{r}}}$, where C_p = 1004.9 J kg⁻¹ K⁻¹ and T_r = 270 K. In the formula for RM-DTE, k_max is the maximum vertical level of the model domain, and N_e is the number of ensemble members.

In Fig. 1, we show the spatial distribution of RM-DTE in the analysis field (0 h), and in the 12-h, 24-h, 36-h and 48-h forecasts. Maximum RM-DTE increases from ~0.9 m s⁻¹ at the analysis time to over 13 m s⁻¹ at 24 h and to over 30 m s⁻¹ at 48 h. In the analysis, the maximum RM-DTE is not in the vicinity of the TC center, whereas at other times it is located to the right of the TC center. The evolution of RM-DTE shows that the uncertainty in initial conditions leads to uncertainty in the forecasts of the TC, and that this uncertainty increases with forecast time, with a significant and rapid increase from 24 h to 48 h.

Figure 1.  Spatial distribution of RM-DTE (units: m s⁻¹) in (a) the analysis (0 h), and in the (b) 12-h, (c) 24-h, (d) 36-h, and (e) 48-h forecast. The typhoon symbol indicates Rammasun’s position

Figure 2 shows the observed and forecasted intensity and track of Rammasun over the 48-h forecast period. As expected, initial condition errors have a large impact on predicted TC intensity but little effect on predicted TC track. In other words, uncertainty in the TC forecasts due to errors in initial conditions is mainly reflected in TC intensity forecast uncertainty. From the predicted evolution of minimum sea level pressure (SLP) in each member (Fig. 2a), it can be seen that the intensity difference among members increases with forecast time, with rapid increase after 24 h. At 48 h, member 11 has the strongest TC intensity and member 2 the weakest, with the difference in predicted intensity reaching 37 hPa. Comparing the observed and predicted track over the 48-h forecast period (Fig. 2b), all members predict a track close to the observed track; at 48 h, the predicted position of the TC in almost all members is within 100 km of the observed location.

Figure 2.  Observed (black line) and ensemble mean forecast (yellow line) of Rammasun’s (a) intensity (minimum SLP; units: hPa) and (b) track over the 48-h forecast period. Also shown are results for individual members (gray lines), with certain members highlighted (member 2, green line; member 4, blue line; and member 11, red line)

The two most extreme ensemble members are members 2 and 11, with the weakest and strongest TC intensity at 48 h, respectively. The intensity results for member 4 lie somewhere between these two extremes (Fig. 2a). The following analyses focus on these three members.

Figure 3 shows the predicted 10-m horizontal wind speed (HWS10) and SLP fields in the three members at 12 h, 24 h, 36 h, and 48 h. All members predict a TC that develops slowly before 24 h. However, the TC in member 4 actually weakens between 12 h and 24 h, ensuring that the range of intensity among the members is small. All three members predict a TC that develops rapidly after 24 h, though development is stronger in member 11 compared with members 4 and 2. This increases the diversity in the intensity predictions among the members; at 48 h, the intensity difference between members 11 and 2 reaches 14 m s⁻¹ and 37 hPa, based on HWS10 and minimum SLP, respectively.

Figure 3.  12-h (row 1), 24-h (row 2), 36-h (row 3), and 48-h (row 4) HWS10 (color shading; units: m s⁻¹) and SLP (contours; units: hPa) forecasts from member 11 (column 1), member 4 (column 2), and member 2 (column 3). The numbers in the lower left corner of each panel are the maximum HWS10 and minimum SLP in the corresponding member at the corresponding forecast time

This variation in TC intensity can be more clearly seen in the evolution of maximum HWS10 in the three members (Fig. 4). Before 24 h, the TC in member 11 is not the strongest; however, after 24 h, the TC in member 11 strengthens so much and so rapidly that its intensity after 30 h becomes the strongest among ensemble members. Member 11 predicts a TC intensity closest to the observations (Figs. 2a and 4).

Figure 4.  The evolution of observed and predicted maximum HWS10 (units: m s⁻¹), together with area-averaged HWS10 (units: m s⁻¹) from the three members. A circle of 100 km radius centered on Rammasun is used to obtain area-averaged quantities

We next investigate uncertainty in the initial conditions, which is necessary for a comprehensive study of TC intensity forecast uncertainty. The overall uncertainty in the initial conditions is characterized by the ensemble spread (ES) and the root-mean-square difference (RMSD). The ES and RMSD are calculated as:

where b' denotes the variable difference between a member and the ensemble mean and b^* denotes the variable difference between two members. N_e is the number of ensemble members and n_g is the number of grid points of the model inner domain.

In Fig. 5, we show the vertical distribution of the initial ES in zonal wind, temperature, and relative humidity, together with the vertical distributions of the RMSD in these variables between members 11 and 4 and between members 11 and 2 (Fig. 5). The RMSDs are larger than the spread for each of the three variables. The RMSD between members 11 and 2 is consistently larger than that between members 11 and 4, especially for temperature. For temperature below the 700-hPa level, the difference between members 11 and 2 is significantly larger than that between members 11 and 4; a similar pattern is found for zonal wind below 850 hPa. For temperature, zonal wind and relative humidity, the largest discrepancy between the two RMSDs occurs at 900 hPa; for temperature and zonal wind at this height, the difference between members 11 and 2 is approximately twice that between members 11 and 4.

Figure 5.  Vertical distribution of the initial ensemble spread (blue) and RMSDs between members 11 and 4 (orange) and between members 11 and 2 (green) for (a) zonal wind (units: m s⁻¹), (b) temperature (units: K) and (c) relative humidity (units: %)

Figure 5 shows that the difference between members 11 and 2 is 0.62−2.22 m s⁻¹ for zonal wind, 0.21−0.66 K for temperature, and 2.21%−12.62% for relative humidity. These amounts are smaller than the corresponding NCEP-assumed errors for rawinsondes, which is about 1.4−3.1 m s⁻¹ for zonal wind, 0.8−1.2 K for temperature, and 5%−20% for relative humidity (Sippel and Zhang, 2008). The analysis error consists of the observational error and the error caused by the analysis method. The analysis error is generally close in magnitude to the observational error. Therefore, although the RMSD between members 11 and 2 is clearly larger than the spread, it is smaller than the analysis error.

Figure 6 shows initial difference fields of mixing ratio, temperature and horizontal wind vectors between members 11 and 2 and between members 11 and 4, at four different heights. For these three variables, the differences between members 11 and 2 and between members 11 and 4 are largest at 900 hPa and smallest at 700 hPa. The difference between members 11 and 2 is larger than the difference between members 11 and 4 at all heights except 700 hPa. The RMSDs in relative humidity are larger at 700 hPa than at the surface (Fig. 5c), whereas the differences in mixing ratio are larger at the surface than at 700 hPa. Furthermore, considering the difference between members 11 and 2, at lower heights (especially at the surface) the area through which the TC will pass shows a large positive difference in humidity. The initial difference fields of equivalent potential temperature (θ_e) between members 11 and 2 and between members 11 and 4 at three heights (Fig. 7) have similar characteristics to the humidity results. It is evident that humidity at lower heights contributed to the TC’s development.

Figure 6.  Initial difference fields of 700-hPa (row 1), 850-hPa (row 2), 900-hPa (row 3) and surface (row 4) mixing ratio (color shading; units: g kg⁻¹), temperature (contours; units: K) and horizontal wind vectors (units: m s⁻¹) between members 11 and 2 (column 1) and between members 11 and 4 (column 2). In each panel, the five typhoon symbols from the lower right corner to the upper left corner indicate the position of Rammasun at 0 h, 12 h, 24 h, 36 h and 48 h, respectively

Figure 7.  Initial difference fields of 700-hPa (column 1), 850-hPa (column 2) and 950-hPa (column 3) θ_e (units: K) between members 11 and 2 (row 1) and between members 11 and 4 (row 2). In each panel, the five typhoon symbols from the lower right corner to the upper left corner indicate the position of Rammasun at 0 h, 12 h, 24 h, 36 h and 48 h, respectively

Equivalent potential temperature θ_e is used as an indicator of convective instability. In Fig. 8 we show the forecasted spatial distribution of 900-hPa θ_e, 900-hPa horizontal wind vectors, and 700-hPa vertical wind speed for members 11, 4, and 2, at four stages during the forecast. Convective instability increases with forecast time in all members, increasing more rapidly and significantly after 24 h; this is especially the case in member 11. This means that the disparity among members is larger after 24 h than it is before. Member 11 has the largest convective instability, so it has the strongest convection development after 24 h. The variation in convective instability is similar to the variation in intensity (Fig. 3). Analysis of the evolution of area-averaged 900-hPa θ_e (Fig. 9a) reveals variation in convective instability over the entire region. After ~18 h, convective instability over the entire region is largest in member 11 and after ~24 h is smallest in member 2.

Figure 8.  Spatial distribution of 900-hPa θ_e (color shading; units: K), 900-hPa horizontal wind vectors (units: m s⁻¹) and 700-hPa vertical wind speed (contours; units: m s⁻¹) at a forecast time of 9 h (row 1), 21 h (row 2), 33 h (row 3), and 45 h (row 4) for member 11 (column 1), member 4 (column 2), and member 2 (column 3)

Figure 9.  Evolution of area-averaged (a) 900-hPa θ_e (units: K) and (b) 700-hPa vertical wind speed (units: m s⁻¹) for members 11 (purple), 4 (green), and 2 (blue). A circle of 100 km radius centered on Rammasun is used to obtain area-averaged quantities

Figure 8 also shows that the ascending motion strengthens as convective instability increases. The larger the convective instability, the stronger the ascending motion. The strong ascending motions are distributed around the region of high θ_e. Before 24 h, strengthening of ascending motions is not apparent, with the maximum upward velocity even decreasing. Based on the evolution of area-averaged 700-hPa vertical wind speed for the three members (Fig. 9b), it is clear that vertical velocity decreases before ~18 h but later intensifies rapidly after ~33 h. Ascending motions strengthen significantly from 33 h to 45 h, especially in member 11 (Fig. 8).

We next show forecasted 3-h accumulated precipitation results from the three members (Fig. 10). The strength of precipitation reflects the strength of convective activity. The precipitation forecasts diverge as precipitation increases; after 24 h, accumulated precipitation is much larger in member 11 compared with members 4 and 2, with member 2 forecasting the lowest accumulated precipitation. This variation in forecasted precipitation is consistent with the variation in convective instability (Fig. 8). An increase in convective instability leads to more convection and reinforcement of ascending movement, which generates more precipitation.

Figure 10.  Forecasted accumulated precipitation (units: mm) at 9−12 h (row 1), 21−24 h (row 2), 33−36 h (row 3), and 45−48 h (row 4) in member 11 (column 1), member 4 (column 2), and member 2 (column 3)

TC development requires a continuous supply of water vapor and heat. In Fig. 11, we show the spatial distributions of the surface upward heat flux, water vapor flux and 10-m horizontal wind vectors forecasted by the three members at four forecast times. At 12 h, the heat flux and water vapor flux predicted by member 11 are both less than those in members 4 and 2. However, the fluxes subsequently increase rapidly and continuously in member 11 and after 24 h they exceed those in the other two members. The heat flux and water vapor flux reduce slightly from 12 h to 24 h in member 4, as is the case in member 2 between 24 h and 36 h. However, in members 4 and 2 both quantities increase significantly after 36 h, with larger increases in member 4 than in member 2. The variation in heat flux and water vapor flux is roughly consistent with the variations in convection (Fig. 8 and Fig. 10) and TC intensity (Fig. 3).

Figure 11.  Forecasted evolution of surface upward heat flux (color shading; units: W m⁻²), water vapor flux (contours; units: g m⁻² s⁻¹) and 10-m horizontal wind vectors (units: m s⁻¹) at a forecast time of 12 h (row 1), 24 h (row 2), 36 h (row 3), and 48 h (row 4) for member 11 (column 1), member 4 (column 2), and member 2 (column 3). The numbers in the upper left corner of each panel are the corresponding maximum heat flux and water vapor flux, respectively

The forecasted evolution of area-averaged surface upward heat flux and water vapor flux in the three members is shown in Fig. 12. The variation of heat flux is similar to that of water vapor flux in all three members. The evolution of the whole-region heat flux and water vapor flux is similar to that of the whole-region convective instability (Fig. 9a) and TC intensity (Fig. 4).

Figure 12.  Evolution of three members’ area-averaged surface upward (a) heat flux (units: W m⁻²) and (b) water vapor flux (units: g m⁻² s⁻¹) in member 11 (purple), member 4 (green), and member 2 (blue). A circle of 100 km radius centered on Rammasun is used to obtain area-averaged quantities

In this section, we investigate the correlation between TC intensity and factors affecting TC development, to explore the contribution of different factors to the uncertainty in TC intensity forecasts. The correlation coefficient r between two variables a any b is calculated as

where a_en and b_en are the ensemble means of the variables a any b, and N_e is the number of ensemble members.

Figure 13 shows the evolution of correlation coefficients between different area-averaged variables. TC intensity before 24 h has little correlation with the final TC intensity (Fig. 13a); this is particularly the case when HWS10 is used as a measure of TC intensity. After 24 h, correlation between the final TC intensity and instantaneous intensity through the forecast increases rapidly with forecast time, resulting in moderate to strong correlation. Among the forecasted variables related to TC development, the strongest correlation with TC intensity is found for the surface upward heat flux and water vapor flux. Correlation between instantaneous TC intensity and heat flux and between TC intensity and water vapor flux remains strong throughout the forecast (Fig. 13b). The correlation coefficients between TC final intensity and heat flux or water vapor flux over the course of the forecast are lower than 0.4 before 30 h, indicating weak correlation. These values increase to above 0.6 after 30 h, indicating that the final intensity of the TC is determined in part by the strength of heat and water vapor fluxes. This result is consistent with the results shown in Fig. 13a.

Figure 13.  Temporal evolution of the correlation coefficient between (a) TC intensity at 48 h and TC intensity at each forecast time (red: HWS10; green: SLP) and (b) TC intensity (HWS10) at 48 h and surface upward heat flux (purple) and water vapor flux (blue) at each forecast time, and between TC intensity (HWS10) at each forecast time and the corresponding surface upward heat flux (red) and water vapor flux (green) at that time. A circle of 100 km radius centered on Rammasun is used to obtain area-averaged quantities

We have used a WRF initial perturbation ensemble forecast of the super-strong typhoon Rammasun (1409) to investigate uncertainty in forecasts of TC intensity and the role of initial condition uncertainty. The ensemble forecasts comprised a control forecast and 20 initial perturbation members.

Initial condition uncertainty led to TC forecast uncertainty, particularly in forecasts of Rammasun’s intensity. This uncertainty increased with forecast time with a more rapid and significant increase after 24 h. At 48 h, the maximum intensity difference between different ensemble members reached 14 m s⁻¹ (HWS10) and 37 hPa (SLP). Before 24 h, Rammasun developed slowly and the member with the strongest final intensity was not at that stage the member with the strongest TC. However, after 24 h the TC in this member strengthened much more than other members, and its intensity became closest to the observed value.

Variation in convective instability, precipitation, surface upward heat flux and water vapor flux was similar to the variation in TC intensity. These variables increased with forecast time, with rapid and significant increase after 24 h. At 48 h, these quantities were larger in the strongest member than in other members. Over the course of the forecast, there was strong correlation between instantaneous TC intensity and the surface upward heat flux and water vapor flux at that instance. After 30 h, there was strong correlation between TC final intensity and the instantaneous surface upward heat flux and water vapor flux.

Divergence in the TC intensity forecasts stems from initial condition uncertainty. Corresponding to the maximum intensity difference at 48 h, differences in the initial conditions are much larger than the initial ES and larger than other initial condition differences. Even so, this difference is smaller than the analysis error on the whole, whether for zonal wind, temperature or relative humidity. Initial condition errors smaller than the analysis errors can result in large error in the TC intensity forecasts after 48 h.

Initial humidity differences, and to a lesser extent temperature differences, at the surface and at lower heights are the key factors that lead to differences in the forecasts of TC intensity among ensemble members. These differences in initial humidity and temperature relate to both the overall values and distribution of these parameters. Compared with the member with the weakest final TC intensity, the member with the strongest final TC intensity had greater initial humidity at the surface and at lower heights in regions through which the TC subsequently passed, as well as higher initial surface air temperatures in these regions. The enhanced humidity and temperature in this member drives greater convective instability. Along with stronger ascending motions, the greater convective instability led to enhanced convective activity, which generated more precipitation. With stronger convective activity, more heat and water vapor were transported upward from the surface. All of these factors made the TC increasingly strong, ensuring that by 48 h the TC in this member was the strongest among all members.

Initial condition uncertainty results in uncertainty in TC intensity forecasts. The TC intensity forecast error is mainly dominated by initial humidity errors and to a lesser extent temperature errors. This result is consistent with the results of other studies (Sippel and Zhang, 2008; van Sang et al., 2008; Zhang and Sippel, 2009; Emanuel and Zhang, 2017), which all emphasized the role of humidity in TC genesis and development. In this work, we have focused on a singular TC and therefore the results are only applicable to Rammasun. To further understand uncertainty in TC intensity forecasts, future studies should focus on a range of TCs with different characteristics.

Acknowledgements. This study received financial support from the National Natural Science Foundation of China (Grant Nos. 41575108 and 41475082).