Abstract

In this paper, the impact of outer loop and partial cycling with the Weather Research and Forecasting Model’s (WRF) three-dimensional variational data assimilation system (3DVAR) is evaluated by analyzing 78 forecasts for three typhoons during 2008 for which Taiwan’s Central Weather Bureau (CWB) issued typhoon warnings, including Sinlaku, Hagupit, and Jangmi. The use of both the outer loop and the partial cycling approaches in WRF 3DVAR are found to reduce typhoon track forecast errors by more than 30%, averaged over a 72-h period. The improvement due to the outer loop approach, which can be more than 42%, was particularly significant in the early phase of the forecast. The use of the outer loop allows more observations to be assimilated and produces more accurate analyses. The assimilation of additional nonlinear GPS radio occultation (RO) observations over the western North Pacific Ocean, where traditional observational data are lacking, is particularly useful. With the lack of observations over the tropical and subtropical oceans, the error in the first-guess field (which is based on a 6-h forecast of the previous cycle) will continue to grow in a full-cycling limited-area data assimilation system. Even though the use of partial cycling only shows a slight improvement in typhoon track forecast after 12 h, it has the benefit of suppressing the growth of the systematic model error. A typhoon prediction model using the Advanced Research core of the WRF (WRF-ARW) and the WRF 3DVAR system with outer loop and partial cycling substantially improves the typhoon track forecast. This system, known as Typhoon WRF (TWRF), has been in use by CWB since 2010 for operational typhoon predictions.

1. Introduction

Typhoons are among the most significant weather systems in Taiwan (Wu and Kuo 1999). Taiwan is an island that is roughly 200 km wide and 400 km long centered at 23.5°N, 121°E. On average, three to four typhoons make landfall there each year. Although the heavy rainfall and gusty winds from typhoons routinely produce significant damage to the island, the rainfall is also one of the most important water resources in Taiwan. Therefore, providing accurate typhoon forecasts is the highest priority of Taiwan’s Central Weather Bureau (CWB). Operational typhoon forecasts are primarily based on the guidance of numerical weather prediction models. Typhoon prediction in the vicinity of Taiwan faces two major challenges: (i) the lack of traditional observational data over the western North Pacific and (ii) the influence of Taiwan’s Central Mountain Range (CMR) on the typhoon circulation patterns. Because of the influence of the CMR, the distribution of winds and rainfall for a typhoon approaching Taiwan is highly sensitive to the track and direction of the approaching storm.

The accuracy of winds, rainfall, ensuing storm surge, and severe weather forecasts for typhoons affecting Taiwan critically depend on the performance of track forecasts. Research on tropical cyclone track and intensity forecasts has suggested this to be a multiscale problem, from the synoptic scale to the convective scale (Holland and Merrill 1984; Molinari and Vollaro 1990; Liu et al. 1997; Ross and Kurihara 1995; Marks and Shay 1998). Therefore, having an accurate initial state for a regional model with realistic representation of multiple-scale processes is a key element in the system’s success. The performance of a data assimilation system that combines observations with short-term model forecasts, analysis update cycles with a well-developed spun-up vortex, and/or a proper method for vortex implantation can significantly impact the skill of a tropical cyclone prediction model (Kurihara et al. 1995; Hsiao et al. 2009).

The data assimilation update cycle is an important component of a typhoon prediction model. In addition to improving the model’s initial conditions, a cycling data assimilation can mitigate model imbalance resulting from interpolation of global model forecasts. Moreover, the data assimilation update-cycling procedure allows higher-resolution forecasts to be used as the background (Hodur 1987; Xiao et al. 2005; Zhao et al. 2006). A data assimilation system relies on specified error statistics to obtain an analysis with an optimal combination of observations and background (which is derived from a short-range forecast).

Over the tropical and subtropical oceans, there are a limited number of traditional observations. Under such conditions, the errors associated with the background field can continue to grow, as there are insufficient observations to constrain the model error growth. As a result, a data assimilation system with continuous (full) cycling can produce inferior analyses. Recently, the National Centers for Environmental Prediction (NCEP) adopted a partial cycling strategy for the North American Mesoscale (NAM) model’s Data Assimilation System (NDAS) with a goal of improving its forecast skills. NDAS partial cycling consists of a cold start from 12 h prior to the analysis time and four subsequent data assimilation cycles at 3-h intervals. Rogers et al. (2009) compared the forecast differences between full and partial cycling from 17 December 2007 to 25 February 2008. The results show a significant improvement from partial cycling for height and wind fields at 48- and 72-h forecasts particularly in the 250–300-hPa layer. Moreover, the largest improvement in temperature occurs in the middle troposphere between 700 and 400 hPa. The partial cycling reduces error growth with time, which results from the accumulation of forecast error in the full cycling run. As a result, partial cycling has been implemented in NDAS operations since December 2008.

In variational data assimilation, the observation operators are linearized in the form of a first-order Taylor series expansion and then the inner loop procedure solves a minimization algorithm for a quadratic (linearized) problem (Courtier et al. 1998). The outer loop procedure solves the nonlinear aspects of the assimilation problem with the first guesses updated for each minimization process (inner loop). Because of the outer loop procedure and the Weather Research and Forecasting Model (WRF) three-dimensional variational data assimilation system’s (3DVAR) incremental formulation, the approximate solution, which is derived from a minimization algorithm, can be improved by solving Taylor series expansion recursively. With the outer loop procedure, the analysis should be as accurate as possible through the linearized Taylor series terms in the calculation of the observation innovation; meanwhile, more observations would be used during the subsequent outer loops. For observations such as global positioning system (GPS) radio occultation (RO) refractivity, which are nonlinearly related to the analysis control variables, more information can be extracted from outer loops (Rizvi et al. 2008). In WRF 3DVAR, the background check quality control (BgQC) procedure, as described by Andersson and Järvinen (1999), was implemented. The background, or first guess, is obtained from the analysis of the previous outer loop. Because the first guess is closer to the observations, some observations that were rejected by the previous outer loop can now pass the BgQC in the current outer loop and influence the analysis. Consequently, the number of observations assimilated is increased with more outer loops. Guo et al. (2006) applied the outer loop strategy for the analysis of Typhoon Haitang (2005) and obtained a positive impact on track forecasts. In addition, Massart et al. (2010) demonstrated the benefits of the outer loop in the 3DVAR approach with the first guess at appropriate time (3D-FGAT) assimilation scheme for atmospheric transport. When the assimilation operators are nonlinear, the outer loop is able to update the linearization around a better state and thus improve the quality of the analysis. Rosmond and Xu (2006) showed that the outer loop strategy produced better assimilation of surface wind speeds and precipitable water in the Naval Research Laboratory Atmospheric Variational Data Assimilation System (NAVDAS). Gridpoint statistical interpolation (GSI) in the National Centers for Environmental Prediction’s (NCEP) Global Data Assimilation System (GDAS), which has been in operation since 2007 (Kleist et al. 2009), included two outer loops. With the lack of observations over the tropical and subtropical oceans, an effective assimilation of limited observations is crucial for typhoon prediction.

