Assimilating AMSU-A Radiance Data with the WRF Hybrid En3DVAR System for Track Predictions of Typhoon Megi (2010)

Received 28 October 2014; revised 12 December 2014; accepted 29 December 2014)

The accuracy of tropical cyclone (TC) track prediction is one of the most important tasks in weather forecasting, crucial for saving lives and property. Over the past two decades, significant improvements have been made in TC track forecasts due to the increased use of remote sensing data, other advanced observations, and improved numerical weather prediction (NWP) models (Houze et al., 2007; Singh et al., 2008; Rappaport et al., 2009; Cangialosi and Franklin, 2011). However, great challenges remain, especially for recurving tracks (George and Gray, 1977; Thu and Krishnamurti, 1992; Holland and Wang, 1995; Zhang et al., 2013b), because a recurving track of a TC involves more uncertainties than an non-recurving one. As one of the most important types of observations for NWP, satellite radiance data can provide valuable temperature and humidity information, especially over areas where conventional observations are limited (Derber and Wu, 1998; English et al., 2000; McNally et al., 2000; Bouttier and Kelly, 2001; Le Marshall et al., 2006; McNally et al., 2006). Satellite radiances can be assimilated into NWP models with retrieval assimilation or direct assimilation (DA) methods. With the retrieval assimilation method, the retrieved temperature and humidity information from radiances using physical or statistical retrieval methods (Goldberg, 1999) are applied in a similar way as with conventional data. The DA method involves assimilating radiance observations directly into NWP models, and is considered superior to retrieval assimilation because the observational error statistics are more accurate (Eyre, 1989; Andersson et al., 1994; Derber and Wu, 1998; English et al., 2000; Bouttier and Kelly, 2001). Moreover, its rapid, real-time processing without the retrieval steps is an advantage for operational data usage, enabling the early introduction of the radiance data into operational systems.

For most operational centers, a three-dimensional variational (3DVAR; e.g., Parrish and Derber, 1992; Lorenc et al., 2000; Barker et al., 2004) DA scheme is employed to obtain the initial conditions. In 3DVAR, observational information is spread to model grid points to correct the background field with the assumption of an isotropic, nearly homogeneous, static, time-invariant with flow-independent background error covariances (BECs). However, these assumptions cannot generally be used for TC assimilation applications due to the strong vortical and nongeostrophic motions of TCs (Hamill and Snyder, 2000); the true flow-dependent nature of the background errors is not captured. Recently, several studies have demonstrated that forecasts from ensemble-based DA can produce comparable or better forecasts than those from 3DVAR in a variety of weather applications (e.g., Li and Liu, 2009; Torn and Hakim, 2009; Hamill et al., 2011; Weng et al., 2011; Zhang et al., 2011; Dong and Xue, 2012).

Ensemble Kalman filter (EnKF) methods benefit from the use of flow-dependent BEC but suffer from rank deficiency due to the smaller ensemble members used. The covariance localization method is often adopted to alleviate this problem. The variational technique is quite efficient and able to process complex nonlinear observations and apply physical constraints, but it is lacking in terms of the inclusion of flow-dependent BECs. To reconcile the advantages and disadvantages of the variational and EnKF methods, increasing effort is being made to hybridize the two approaches, rather than settling for one particular method. A hybrid DA approach that couples the ensemble DA technique into the variational framework (e.g., Hamill and Snyder, 2000; Lorenc, 2003; Wang et al., 2008a; Zhang et al., 2009; Wang, 2010; Schwartz et al., 2013) has shown great promise for weather applications. Many previous studies have shown that hybrid approaches yield comparable or better forecasts than those based purely on 3DVAR that do not incorporate ensemble BECs, and can outperform forecasts initialized by standalone EnKFs (e.g., Buehner, 2005; Buehner et al., 2010; Hamill et al., 2011; Wang, 2011; Li et al., 2012; Zhang and Zhang, 2012; Schwartz et al., 2013; Wang et al., 2013; Zhang et al., 2013a; Pan et al., 2014; Schwartz and Liu, 2014).

The hybrid approach is popular for the following three reasons. First, the hybrid technique can be easily implemented in pre-existing variational DA frameworks (e.g., Wang et al., 2007b; Zhang et al., 2013a; Pan et al., 2014). Second, with a model-space covariance localization technique, the assimilation of non-local observations, such as satellite radiance data, may be more effective in hybrid frameworks than in EnKFs that use observation-space localization (Campbell et al., 2010). And third, the hybrid method can save computational cost by using the ensembles at a coarser resolution than deterministic hybrid analysis (e.g., Rainwater and Hunt, 2013), and by producing similar results as EnKFs but with a smaller ensemble compared to the traditional EnKF.