More recently, the WRF has increasingly been applied to tropical cyclone forecasts by the science community (Davis et al. 2008; Hsiao et al. 2010; Lin et al. 2010; Torn 2010). Davis et al. (2008) showed that the use of a mixed layer ocean model and improved surface momentum exchange in the WRF improved the intensity forecast for Hurricane Katrina (2005). Hsiao et al. (2010) showed that a vortex relocation scheme in the WRF improved the Advanced Research core of the WRF’s (WRF-ARW) typhoon initialization and forecasts, particularly when it was coupled with data assimilation update cycling. Meanwhile, Lin et al. (2010) simulated Hurricane Isabel (2003) on extreme rainfall, winds, and surge using the WRF coupled with the 2D, depth-averaged hydrodynamic Advanced Circulation Model (ADCIRC). A WRF-ARW based ensemble Kalman filter (EnKF) system was shown to give superior performance in track and intensity prediction for 10 Atlantic basin tropical cyclones (Torn 2010). An important advantage of EnKF data assimilation is its ability to use flow-dependent error statistics estimated from short-term ensemble forecasts and does not require a bogus vortex (Chen and Snyder 2007; Torn and Hakim 2009). These results suggest that the ARW can be used as a tropical cyclone prediction model with respectable skills.

In this study, we assess the impact of the partial cycling and outer loop approaches in the WRF 3DVAR system on tropical cyclone prediction. The model description and experiment design are presented in section 2. Results of forecasts for Typhoons Sinlaku, Hagupit, and Jangmi (2008), for which CWB issued typhoon warnings, are discussed in section 3. Finally, conclusions are given in section 4.

2. WRF description and the experiment design

The WRF-ARW uses a terrain-following, hydrostatic-pressure vertical coordinate with the model top being a constant pressure surface of 30 hPa in this study. The CWB operational configuration used here is based on version 3.1.1 with a 221 × 127 grid and 45-km grid spacing. There are 45 vertical levels with higher resolution in the planetary boundary layer.1 The model physics options include the Goddard microphysics scheme (Tao et al. 2003), the Kain–Fritsch cumulus parameterization scheme (Kain and Fritsch 1990), the Yonsei University (YSU) planetary boundary layer scheme (Hong et al. 2006), the Noah land surface model (Chen and Dudhia 2001), and the Rapid Radiative Transfer Model (RRTM) longwave (Mlawer et al. 1997) and Goddard shortwave (Chou and Suarez 1994) radiation schemes. For this version of WRF, the RRTM longwave radiation physics scheme was noted to produce an upper-level cold bias in the forecasts (Cavallo et al. 2011). Detailed descriptions of the WRF can be obtained from Skamarock et al. (2008).

For WRF 3DVAR, we adopt a partial cycling approach similar to that in NDAS. In NDAS, the partial cycling includes a cold start from 12 h prior to the analysis time and subsequent data assimilation cycling every 3 h. For our study, the partial cycling begins with a cold start that is initialized from NCEP Global Forecast System (GFS) analyses, followed by two update cycles using the 6-h WRF forecasts from the previous cycle. We do not follow the partial cycling strategy of NDAS exactly in WRF 3DVAR due to concerns about the imbalances between mass and wind fields in a 3-h forecast. We investigate the initial imbalance characteristics in short-term WRF forecasts based on WRF 3DVAR analysis using the tendency of surface pressure (Lynch and Huang 1992). The mean absolute surface pressure tendency, N, is defined as

 
formula

where ps is the surface pressure and the summation denotes calculation over the whole model domain. In addition, N is a characteristic quantity that reflects the overall balance of the model state. The forecast pressure tendency (from a 6-h data assimilation cycle) starting at 1800 UTC 26 September 2008 is selected as an example to illustrate the evolutions of N in the first 9-h forecasts (Fig. 1). The pressure tendency is reduced from 1.75 to 1 hPa (3 h)−1 in the first 4 h. Afterward, it is stabilized with the value of about 1 hPa (3 h)−1. Because of the balance constraints built into the preconditioning of the cost-function minimization of the WRF 3DVAR system, the initial imbalance tends to disappear after 4 h. Therefore, we use a 6-h assimilation update cycle for partial cycling, instead of the 3-h assimilation update cycle, to minimize the impact of the imbalance in the background fields.

Fig. 1.

The 9-h evolutions of the mean absolute surface pressure tendency N [hPa (3 h)−1] averaged over the forecast area of WRF forecasts.

Fig. 1.

The 9-h evolutions of the mean absolute surface pressure tendency N [hPa (3 h)−1] averaged over the forecast area of WRF forecasts.

The basic configuration of WRF 3DVAR (Skamarock et al. 2008) is based on the iterative minimization of the cost function:

 
formula

where x, xb, and yo are the analysis state, first guess (or background), and observation, respectively. In addition, and are background and observation error covariance matrices, respectively. Comparison against observations (yo) is made possible by y = H(x), which transforms the gridded analysis x to observation space using the nonlinear observation operator (H). The “cv5” option (for background error statistics), which is adopted here, includes the following physical-space control variables: streamfunction, unbalanced velocity potential, unbalanced surface pressure, unbalanced temperature, and “pseudo–” relative humidity. Under the assumption that the errors cross correlations are negligible, the control variable space v, which is used by WRF 3DVAR, is related to model space through the control variable transforms :

 
formula

where = T. In terms of control variable and incremental formulation (Courtier et al. 1994) via quadratic approximation, the cost function (1) can be written as

 
formula

where d = yoH(xb) is the departure from observation (innovation) and is the linearized nonlinear observation operator (H) in the form of a first-order Taylor series expansion. The gradient of cost function (3) is used to solve the minimization problem as

 
formula

Thus, the following analysis can be obtained when the first derivative (4) is equal to zero:

 
formula

With the inner loop minimization algorithm, the increment is computed from (5) and added back to the control variable (v). As a result, the control variable can be modified with each minimization process (or outer loop). In addition, nonlinearities in the observation operator can be included and observation quality control is repeated in each outer loop. Observations that were rejected in the previous outer loop may be accepted in the subsequent outer loop.

To illustrate the effect of the outer loop with the use of the nonlinear observation operator, we perform a single observation test. Figure 2 gives an example of the evolution of the cost function and gradient norm as a function of minimization iterations. The case is the single GPS RO refractivity observation test with a value of 10 units on the 20th model level at 1800 UTC 26 September 2008. It clearly shows the large reduction in the cost function in the first outer loop, as well as the gradient norm, and the WRF 3DVAR minimization converges within the third outer loop. The single observation example indicates that three outer loops may be sufficient for a nonlinear observation operator, such as GPS RO (which is the only highly nonlinear operator). It is noted that the original proposers for the outer loop strategy (Courtier et al. 1994; Veerse and Thépaut 1998) tested four outer loops in their pioneering work, while the operational GSI in the NCEP GFS implemented two outer loops (Kleist et al. 2009). Massart et al. (2010) showed that more than one outer loop can significantly improve the analysis when the background state is inaccurate. Based on these previous studies, which proposed two to four outer loops, we decided to use three outer loops in our WRF 3DVAR simulations.

Fig. 2.

The (a) cost functions and (b) gradient norm as functions of minimization iterations from Eq. (5) with outer loops. Numbers 1–5 denote the time of the outer loop.

Fig. 2.

The (a) cost functions and (b) gradient norm as functions of minimization iterations from Eq. (5) with outer loops. Numbers 1–5 denote the time of the outer loop.

Typhoon vortex initialization is an important component for numerical typhoon prediction. The typhoon initialization procedure in this study is based on the work of Hsiao et al. (2010), in which a vortex relocation scheme is utilized to minimize the typhoon position error in the first guess. Meanwhile, synthetic data from a bogus vortex are assimilated as part of the 3DVAR analysis to better define the typhoon structure.

Four sets of data assimilation experiments are conducted, each with a 72-h forecast integration, to assess the impact of partial cycling and the outer loop approaches on typhoon track forecasts. For the full cycling experiment, the initial analysis is cold started (using the NCEP GFS as the first guess) from 0600 UTC 8 September 2008. Subsequently, it carries its own 6-h data assimilation cycle with the first guess for the 3DVAR analysis derived from the 6-h forecast of the previous run. For partial cycling, the analysis is generated by a cold start at 12 h prior to the analysis time (using NCEP GFS analysis as the first guess), followed by two data assimilation cycles at 6-h intervals. Therefore, there are four sets of numerical experiments: full cycling runs with (FWO) and without (FNO) outer loops; partial cycling runs with (PWO) and without (PNO) outer loops, for a total of 78 cases (1800 UTC 8 September–0000 UTC 28 September 2008) in each experiment set. Furthermore, the full cycling experiment without an outer loop (FNO) serves as the control run to assess the impact of partial cycling and outer loops on typhoon track forecasting.

3. Experimental results

a. Synopsis of typhoon cases

During the experimental period, there were three typhoons over the western North Pacific: Sinlaku, Hagupit, and Jangmi. In addition, both Sinlaku and Jangmi reached supertyphoon intensity before making landfall on Taiwan. Typhoon Hagupit only reached typhoon intensity. Typhoon Sinlaku (Fig. 3a) developed on 8 September 2008 east of the Philippines and moved slowly north-northwestward to the west of the subtropical high. It made landfall on 14 September and then recurved between Taiwan and China. Finally, Sinlaku sped up gradually toward Japan along the western edge of the subtropical high. On 19 September, Hagupit was upgraded from a tropical depression to a tropical storm east of the Philippines (Fig. 3b). It moved steadily northwestward and intensified into a typhoon 2 days later. Hagupit continued to intensify and approached southern China; it finally struck Guangdong Province and dissipated quickly on 24 September. At the same time, Typhoon Jangmi (Fig. 3c) formed between Guam and the Philippines, then moved to the northwest and made landfall on 28 September. After passing through Taiwan, Jangmi experienced a rapid weakening, recurved north-northeastward and later dissipated to the south of Japan.

Fig. 3.

CWB best-track positions for Typhoons (a) Sinlaku, (b) Hagupit, and (c) Jangmi plotted every 6 h with labels indicating the September date of the 0000 UTC position in 2008. The maximum surface wind speed (m s−1; from CWB) is also indicated.

Fig. 3.

CWB best-track positions for Typhoons (a) Sinlaku, (b) Hagupit, and (c) Jangmi plotted every 6 h with labels indicating the September date of the 0000 UTC position in 2008. The maximum surface wind speed (m s−1; from CWB) is also indicated.

b. Analysis of the track forecasts

Figure 4 shows the mean track errors of the four sets of experiments and the impact of partial cycling and the outer loops. A Student’s t test was used to measure the statistical significance of the mean track error, and the 95% confidence level was considered the significance delimiter (Table 1). Note that tests would have to result in a probability value (p value) of 0.05 or less before the results could be considered significant. The results show that all the statistics are significant at the 95% confidence level except for the differences between full and partial cycling (FNO versus PNO) during the early forecast hours.

Fig. 4.

Mean typhoon track errors for FNO (black bars), FWO (white bars), PNO (dark gray bars), and PWO (light gray bars) experiments and the percentage of the mean track error improvement seen by the outer loops (FWO vs FNO with gray closed circle), partial cycling (PNO vs FNO with open circle), and both of their influences (PWO vs FNO with black closed circle) from 78 WRF-ARW forecasts.

Fig. 4.

Mean typhoon track errors for FNO (black bars), FWO (white bars), PNO (dark gray bars), and PWO (light gray bars) experiments and the percentage of the mean track error improvement seen by the outer loops (FWO vs FNO with gray closed circle), partial cycling (PNO vs FNO with open circle), and both of their influences (PWO vs FNO with black closed circle) from 78 WRF-ARW forecasts.