However, to the best of our knowledge, no study has yet been published that applies the Weather Research and Forecasting (WRF) Hybrid En3DVAR DA to the assimilation of radiance data from the Advanced Microwave Sounding Unit-A (AMSU-A) for improving TC track prediction. Accordingly, this preliminary study examines the potential benefits of this computationally efficient procedure for improving TC track forecasting when assimilating AMSU-A radiance data.

More specifically, this study investigates the impacts of DA on TC track forecasting within the WRF Hybrid DA system (Barker et al., 2012), similar to in (Wang, 2011). However, this work differs from that of (Wang, 2011) in several important ways. First, this work investigates the effect of assimilating AMSU-A radiance observations on track predictions, while (Wang, 2011) did not assimilate any satellite radiance data. Second, in this study, each ensemble member is updated by running the hybrid analysis system multiple times with perturbed observations, whereas (Wang, 2011) used an ensemble transformation Kalman filter (ETKF) to generate the analysis ensemble. Moreover, we systematically examined the influence of the static BEC related length scale factor and variance scale factor on the prediction of the track of a typhoon for radiance DA.

The remainder of this paper is organized as follows. In section 2, we provide a brief introduction to the WRF Hybrid DA system system, together with aspects on radiance data assimilation. An overview of Typhoon Megi (2010) and the experimental settings are described in section 3. Section 4 presents the main results. Conclusions and further discussion are provided in the last section.

The WRF Hybrid DA system works by including the extended control variables originally proposed by (Lorenc, 2003). The theoretical details can be found in Wang et al. (2007a, 2008b). We list the most important mathematical formulations in this section. The final analysis increment x' is a sum of two terms, described as

\begin{equation} \label{eq1} \bm{x}'=\bm{x}'_1+\sum_{k=1}^K(\bm{a}_k\circ \bm{x}_{{\rm e},k}) . (1)\end{equation}

On the right hand side of Eq. (2), the first term x'₁ is the increment associated with the static BEC in 3DVAR. The second term of Eq. (2) represents a local linear combination of ensemble perturbations, and a is the increment associated with the flow-dependent ensemble covariance. Vectors a_k(k=1,…,K) are the extended control variables (Lorenc, 2003) for K ensemble members, and variable x _e,k is the kth ensemble perturbation normalized by \(\sqrtK-1\). Symbol ° denotes the element-by-element product of the vectors a_k and x _e,k. The coefficients of a_k for each ensemble member vary in space due to the necessary ensemble covariance localization, which is conducted in the model state variable space. Without localization, the coefficients of a_k will be the same and can be represented by scalars (Lorenc, 2003). The analysis increment x'₁ and the coefficients of a_k are obtained by minimizing the following hybrid cost function:

\begin{eqnarray} \label{eq2} J(\bm{x}'_1,\bm{a})&=&\beta_1J_1+\beta_2J_{\rm e}+J_{\rm o}\nonumber\\[0.5mm] &=&\beta_1\dfrac{1}{2}(\bm{x}'_1)^{\rm T}\bm{B}^{-1}(\bm{x}'_1)+\beta_2\dfrac{1}{2}(\bm{a})^{\rm T}\bm{A}^{-1}(\bm{a})+\nonumber\\[0.5mm] &&\dfrac{1}{2}(\bm{y}'_{\circ}-\bm{Hx}')^{\rm T}\bm{R}^{-1}(\bm{y}'_{\circ}-\bm{Hx}') , (2)\end{eqnarray}

where a sum of J₁ and J _e with different weight replaces the usual WRF 3DVAR background term J₁. Here, J₁ is associated with the static covariance matrix B. J _e is associated with the ensemble covariance for the extended control vector a, which is formed by concatenating K vectors a_k, k=1,…, K as a ^T=(a₁ ^T,a₂ ^T,…,a _K ^T), which in turn is controlled by the block diagonal matrix A for spatial correlation. J _o is the observation term associated with the observation error covariance matrix R. As in 3DVAR, y'_°=y_°-H(x _b) is the innovation vector. Here, y_° is the observation vector, x _b is the background state vector, and H is the observation operator. In Eq. (3), H is the linearized observation operator. Two coefficients, β₁ and β₂, determine respective weights given to the flow-dependent ensemble covariance and static covariance (Hamill and Snyder, 2000; Etherton and Bishop, 2004; Wang et al., 2007b), which are constrained as

\begin{equation} \label{eq3} \dfrac{1}{\beta_1}+\dfrac{1}{\beta_2}=1 . (3)\end{equation}

The AMSU-A 1b radiance data are ingested into the WRF 3DVAR and WRF Hybrid DA system in this study. AMSU-A is a line-scanned microwave sensor with 15 sensitive channels, each with a 2343 km swath width. It measures 30 pixels in each swath, with an approximate 48 km diameter footprint at nadir. In this study, a subset of AMSU-A channels is chosen to be assimilated. Channels 1-2 and 15 are located in window regions and are thus not assimilated since they are sensitive to uncertain surface parameters, cloud and precipitation. Channels 3-14 are sensitive to temperature, among which only channels 5-7 are assimilated because they peak under the model top (20 hPa).

The Community Radiative Transfer Model (CRTM; Han et al., 2006; Liu and Weng, 2006) is coupled within the WRFDA system (Barker et al., 2012) described in section 2.1 as the observation operator for AMSU-A radiance. This is then used to calculate the simulated radiances with the temperature and moisture information from model states. Radiance data over mixed surfaces and observations with large scan angles are rejected. A radiance observation with large bias (the bias-corrected observation minus the CRTM modeled radiance), exceeding either 15 K or 5r is rejected, where r is the specified observation error standard deviation for brightness temperature. Following (Liu et al., 2012), we use the full-sky AMSU-A radiances in this study without any special cloud detection procedure. Better results for track forecasts are obtained when the thinning mesh is set to roughly 6-8 times the grid resolution for the Typhoon Megi (2010) case. Thus, we determine a 120 km thinning mesh as the first attempt to study assimilating AMSU-A radiances using hybrid methods. We correct the systematic biases from observed radiances before assimilation using the same method as in (Liu et al., 2012) and (Xu et al., 2013). The radiance bias is expressed inside a modified forward operator with a linear combination of several predictors (the scan position, the square and cube of the scan position, the 1000-300 hPa and 200-50 hPa layer thicknesses, surface skin temperature, and total column water vapor) and their coefficients. The coefficients are updated via a variational minimization process by including them in the cost function [Eq. (3)] as control variables (Derber and Wu, 1998; Aulign\'e et al., 2007; Dee and Uppala, 2009).

Super Typhoon Megi (2010)——one of the most destructive TCs over the western North Pacific and South China Sea in 2010——was chosen for our experiments in this study. Megi (2010) was identified as a tropical disturbance by the Joint Typhoon Warning Center (JTWC) when it was about 600 km to the east of the Philippines at 0000 UTC 12 October 2010. Megi (2010) then developed quickly that same day. The JTWC classified the vortex as a tropical depression before 0900 UTC 13 October. Then, later, on 14 October, the Japan Meteorological Agency (JMA) upgraded Megi (2010) to a severe tropical storm and the JTWC upgraded it to a category-1 typhoon. On 15 October, Megi (2010) moved northwestwards and gradually intensified into a typhoon over the Pacific to the east of the Philippines. As shown in Figs. 1 and 2, Megi (2010) initially moved northwestward, and then turned west-southwestward. It experienced a rapid intensification during 16-18 October, reaching its peak intensity at 1200 UTC 17 October, with minimum sea level pressure (MSLP) of 895 hPa and maximum surface wind (MSW) of up to 72 m s^-1. Megi (2010) made landfall over Luzon Island as a super typhoon at 0425 UTC 18 October and became weaker after crossing the island. After that, it re-intensified rapidly from category-2 to category-4 (MSLP: 935 hPa) early on 19 October. At 0000 UTC 20 October, Megi (2010) experienced a sharp turn from westward to northward, an unusual track change that was not forecast by any of the leading operational centers. On 23 October, it made a second landfall as a tropical storm at Zhangpu in Fujian Province, and finally dissipated gradually the next day.

Figure 1.  The WRF model domain and the 6 hourly track of Typhoon Megi (2010) (black typhoon symbols) during 0000 UTC 13 October to 0600 UTC 24 October 2010. Also shown is a snapshot of AMSU-A radiance observations (green dots from NOAA-15, red dots from NOAA-16, and blue dots from METOP-2) assimilated at 0000 UTC 18 October 2010.

Figure 2.  Time series of the minimum sea level pressure (unit: hPa) and the maximum surface winds (units: m s^-1) of Typhoon Megi (2010) from the China Meteorological Administration best-track data.

The WRF model (Skamarock et al., 2008) was used to conduct all the experiments. WRF is a three-dimensional, compressible, non-hydrostatic atmospheric model using terrain-following, mass-based sigma coordinate levels; its governing equations are written in flux form. All experiments were conducted over the single domain with a 15 km horizontal grid spacing. There were 450× 400 grid points horizontally and 43 vertical levels with the model top at 20 hPa. The following parameterizations were used: the WRF Single-Moment 6-Class scheme (Hong et al., 2004); the Kain-Fritsch cumulus parameterization (Kain and Fritsch, 1990, 1993; Kain, 2004) with a modified trigger function (Ma and Tan, 2009); the Yonsei University (YSU) boundary layer scheme (Noh et al., 2003); the 5-layer thermal diffusion model for land surface processes scheme; the RRTM longwave radiation scheme (Mlawer et al., 1997); and the MM5 shortwave radiation scheme (Dudhia, 1989).

Several initial experiments were configured to compare the fundamental difference in using flow-dependent and static BEC and to evaluate the impact of DA when assimilating AMSU-A data on the subsequent forecast of Typhoon Megi (2010). Table 1 summarizes the design of the experiments.

The NOAMA experiment assimilated conventional observations only from the operational Global Telecommunication System dataset of the National Centers for Environmental Prediction (NCEP) with the traditional 3DVAR method. Similar to the NOAMA experiment, the 3DVAR_AMA experiment also assimilated AMSU-A radiances from the NOAA-15, 16, 18 and METOP-2 satellites, besides the conventional observations in NOAMA. Figure 1 shows a snapshot of AMSU-A radiance observations assimilated at 0000 UTC 18 October 2010, for which time AMSU-A radiance data from NOAA-18 are not available. The other two experiments (HYBRID_AMAS and HYBRID_AMAF) employed the hybrid DA method with 40 ensemble members using the mean of the ensemble forecasts as the background. They were identical except with different weighting factors, β₁ and β₂, listed in Table 1. It should also be noted that both 3DVAR_AMA and HYBRID_AMAS used the static background error and the only difference between 3DVAR_AMA and HYBRID_AMAS was the data assimilation background. 3DVAR_AMA used the deterministic forecast as the background, whereas HYBRID_AMAS used the ensemble mean instead for the same purpose. If the difference in the results of 3DVAR_AMA and HYBRID_AMAF are found to be large, while those of 3DVAR_AMA and HYBRID_AMAS are quite similar to each other, it will suggest that this difference is mainly contributed to by the flow-dependent ensemble covariance, rather than the use of the ensemble forecast mean as the background. In addition to the above four experiments, we also conducted experiments with different weighting factors of 75%, 50%, and 25% toward the ensemble contribution. The results did not significantly improve upon the results of HYBRID_AMAF, similar to the findings of (Li et al., 2012). Therefore, we focus on comparing HYBRID_AMAF and 3DVAR in this paper.

Figure 3.  Experimental design of NOAMA, 3DVAR_AMA, HYBRID_AMAS, and HYBRID_AMAF. The 6 h spin-up run was conducted from 1800 UTC 17 through 0000 UTC 18 October 2010. The 6 h cycling forecast-analysis experiments were carried out starting at 1800 UTC and ending at 2300 UTC.

To reduce spurious correlations caused by sampling error due to the limited ensemble member, a 750 km horizontal localization radius was applied to the ensemble covariance, and a standard vertical localization described in Wang et al. [2014, Eq. (12)] was applied. In all the DA experiments, all observations within 1.5 h were assumed to be valid at the analysis time. The static BEC statistics in the 3DVAR system were derived from the differences between the 24 h and 12 h forecasts valid at the same time, using the National Meteorological Center (NMC) method (Parrish and Derber, 1992). These forecasts were generated using the WRF model on the same model grid for each day of September 2010, using the GFS (Global Forecast System) analyses as the initial and boundary conditions. The static BEC matrix was obtained using the WRFDA utility called CV5 (Barker et al., 2012) for five control variables (stream function, unbalanced temperature, surface pressure and velocity potential, and relative humidity). All experiments were initialized at 1800 UTC 17 October 2010 (Fig. 3), using the NCEP operational 0.5^{°× 0.5°} GFS analysis data as the initial and lateral boundary conditions. With the initial and lateral boundary conditions at this time, the initial ensemble was generated by taking Gaussian random draws with the static BECs and zero mean (Torn et al., 2006) and adding the perturbations to the GFS analysis. The deterministic and ensemble fields produced at 1800 UTC 17 October initialized the 6-h WRF forecasts, which served as backgrounds for the first 3DVAR, and hybrid analyses at 0000 UTC 18 October. A 120 h deterministic forecast was separately initialized at 0000 UTC 18 October by the analysis with 3DVAR and the hybrid configurations. The 6 h cycling forecast-analysis experiments were carried out for all experiments until 0000 UTC 23 October. For the following cycles, the background was the previous cycle's 6 h WRF forecast (initialized from the hybrid mean analysis). For those cycles, 24 h deterministic forecasts were carried out to evaluate the impact of cycling assimilation on forecasts. Digital filter initialization (DFI; Lynch and Huang, 1992; Huang and Lynch, 1993) using a twice-DFI scheme and the Dolph filter (Lynch, 1997) with a 2 h backward integration were applied to all 120 h forecasts. Each ensemble member was updated by running the hybrid analysis system multiple times with perturbed observations. The covariance relaxation method of (Zhang et al., 2004) was employed to maintain ensemble spread, where the final inflated ensemble perturbation is a weighted average of prior perturbation (as the background ensemble) and posterior perturbation (the ensemble analysis). In this study, the weight for the posterior perturbation was set to 0.8.

In this section, we evaluate the analyses and forecasts of each of the DA experiments. The model predictions of Megi (2010)'s track were verified against the best-track analyses of the China Meteorological Administration (CMA) (Yu et al., 2007). Aspects of the ensemble spread performance were also examined since the key with ensemble-based DA is the use of an ensemble to estimate the flow-dependent forecast error. The differences in the analyses increments from 3DVAR and the hybrid DA are further diagnosed to evaluate how these differences contribute to the difference in the track forecasts. RMSE profiles of the 24-h forecasts are also displayed when using a set of conventional observations as reference.

It is well known that WRF 3DVAR CV5 BE modeling using the NMC method tends to overestimate the error covariance of winds; and therefore, the spatial correlation scales are excessively large (Lee et al., 2006). In this case, small-scale observed information can be filtered out in the analysis step, and the locally observed information is easily spread out to large spatial distances (Daley, 1991). The impact on typhoon forecasts of tuning the background error has been revealed by several numerical studies using conventional observations (Gu et al., 2005; Guo et al., 2006). However, to the best of our knowledge, there are currently no specific studies published on tuning the BECs to optimize the impact of using AMSU-A radiance observations on the prediction of typhoon tracks. In order to compare with the hybrid DA method as objectively as possible, we conducted a set of sensitivity experiments to investigate the impact of the tuned background-error covariance by tuning the LEN_SCALING parameters defined as the final percentage of the length scale of the background error in use for the five variables (perturbation stream function, velocity potential, temperature, humidity, and surface pressure) and VAR_SCALING parameters as the final percentage of the variance of background error for the same five variables on the prediction of Typhoon Megi (2010). For both LEN_SCALING and VAR_SCALING, 1.0 is the default value in WRF 3DVAR, which means that no reduction is conducted for the length scale or the background error variance.

First, how the physical length scale corresponds to different LEN_SCALING settings is revealed by the results of several single observation tests. A 1 m s^-1 u-wind increment was added in at the 19th model level at 0000 UTC 18 October. Figure 4 shows the temperature and horizontal wind increments obtained in response to the u-wind difference for LEN_SCALING = 1.0, 0.5 and 0.1 with a fixed VAR_SCALING (1.0) in DA. We can see that the radius of influence of the increments was quite large when we set LEN_SCALING to 1.0, whereas it was much smaller for LEN_SCALING = 0.1. The results indicate that tuning of the length scale factor could be necessary to optimize the performance of radiance data assimilation, even with the background error calculated locally using the NMC method.

Figure 4.  Increment\;distributions\;of\;wind(vectors;\;units:\;m\;s^-1) and temperature (unit: K) at the 19th model level for LEN_ SCALING = (a) 1.0 (b) 0.5 and (c) 0.1.

For the case of Typhoon Megi (2010), we set VAR_ SCALING to 1.0, 0.5 and 0.1. For each chosen value of VAR_SCALING, LEN_SCALING was varied from 0.1 to 1.0 in steps of 0.1. In total, there were 30 different combinations of LEN_SCALING and VAR_SCALING. Figure 5 shows the predicted typhoon tracks and their corresponding forecast errors as a function of forecast time from the 30 forecast runs. For most cases, AMSU-A radiance assimilation reduced the forecast track errors compared to the NOAMA experiment. The track from the NOAMA experiment had an east bias that caused the typhoon to land unexpectedly in Taiwan. With the AMSU-A observations assimilated, it seems that with a fixed value of VAR_SCALING, the predicted typhoon tracks shifted gradually westward when we increased the value of LEN_SCALING. This indicates that with larger LEN_SCALING, the NMC method tends to overestimate the length scale, which determines the pattern and the extent to which the observed radiance information spreads out.

In Figs. 5a and b, for VAR_SCALING = 1.0, and when LEN_SCALING was larger than 0.5, the predicted tracks shifted too far west, with a track error even larger than that in the NOAMA experiment. With the smallest value of LEN_SCALING (0.1), the forecast track stayed almost the same as NOAMA. With reduced VAR_SCALING (0.5), predicted tracks with westward bias were still observed from the AMSU-A assimilation experiments for most cases, and the shift was also obvious with an increasing LEN_SCALING, but not as exaggerated as those cases with VAR_SCALING = 1.0 (Figs. 5c and d). The assimilations with moderate values (0.5 and 0.6) for LEN_SCALING produced a smaller forecast error than the other runs. When VAR_SCALING was further reduced to 0.1 (Figs. 5e and f), the predicted tracks for most values of LEN_SCALING were closer to the best track, compared to those in Figs. 5a and b. Nevertheless, the results still showed a westward bias tendency of the track with increasing LEN_SCALING (Fig. 5e). From Fig. 5f, it seems that the forecasts with LEN_SCALING = 0.5 and 0.6 had smaller track error than those from others forecasts. Based on the above results, we conclude that the tuning of LEN_SCALING in the recursive filter process of the static BEC is meaningful to the analysis from data assimilation and subsequent numerical prediction. On the other hand, it seems there is no optimal setting in terms of tuning VAR_SCALING for radiance data assimilation. Accordingly, for all the experiments discussed in the following section, we used VAR_SCALING = 1 and LEN_SCALING = 0.5.

Figure 5.  The (a, c, e) tracks and (b, d, f) forecast track error initialized at 0000 UTC 18 October 2010 using the settings (a, b) VAR_SCALING = 1.0 and LEN_SCALING = 0.1-1.0 in steps of 0.1; (c, d) VAR_SCALING = 0.5 and LEN_SCALING = 0.1-1.0 in steps of 0.1; and (e, f) VAR_SCALING = 0.1 and LEN_SCALING = 0.1-1.0 in steps of 0.1. The gray curve is the control run with no AMSU-A data and the black curve is the best track.

Figure 6.  Ensemble spread for (a) wind speed (units: m s^-1) and (b) temperature (unit: K) valid at 0000 UTC 18 October 2010 at the 19th model level.

Figure 7.  Static background error forecast RMSEs against GFS analyses, and ensemble spread for (a) wind speed (units: m s^-1) and (b) temperature (unit: K) at 0000 UTC 18 October 2010.

Since a high-quality prior ensemble is the key to successful hybrid analyses, it is important to evaluate the ensemble performance. Figure 6 shows the ensemble spread of wind and temperature on the 19th model level after the 6 h forecast valid at 0000 UTC 18 October, when Typhoon Megi (2010) intensified. The 6 h forecast ensemble directly provides the flow-dependent background error in the hybrid system. Since there are more forecast uncertainties where the spread is large, observations in areas with large ensemble spread are most likely to have a large impact in the hybrid system. Likewise, observations are less likely to influence the analyses in areas where the spread is smaller. The patterns reflecting the features of the observation locations and meteorological conditions can be seen in Figs. 6a and b. The spread was large over the Tibetan Plateau in the west, where few observations were available to constrain the model. Conversely, spread was smallest over Eastern China, where observations were plentiful. A local maximum of spread was obvious for wind and potential temperature in the northeast of the Philippines, where Typhoon Megi (2010) moved, reflecting the forecast uncertainty for a TC. The ensemble spread also suggests larger forecast uncertainties around Megi (2010) than in its environment.

The ensemble spread can represent the ensemble mean forecast error in a well-calibrated system (Houtekamer et al., 2005). The forecast RMSEs, ensemble spread, and the static background error from the WRFDA (Wang et al., 2014) using the NMC method and the WRFDA modeling are shown in Figs. 7a and b. The forecast RMSEs were obtained by comparing the forecast ensemble mean with the GFS analyses. The static background error from the WRFDA modeling was estimated based on the ensemble perturbations using Gaussian random draws with the static BECs and zero mean. The ensemble mean RMSEs of wind were less than 3 m s^-1 and the temperature was less than 1 K at most levels. For winds, the static background error from the WRFDA modeling was largely underestimated, consistent with the results in (Wang et al., 2014); whereas, the ensemble spread was much closer to the RMSEs compared to the static background error. For temperature, both the ensemble spread and the static background error were greater than the corresponding RMSEs between model level 34 and model level 41. For other levels, the ensemble spread was much closer to the RMSEs, while the static background error underestimated the forecast errors. Overall, the ensembles were reasonably well calibrated. The final background error from the hybrid system, as a mix of the flow-dependent and the static background error, plays an important role in the data assimilation procedure.

Figure 8.  The predicted (a) tracks and (b) track errors from 0000 UTC 18 to 0000 UTC 23 October 2010, along with the best track.

The predicted typhoon tracks and track errors from NOAMA, 3DVAR_AMA, HYBRID_AMAS and HYBRID_ AMAF are shown in Figs. 8a and 8b for the 120 h forecast from 0000 UTC 18 October through 0000 UTC 23 October 2010. The best track is also shown in Fig. 8a. The center of the typhoon is defined as the location of MSLP. All experiments had similar track forecasts for the first 24 hours for model integration. In NOAMA, the predicted typhoon tracks had a substantially large eastward bias compared to the best track, with an average track error close to 238 km within the 120 h. In fact, such a strong bias has been observed by several other global forecasting models initialized in a similar time period. This phenomenon could be partly attributed to the abnormal heavy northwesterly flows associated with the trough over mainland China (Kieu et al., 2012). The forecasts with radiance assimilation produced a much improved position and subsequent westward movement that occurred between 24 h and 48 h. Of note is that the tracks from all the radiance DA experiments were rather closer, with a track error of less than 100 km during the first 84 h forecast. After that, the vortices of 3DVAR_AMA and HYBRID_AMAS tended to move too quickly to the northwest and the track errors increased with forecast time to almost 300 km; whereas, the track forecast from HYBRID_AMAF agreed much better with the best track. Note that the tracks of 3DVAR_AMA and HYBRID_AMAS were very close, with mean track errors of 110 km and 106 km, respectively. These results suggest that the improvement in HYBRID_AMAF compared to 3DVAR can be mainly attributed to the use of the flow-dependent ensemble covariance, more so than the use of the ensemble mean as the background forecast. HYBRID_AMAF had the smallest mean track error (48 km), especially for forecast lead times longer than 24-h with the recurving track predicted. In general, the forecast track from the HYBRID_AMAF experiment aligned more closely with the best track than those from NOAMA, 3DVAR_AMA and HYBRID_AMAS.

Next, the analysis differences among the DA methods are examined to identify how these differences affect the subsequent forecasts. Figure 9 shows the analysis increment (shaded) and the background (contours) for the 500 hPa geopotential height at 0000 UTC 18 October. The analysis increment from 3DVAR_AMA and HYBRID_AMAF were different for both Megi (2010) (the fields around the typhoon) and its environment (the large-scale fields). In HYBRID_ AMAF, the analysis increment for Megi (2010) was larger than that for its environment (Fig. 9b), while in 3DVAR_ AMA, the analysis increment for Megi (2010) was comparable or less than the increment for the environment (Fig. 9a). In 3DVAR, there was an increase in the geopotential height to the northeast of the Philippines and a decrease in the geopotential height in southern China. The increment for Megi (2010) from HYBRID_AMAF was much stronger compared to that from 3DVAR_AMA. The increment in HYBRID_AMAF showed an obvious dipole pattern, which could push Megi (2010) to move to the northeast of the background forecasted position. This suggested that HYBRID_AMAF was able to adjust the position of Megi (2010) in the background forecast after the DA with the typhoon-specific increment, which may have contributed to the better prediction of the subsequent recurving track. The results indicate that HYBRID_AMAF can systematically correct the position of the typhoon by the DA with a flow-dependent ensemble covariance, even without other techniques such as the extra vortex relocation procedure, which is widely employed by the operational 3DVAR system.

Figure 9.  The 500 hPa geopotential height increments (color shading) from (a) 3DVAR_AMA and (b) HYBRID_AMAF valid at 0000 UTC 18 October 2010 for Typhoon Megi (2010). The black contours are the corresponding background 500 hPa geopotential height fields valid at 0000 UTC 18 October 2010.

Figure 10.  Vertical profiles of 24 h forecast RMSEs for (a) u-wind (units: m s^-1), (b) v-wind (units: m s^-1), (c) temperature (unit: K), and (d) specific humidity (units: g kg^-1), against the conventional observations. The numbers of conventional observations are shown on the right of each panel.

The performances of forecasts with 24 h lead times from the different DA experiments are evaluated in this sub-section. The RMSE profiles of horizontal wind, temperature, and specific humidity were calculated based on the difference of the model forecasts initialized from 0000 UTC 1800 to 0000 UTC 2300 and a set of conventional observations (Fig. 10); for example, radiosondes and the atmospheric motion vector winds from geostationary satellites (GeoAMV). The largest RMSE of horizontal wind and temperature in each experiment was observed near the upper troposphere, while the largest RMSE of specific humidity was around the lower troposphere. HYBRID_AMAS had slightly smaller RMSEs than 3DVAR_AMA for almost all the variables at most levels, except for moisture near the surface. HYBRID_AMAF substantially outperformed 3DVAR_AMA in terms of the horizontal winds and temperature in the lower and mid-troposphere and for moisture at all levels, though we tuned the static BEC to optimize 3DVAR_AMA's performance. HYBRID_AMAF had smaller RMSEs than HYBRID_AMAS for most variables at most times. The improvement for specific humidity can be attributed to the multivariate correlations between moisture variables and other variables from the ensemble forecasts, which did not exist in the static BEC in 3DVAR_AMA.

In this study, several experiments involving the assimilation of AMSU-A radiance data for Super Typhoon Megi (2010) were performed to investigate the effect of the hybrid DA approach on TC track forecasts. Instead of using an EnKF to update the ensemble, as in earlier studies, in this study each ensemble member was updated by running the hybrid analysis system multiple times with perturbed observations. Model predictions of Megi (2010)'s track were verified against the best-track analyses of the CMA. Further diagnostics were conducted to evaluate how the analyses increment differences contributed to the subsequent forecast differences. Aspects of the flow-dependent error covariance represented by the ensemble spread from the ensemble forecasts were examined. The 24 h forecasts were also verified against a set of conventional observations. The main findings from the experiments can be summarized as follows:

(1) The tuning of length scale factor is a necessary task to optimize the performance of radiance data assimilation, even with the background error calculated locally using the NMC method. In the case of Typhoon Megi (2010), a moderate length scale factor with roughly 50% of the default value improved the prediction of the typhoon track, but how the forecast error depended on the variance scale factor was more mixed.

(2) AMSU-A radiance data assimilation had a slight positive impact on the track forecast from 3DVAR, even when the static BEC was carefully tuned to optimize its performance. When the hybrid DA was employed, the track forecast was significantly improved, especially with respect to the sharp recurving from west-northwest to north after crossing the Philippines. In general, the forecast track from the HYBRID_AMAF experiment aligned more closely with the best track than NOAMA, 3DVAR_AMA and HYBRID_AMAS Megi (2010).

(3) The increments for both the fields around the typhoon and the large-scale environment produced by the hybrid DA and 3DVAR were different. The hybrid DA is capable of systematically adjusting the position of the typhoon with the typhoon-specific BEC used. (4) For 24 h forecasts, the hybrid DA experiment with the use of full flow-dependent background error substantially outperformed 3DVAR in terms of the horizontal winds and temperature in the lower and mid-troposphere and for moisture at all levels, though we tuned the static BEC to optimize its performance.

The findings of this study are encouraging and suggest that hybrid DA can initialize higher quality forecasts than the 3DVAR system. Given that any improvements in the typhoon intensity forecast were not obvious due to the limitations of NWP modeling related to aspects such as dynamics, physical parameterizations, and spatial resolution, only track forecasts are emphasized in this paper. Further investigations into the use of the hybrid radiance data assimilation system for improving typhoon intensity forecasts are ongoing. In this study, we used just one outer loop during the variational minimization process. In the future, we intend to focus on employing multiple outer loops in 3DVAR frameworks and retuning the scale length for the 3DVAR recursive filter to assess the relative contributions of the static and ensemble BECs. To understand the relative advantages and disadvantages of different techniques in typhoon forecasts, direct and thorough comparisons with other DA techniques such as the EnKF, 4DVAR, and En-4DVAR are also planned.