Table 1.

The probability values calculated from Fig. 4 for comparison with FNO using a Student’s t test.

The probability values calculated from Fig. 4 for comparison with FNO using a Student’s t test.
The probability values calculated from Fig. 4 for comparison with FNO using a Student’s t test.

The track error from PWO (with partial cycling and outer loop) is smaller than the other three experiments throughout the 72-h forecasts. The impacts of outer loops, partial cycling, and the combination of outer loops and partial cycling are assessed individually from the comparison with the control run (FNO). With the use of outer loops, the track error is reduced by more than 60% at the initial time. This improvement gradually decreases with the forecast time and reaches a minimum of 28% at the 60-h forecast. On the other hand, the improvement in track forecast due to partial cycling is statistically significant after the 24-h forecast, and peaks at 22% at the 36-h forecast. It is clearly shown that the impact of outer loops on typhoon track forecast is statistically larger than that of partial cycling. The combination of outer loops and partial cycling (PWO) shows the largest improvement, with more than 30% averaged over the 72-h forecast period and with a 64% improvement at the initial time. The use of both outer loops with partial cycling yields the greatest reduction in track forecast error.

The performance of individual experiments relative to the control run (FNO) from 0- to 72-h forecasts is summarized in Table 2. In most cases the track errors from the three experiments are lower than those of FNO, indicating that both outer loops and partial cycling contribute to improved track forecasts. For the FWO experiment, 333 cases (out of a total of 488) perform better than the control run. This is consistent with the average forecast improvement shown earlier, indicating the positive impact of the outer loop method. For the PNO experiment, the influence of the partial cycling in itself (without outer loop) is mostly neutral with a slight positive impact. The impact of partial cycling on track forecast increases after the 24-h forecasts, since the percentages of cases with lower track errors are 53%, 62%, and 61% during the 0–24-, 36–48-, and 60–72-h forecast periods, respectively. Combining both approaches in WRF 3DVAR shows greater improvement in the track forecasts than does using either approach alone. The results are consistent with the track forecast errors as a function of forecast time shown in Fig. 4.

Table 2.

The numbers of cases with track forecast errors smaller and larger than the errors in the control run (FNO).

The numbers of cases with track forecast errors smaller and larger than the errors in the control run (FNO).
The numbers of cases with track forecast errors smaller and larger than the errors in the control run (FNO).

c. Verification of the analysis and forecasts

The performance of the four experiments during the three typhoons’ life cycles is examined both in terms of the analysis and forecasts. The root-mean-square error (RMSE) of horizontal winds (U, V), temperature (T), and the mixing ratio of water vapor (Q) are calculated between model forecasts and NCEP GFS analysis over the model domain. A total of 78 individual 72-h WRF deterministic forecasts are verified against the GFS analyses at 0000, 0600, 1200, and 1800 UTC each day.

The statistical significance of all results is assessed using a bootstrap resampling method (Wilks 1997) where random samples (with replacement) are generated from the distribution of the RMSE. In this implementation, the entire distribution is resampled 10 000 times. Only significance levels greater than 90% are considered meaningful herein. Figure 5 shows the vertical distribution of the mean RMSE of the analysis of U, V, T, and Q averaged over all cases for each experiment. It is clear that the analysis RMSEs of PWO and FWO are significantly smaller than that of FNO for all variables. In general, PNO has smaller RMSEs than FNO in the middle and upper levels for temperature and winds. However, the moisture shows little difference between PNO and FNO, except in the middle level.

Fig. 5.

Vertical profile of the analysis RMSEs of (a) U, (b) V, (c) T, and (d) Q for all experiments. Error bars denote the 5% and 95% percentiles determined from bootstrap resampling.

Fig. 5.

Vertical profile of the analysis RMSEs of (a) U, (b) V, (c) T, and (d) Q for all experiments. Error bars denote the 5% and 95% percentiles determined from bootstrap resampling.

The PWO experiment consistently has the lowest RMSE for both dynamic and thermodynamic variables, and is followed by FWO, PNO, and FNO. In addition, the maximum difference between the FNO and PWO RMSEs is about 2 m s−1 for U and V, 0.2 K for T, and 0.7 g kg−1 for Q. The RMSEs of PWO for U, V, T, and Q are 28%, 28%, 21%, and 24% smaller than in the control run, respectively (Figs. 5a–d). Note that the analysis of all variables for experiments with outer loops in WRF 3DVAR almost always has a lower RMSE regardless of whether partial cycling is used or not. In other words, for the initial analysis, the use of outer loops has a greater impact than the use of partial cycling. However, the best analysis is obtained when both outer loop and partial cycling approaches are used. These results are consistent with the typhoon position error at the initial time shown earlier. The greatest improvement in temperature due to partial cycling is found in the middle levels between 900 and 400 hPa, similar to the results of Rogers et al. (2009).

Figure 6 summarizes the vertically averaged RMSEs for each experiment averaged over a total of 78 cases from 0- to 72-h forecast at 24-h intervals. The evaluation of statistical significance of domain-averaged forecast RMSEs, as well as the analysis RMSEs (Fig. 5), based on bootstrap resampling, indicates that differences due to the outer loop scheme (FWO), as well as the combined outer loop and partial cycling approach (PWO), are statistically significant. Except for the moisture field, results of the partial cycling (PNO) experiment are statistically different from those of the full cycling experiment (FNO). In addition, PWO consistently has smaller RMSEs than do the other three experiments for horizontal winds, temperature, and moisture fields at the initial time as well as the subsequent forecast period. In general, the errors increase with forecast time for all variables, with the differences among the four experiments decreasing with time. These results indicate that the outer loop significantly improves the analysis, with its impacts gradually decreasing with forecast time, and that partial cycling only slightly improves (less than 10%) the analysis and forecasts for the wind and temperature fields. However, there is an RMSE reduction in the horizontal wind field throughout the forecasts, which is particularly apparent in the 48- and 72-h forecasts and is attributable to the use of partial cycling. In summary, the use of outer loops improves the analysis and forecasts, particularly during the early stages of the forecast. On the contrary, partial cycling leads to reduced forecast errors for wind and temperature fields, especially after the 48-h forecast, despite the fact that the initial improvement is relatively modest.

Fig. 6.

Domain-averaged (in both horizontal and vertical grids) forecast RMSEs of (a) U, (b) V, (c) T, and (d) Q for all experiments. Error bars denote the 5% and 95% percentiles determined from bootstrap resampling.

Fig. 6.

Domain-averaged (in both horizontal and vertical grids) forecast RMSEs of (a) U, (b) V, (c) T, and (d) Q for all experiments. Error bars denote the 5% and 95% percentiles determined from bootstrap resampling.

d. Observations assimilated

To examine the impact of outer loops on the analyses and forecasts, the usage of various types of observations by the WRF 3DVAR system is examined. These data types include 12 categories: surface synoptic observations (SYNOP—surface pressure, temperature, wind, relative humidity), buoys (BUOY—surface pressure, temperature, wind, relative humidity), aviation routine weather reports (METAR—wind, temperature), Quick Scatterometer (QSCAT—surface wind), ships (SHIP—surface pressure, temperature, wind, relative humidity), aviation reports (AIREP—wind, temperature), pilot reports (PILOT—wind), soundings (SOUND—surface pressure, temperature, wind, relative humidity), bogussing (BOGUS—temperature, wind, relative humidity), GPS (refractivity), geosynchronous atmospheric motion vectors (GEOAMV—wind), and satellite temperature (SATEM—temperature). The quality control checks are applied to these data prior to being assimilated by the WRF 3DVAR. These checks include removing observations outside the domain, excluding location/time duplicates and incomplete observations (e.g., no location), and ensuring the vertical consistency of upper-air profiles (such as removing superadiabatic lapse rates from temperature profiles, i.e., not allowing the potential temperature to decrease with the height). Observations whose innovations are greater than 5 times the assumed observation error standard deviation are rejected (Skamarock et al. 2008).

Figure 7 displays the percent usage of all the observation types mentioned earlier for all experiments. Note that, since the bogus vortex data are always assimilated (in all experiments) to better define the typhoon structure with dynamic consistency, they are excluded from the statistics presented in Fig. 7. For other types of observations, more than 90% of the data are used in each experiment. In general, more than 95% of the observations are assimilated, and the upper (lower) quartile in all experiments is more than 97% (94%). The experiment that has the lowest data usage is the control run (FNO) with a 96.3% mean value for all data types; data usage slightly increases in the PNO experiment (96.6%) with the use of partial cycling. With the use of outer loops there is an increase of about 2% in data usage. Similar differences are also reflected in the minimum and median values of data usage. Although the overall observation increase from 96% to 98% may not appear to be significant, this represents more than 700 additional observations that are used in the assimilation. In addition, the largest difference (about 300) among the observations is GPS RO refractivity. This indicates that close to half of the additional observations assimilated are the GPS RO data.

Fig. 7.

Box-and-whisker plot of the percentage of observations used as derived from 12 categories of observations for each experiment. Values from the mean of the observation percentage are indicated with dots. The box-and-whisker plot is interpreted as follows: the middle line shows the median value; the top and bottom of the box show the upper and lower quartiles (i.e., 75th and 25th percentile values), respectively; and the whiskers show the minimum and maximum values.

Fig. 7.

Box-and-whisker plot of the percentage of observations used as derived from 12 categories of observations for each experiment. Values from the mean of the observation percentage are indicated with dots. The box-and-whisker plot is interpreted as follows: the middle line shows the median value; the top and bottom of the box show the upper and lower quartiles (i.e., 75th and 25th percentile values), respectively; and the whiskers show the minimum and maximum values.

Since GPS RO refractivity is nonlinearly related to the analysis control variables, the outer loop method is able to extract more information from such data. To examine the utilization of GPS RO data, the rejection number of GPS RO refractivity with time is shown in Fig. 8. The number is derived from the count of each observation at each vertical level. There is a significant difference between experiments with and without the outer loops in terms of the utilization of GPS RO data. With outer loops, the rejection number is about 100, which is far less than the experiments without outer loops (~400). In general, in terms of the usage of GPS RO data, all the experiments have similar distributions with time during the study period. Figure 9 shows the loci of the additional observations assimilated in FWO (as compared with FNO) for the initial analysis of 78 forecasts. The additional GPS RO data used in FWO are broadly distributed over the domain, with higher density over the higher latitudes, which is consistent with the distribution of GPS RO soundings from the Constellation Observing System for Meteorology, Ionosphere, and Climate (COSMIC) constellation with a 72° inclination orbit (Anthes et al. 2008). Note that the distribution of additional RO soundings used in FWO, shown in Fig. 9, is almost the same over oceans and land. However, other upper-air observations, such as sounding and pilot, are sparse over the ocean. The assimilation of additional GPS RO observations over the tropical and subtropical oceans (where there is a lack of the traditional observations) would help improve the quality of the analysis.

Fig. 8.

Time series of the rejection numbers of GPS refractivity for all experiments.

Fig. 8.

Time series of the rejection numbers of GPS refractivity for all experiments.

Fig. 9.

The distribution of the difference between FNO and FWO in the RO sounding usage (gray dots). Star and cross symbols denote differences in sounding and pilot observation usage patterns, respectively.

Fig. 9.

The distribution of the difference between FNO and FWO in the RO sounding usage (gray dots). Star and cross symbols denote differences in sounding and pilot observation usage patterns, respectively.

Further investigation indicates that 85% of the additional GPS RO observations (which were rejected in the first outer loop) were recovered during the second outer loop, and the remaining 15% were captured during the third outer loop. This is consistent with Rizvi et al. (2008) in that the outer loop allows observations rejected in previous iterations to be accepted when the updated analysis becomes closer to the observation. The outer loop approach is particularly beneficial to the assimilation of nonlinear variables such as GPS RO refractivity. With more observations, particularly the GPS RO data, used in the WRF 3VAR system, improved analysis and subsequent forecasts are obtained.

e. Verification of first guess

Because of the modest impact of partial cycling on track forecast error and gridpoint analysis and forecast errors, we now examine the time series of vertically averaged RMSEs for PNO versus FNO to investigate how partial cycling influences the quality of the initial state (Fig. 10). A Student’s t test is used to measure the statistical significance of partial cycling, and the 95% confidence level is considered significant herein. The probability value (p value) for the comparison between FNO and PNO is less than 0.001 for the U, V, and T fields, and 0.33 for the Q field. The partial cycling experiment (PNO) has smaller RMSEs than does the full cycling version (FNO), which is statistically significant, for wind and temperature fields at nearly all times. However, the difference between these two experiments does not increase substantially with time. The reason for this might be either the study period is not long enough to accumulate the forecast errors, or the assimilation of observations leads to an improved first guess and consequently the model error of the full cycling approach becomes stabilized. To illustrate the influence of data assimilation on the full cycling run, we conduct an additional full cycling experiment (FND) without the assimilation of observations. As expected, the FND experiment has a larger RMSE and the error increases steadily with time. This suggests that the assimilation of observations plays an important role in controlling the growth of the analysis error.

Fig. 10.

Time series of the vertically averaged first-guess RMSEs of (a) U, (b) V, (c) T, and (d) Q for FNO, PNO, and FND.

Fig. 10.

Time series of the vertically averaged first-guess RMSEs of (a) U, (b) V, (c) T, and (d) Q for FNO, PNO, and FND.

To further examine the impact of partial cycling, we show in Fig. 11 the time evolution of the domain-averaged absolute bias for PNO and the difference between FNO and PNO. The positive differences between these two experiments indicate that PNO has fewer errors than FNO. In addition, the horizontal wind errors increase with forecast time and height. Comparison of these two experiments shows that PNO has smaller bias at all levels throughout the 72-h forecasts. For the temperature field, the largest difference between FNO and PNO is located at the middle levels (from 900 to 400 hPa) at the initial time. Despite the relatively small difference between the partial and full cycling methods, the verification of storm track, analysis, and forecasts all exhibit better performance with the use of the partial cycling approach.

Fig. 11.

Domain-averaged forecast absolute biases of (a) U, (b) V, (c) T, and (d) Q for PNO experiments (shaded). The difference between FNO and PNO (FNO − PNO) is plotted with contours.

Fig. 11.

Domain-averaged forecast absolute biases of (a) U, (b) V, (c) T, and (d) Q for PNO experiments (shaded). The difference between FNO and PNO (FNO − PNO) is plotted with contours.

To gain further insights into the impact of partial cycling, we now examine the spatial distribution of the difference fields. Because of the significant impact on winds at 300 hPa and on temperature at 700 hPa (Fig. 5), we show in Fig. 12 the difference between the composite analysis of FNO (PNO) and the NCEP GFS analysis from 78 cases. Overall, PNO has smaller differences than NCEP GFS analysis in winds and temperature. If we consider the differences from the NCEP GFS analysis to be errors, then it is clear that the FNO has a larger warm bias over the western North Pacific and the Indian Ocean. Both FNO and PNO have excessive easterlies over the Indian Ocean, with a larger magnitude in FNO. The reason for the warm bias over the ocean might be attributable to the lack of satellite radiance assimilation. The operational global data assimilation system makes good use of satellite radiance data and, therefore, will be able to correct the model systematic errors more effectively over the ocean where there are few traditional observations. On the contrary, regional models (such as the CWB regional forecast system) that do not assimilate satellite radiance are expected to have a larger systematic model error over the ocean. The assimilation of GPS RO data would help, but may not be sufficient to overcome this problem in regional models. The fact that partial cycling (which always began with a cold start based on NCEP GFS analysis) still has a warm bias over the ocean (admittedly smaller than that of FNO) suggests that the version of WRF used in the CWB’s operations is responsible for creating the systematic warm bias. One possible source is from the use of the Kain–Fritsch scheme without an explicit shallow convection parameterization (Torn and Davis 2012).

Fig. 12.

The differences between the composite analyses of (left) FNO and (right) PNO and NCEP from 78 cases for (a),(b) 700-hPa temperature, and 300-hPa (c),(d) U and (e),(f) V fields.

Fig. 12.

The differences between the composite analyses of (left) FNO and (right) PNO and NCEP from 78 cases for (a),(b) 700-hPa temperature, and 300-hPa (c),(d) U and (e),(f) V fields.

An important difference between the regional and global data assimilations is the existence of a lateral boundary. For regional data assimilation, observations outside the boundary would not have an influence on the analysis inside the model domain. Global data assimilation does not have such an issue. For a cycling data assimilation, an important source of error is the model systematic error. The fact that the systematic difference decreases toward the lateral boundaries (LBC) for both PNO and FNO (Figs. 12a and 12b), when verified against the NCEP GFS analysis, indicates that the forecast lateral boundary conditions (LBCs) from NCEP GFS forecasts help control the growth of the model systematic error. The smaller systematic error over the interior model domain (particularly over the ocean) in PNO, when compared with FNO, suggests that the WRF used in CWB operations does contain a systematic error larger than that of the NCEP GFS forecast. Partial cycling, which continuously brings in improved analyses from NCEP GFS over the ocean, helps to suppress the growth of the analysis error.

4. Conclusions

In this paper, we examine the performance of WRF 3DVAR with the use of outer loops and partial cycling on 78 cases of analyses and forecasts for Typhoons Sinlaku, Hagupit, and Jangmi (2008), for which CWB issued typhoon warnings. Overall, the use of both the outer loop and partial cycling approaches produces the largest improvement in the typhoon track forecast, analysis, and forecasts that are statistically significant. In addition, the outer loop results in a larger improvement in typhoon track forecasts than does partial cycling in the CWB’s operational configuration. Further analysis indicates that outer loops allow more observations, particularly for the GPS RO refractivity, to be used in the assimilation. These additional GPS RO observations are distributed uniformly over both land and ocean, while the additional upper-air sounding observations are only located over land. As a result, more GPS RO observations are assimilated through outer loops, which improves the analysis over the western North Pacific Ocean, where there is a lack of traditional observations.

Verification of WRF 3DVAR analysis with continuous cycling against NCEP GFS analysis averaged over the 78 cases reveals visible systematic temperature errors over the western North Pacific Ocean and easterly wind errors over the Indian Ocean. These systematic errors are most likely caused by WRF physical parameterizations. Partial cycling is shown to be effective in suppressing the growth of systematic model error, through periodic cold starts with the NCEP GFS global analysis, followed by two analysis update cycles.

The use of the outer loop method in WRF 3DVAR requires little extra computation for the added minimization of the cost function. There is less than a 1-min wall-clock time difference with the use of three outer loops on an IBM supercomputer with 16 “power5+” processors and 32 GB of memory. The use of partial cycling will require the generation of a cold start at 12 h prior to the analysis time, followed by two data assimilation cycles at 6-h intervals. The extra computational cost is about 20 min with a 64-processor IBM supercomputer. However, the extra computations can be performed in advance prior to the analysis time and will not result in a delay in the operation. Therefore, neither the outer loop nor the partial cycling approach will result in significant additional computational burden, and both can be easily implemented in the CWB’s operational forecast system.

In summary, we have shown that improved analyses and forecasts of typhoons can be achieved with the use of both outer loops and partial cycling in the WRF 3DVAR data assimilation system. The operational CWB typhoon forecast system now includes these features, together with the vortex implantation procedure implemented by Hsiao et al. (2010). This system has been found to give superior performance for typhoon track forecasting with the CWB operational forecasting system. Because the configurations of this WRF-3DVAR/WRF model system at CWB are focused on typhoon forecasting, it is designated as TWRF, which stands for the typhoon WRF. The TWRF modeling system became operational at CWB for typhoon forecasting over the western North Pacific region in 2010.

Acknowledgments

This study was supported by the National Science Council of the R.O.C. under Grants NSC98-2625-M-052-008 and NSC99-2625-M-052-003-MY3, American Institute in Taiwan under Award CWB05-071-119, and the U.S. National Science Foundation under Cooperative Agreement AGS-0918398/CSA No. AGS-0939962. We thank Dr. Xingqin Fang and Dr. Wan-Shu Wu for reading an earlier version of this manuscript and offering many thoughtful comments and suggestions. We also thank the three anonymous reviewers for their constructive comments and suggestions on the manuscript, as well as Ms. Rebecca Nowak and Ms. Ingrid Moore for technical editing of the manuscript.

REFERENCES

REFERENCES
Andersson
,
E.
, and
H.
Järvinen
,
1999
:
Variational quality control
.
Quart. J. Roy. Meteor. Soc.
,
125
,
697
722
.
Anthes
,
R. A.
, and
Coauthors
,
2008
:
The COSMIC/FORMOSAT-3 mission: Early results
.
Bull. Amer. Meteor. Soc.
,
89
,
313
333
.
Cavallo
,
S. M.
,
J.
Dudhia
, and
C.
Snyder
,
2011
:
A multilayer upper-boundary condition for longwave radiative flux to correct temperature biases in a mesoscale model
.
Mon. Wea. Rev.
,
139
,
1952
1959
.
Chen
,
F.
, and
J.
Dudhia
,
2001
:
Coupling an advanced land-surface/hydrology model with the Penn State–NCAR MM5 modeling system. Part I: Model description and implementation
.
Mon. Wea. Rev.
,
129
,
569
585
.
Chen
,
Y.
, and
C.
Snyder
,
2007
:
Assimilating vortex position with an ensemble Kalman filter
.
Mon. Wea. Rev.
,
135
,
1828
1845
.
Chou
,
M.-D.
, and
M. J.
Suarez
,
1994
:
An efficient thermal infrared radiation parameterization for use in general circulation models. NASA Tech. Memo. 104606, Vol. 3, 85 pp
.
Courtier
,
P.
,
J.-N.
Thépaut
, and
A.
Hollingsworth
,
1994
:
A strategy for operational implementation of 4D-Var, using an incremental approach
.
Quart. J. Roy. Meteor. Soc.
,
120
,
1367
1387
.
Courtier
,
P.
, and
Coauthors
,
1998
:
The ECMWF implementation of three-dimensional variational (3DVAR) data assimilation. Part I: Formulation
.
Quart. J. Roy. Meteor. Soc.
,
123
,
1
26
.
Davis
,
C. A.
, and
Coauthors
,
2008
:
Prediction of landfalling hurricanes with the Advanced Hurricane WRF model
.
Mon. Wea. Rev.
,
136
,
1990
2005
.
Guo
,
Y.-R.
,
H.-C.
Lin
,
X. X.
Ma
,
X.-Y.
Huang
,
C. T.
Terng
, and
Y.-H.
Kuo
,
2006
:
Impact of WRF-Var (3DVar) background error statistics on typhoon analysis and forecast. WRF Users’ Workshop, Boulder, CO, NCAR, P4.2. [Available online at http://www.mmm.ucar.edu/wrf/users/workshops/WS2006/abstracts/PSession04/P4_2_Guo.pdf.]
Hodur
,
R. M.
,
1987
:
Evaluation of a regional model with an update cycle
.
Mon. Wea. Rev.
,
115
,
2707
2718
.
Holland
,
G. J.
, and
R. T.
Merrill
,
1984
:
On the dynamics of tropical cyclone structure changes
.
Quart. J. Roy. Meteor. Soc.
,
110
,
723
745
.
Hong
,
S.-Y.
,
Y.
Noh
, and
J.
Dudhia
,
2006
:
A new vertical diffusion package with an explicit treatment of entrainment processes
.
Mon. Wea. Rev.
,
134
,
2318
2341
.
Hsiao
,
L. F.
,
M. S.
Peng
,
D. S.
Chen
,
K. N.
Huang
, and
T. C.
Yeh
,
2009
:
Sensitivity of typhoon track predictions in a regional prediction system to initial and lateral boundary conditions
.
J. Appl. Meteor. Climatol.
,
48
,
1913
1928
.
Hsiao
,
L. F.
,
C. S.
Liou
,
T. C.
Yeh
,
Y. R.
Guo
,
D. S.
Chen
,
K. N.
Huang
,
C. T.
Terng
, and
J. H.
Chen
,
2010
:
A vortex relocation scheme for tropical cyclone initialization in Advanced Research WRF
.
Mon. Wea. Rev.
,
138
,
3298
3315
.
Kain
,
J. S.
, and
J. M.
Fritsch
,
1990
:
A one-dimensional entraining/detraining plume model and its application in convective parameterization
.
J. Atmos. Sci.
,
47
,
2784
2802
.
Kleist
,
D. T.
,
D. F.
Parrish
,
J. C.
Derber
,
R.
Treadon
,
W.-S.
Wu
, and
S.
Lord
,
2009
:
Introduction of the GSI into the NCEP Global Data Assimilation System
.
Wea. Forecasting
,
24
,
1691
1705
.
Kurihara
,
Y.
,
M. A.
Bender
,
R. E.
Tuleya
, and
R. J.
Ross
,
1995
:
Improvements in the GFDL hurricane prediction system
.
Mon. Wea. Rev.
,
123
,
2791
2801
.
Lin
,
N.
,
J. A.
Smith
,
G.
Villarini
,
T. P.
Marchok
, and
M. L.
Baeck
,
2010
:
Modeling extreme rainfall, winds, and surge from Hurricane Isabel (2003)
.
Wea. Forecasting
,
25
,
1342
1361
.
Liu
,
Y.
,
D.-L.
Zhang
, and
M. K.
Yau
,
1997
:
A multiscale study of Hurricane Andrew (1992). Part I: Explicit simulation and verification
.
Mon. Wea. Rev.
,
125
,
3073
3093
.
Lynch
,
P.
, and
X.-Y.
Huang
,
1992
:
Initialization of the HIRLAM model using a digital filter
.
Mon. Wea. Rev.
,
120
,
1019
1034
.
Marks
,
F. D.
, and
L. K.
Shay
,
1998
:
Landfalling tropical cyclones: Forecast problems and associated research opportunities
.
Bull. Amer. Meteor. Soc.
,
79
,
305
323
.
Massart
,
S.
,
B.
Pajot
,
A.
Piacentini
, and
O.
Pannekoucke
,
2010
:
On the merits of using a 3D-FGAT assimilation scheme with an outer loop for atmospheric situations governed by transport
.
Mon. Wea. Rev.
,
138
,
4509
4522
.
Mlawer
,
E. J.
,
S. J.
Taubman
,
P. D.
Brown
,
M. J.
Iacono
, and
S. A.
Clough
,
1997
:
Radiative transfer for inhomogeneous atmosphere: RRTM, a validated correlated-k model for the longwave
.
J. Geophys. Res.
,
102
(
D14
),
16 663
16 682
.
Molinari
,
J.
, and
D.
Vollaro
,
1990
:
External influences on hurricane intensity. Part II: Vertical structure and response of the hurricane vortex
.
J. Atmos. Sci.
,
47
,
1902
1918
.
Rizvi
,
S. R. H.
,
Y.-R.
Guo
,
H.
Shao
,
M.
Demirtas
, and
X.-Y.
Huang
,
2008
:
Impact of outer loop for WRF data assimilation system (WRFDA). WRF Users’ Workshop, Boulder, CO, NCAR, P5.3. [Available online at http://www.mmm.ucar.edu/wrf/users/workshops/WS2008/abstracts/P5-03.pdf.]
Rogers
,
E.
, and
Coauthors
,
2009
:
The NCEP North American Mesoscale modeling system: Recent changes and future plans. Preprints, 23rd Conf. on Weather Analysis and Forecasting/19th Conf. on Numerical Weather Prediction, Omaha, NE, Amer. Meteor. Soc., 2A4. [Available online at http://ams.confex.com/ams/pdfpapers/154114.pdf.]
Rosmond
,
T.
, and
L.
Xu
,
2006
:
Development of NAVDAS-AR: Non-linear formulation and outer loop tests
.
Tellus
,
58A
,
45
58
.
Ross
,
R. J.
, and
Y.
Kurihara
,
1995
:
A numerical study on influences of Hurricane Gloria (1985) on the environment
.
Mon. Wea. Rev.
,
123
,
332
346
.
Skamarock
,
W. C.
, and
Coauthors
,
2008
:
A description of the Advanced Research WRF version 3. NCAR Tech. Note NCAR/TN-475+STR, 113 pp
.
Tao
,
W.-K.
, and
Coauthors
,
2003
:
Microphysics, radiation and surface processes in the Goddard Cumulus Ensemble (GCE) model
.
Meteor. Atmos. Phys.
,
82
,
97
137
.
Torn
,
R. D.
,
2010
:
Performance of a mesoscale ensemble Kalman filter (EnKF) during the NOAA High-Resolution Hurricane Test
.
Mon. Wea. Rev.
,
138
,
4375
4392
.
Torn
,
R. D.
, and
G. J.
Hakim
,
2009
:
Ensemble data assimilation applied to RAINEX observations of Hurricane Katrina (2005)
.
Mon. Wea. Rev.
,
137
,
2817
2829
.
Torn
,
R. D.
, and
C.
Davis
,
2012
:
The influence of shallow convection on tropical cyclone track forecasts
.
Mon. Wea. Rev.
,
140
,
2188
2197
.
Veerse
,
F.
, and
J.-N.
Thépaut
,
1998
:
Multiple-truncated incremental approach for four-dimensional variational assimilation
.
Quart. J. Roy. Meteor. Soc.
,
124
,
1889
1908
.
Wilks
,
D. S.
,
1997
:
Resampling hypothesis tests for autocorrelated fields
.
J. Climate
,
10
,
65
82
.
Wu
,
C. C.
, and
Y.-H.
Kuo
,
1999
:
Typhoons affecting Taiwan: Current understanding and future challenges
.
Bull. Amer. Meteor. Soc.
,
80
,
67
80
.
Xiao
,
Q.
,
Y.-H.
Kuo
,
J.
Sun
,
W. C.
Lee
,
E.
Lim
,
Y. R.
Guo
, and
D. M.
Barker
,
2005
:
Assimilation of Doppler radar observations with a regional 3DVAR system: Impact of Doppler velocities on forecasts of a heavy rainfall case
.
J. Appl. Meteor.
,
44
,
768
788
.
Zhao
,
Q.
,
J.
Cook
,
Q.
Xu
, and
P.
Harasti
,
2006
:
Using radar wind observations to improve mesoscale numerical weather prediction
.
Wea. Forecasting
,
21
,
502
522
.

Footnotes

1

The 45 vertical levels with sigma values in the terrain-following coordinate are 1.0, 0.995, 0.988, 0.98, 0.97, 0.96, 0.945, 0.93, 0.91, 0.89, 0.87, 0.85, 0.82, 0.79, 0.76, 0.73, 0.69, 0.65, 0.61, 0.57, 0.53, 0.49, 0.45, 0.41, 0.37, 0.34, 0.31, 0.28, 0.26, 0.24, 0.22, 0.2, 0.18, 0.16, 0.14, 0.12, 0.1, 0.082, 0.066, 0.052, 0.04, 0.03, 0.02, 0.01, and 0.0.