The impact of ocean data assimilation on seasonal predictions based on the National Climate Center climate system model

1. Introduction

Oceans play a key role in the predictability of the climate system due to their tremendous thermal inertia compared to atmosphere or land (Counillon et al., 2014). Accuracy of the ocean initialization during modeling can significantly impact seasonal to decadal climate predictions (Alves et al., 2011; Zheng and Zhu, 2015). A common strategy to obtain the optimal initialization is to assimilate available ocean observations into ocean models, which aim to produce the best estimates of ocean states.

There have been many advances in data assimilation techniques ranging from the relatively simple optimum interpolation (OI) and three-dimensional variational methods (3DVAR) to more sophisticated four-dimensional variational methods (4DVAR) and the ensemble Kalman filter (EnKF) approach. The OI and 3DVAR based schemes are computationally cheap to perform and have been widely used in operational ocean forecasting systems. However, both OI and 3DVAR use the time-invariant background error covariance, which tends to produce inaccurate analyses in areas with highly nonlinear flows. This problem can be partly solved by using the flow-dependent error covariance adopted in EnKF and 4DVAR.

Although EnKF and 4DVAR have been used in many studies, practical problems still exist for realistic ocean applications, especially for operational global ocean data assimilation systems. One disadvantage is that EnKF and 4DVAR are computationally expensive to perform. For example, the computational costs of EnKF increase linearly with the ensemble number N. A value of N>20 is unaffordable for operational forecasting given the current limited computer resources, while EnKF usually requires more than 20 ensemble members (e.g., Miyazawa et al., 2012; Xu et al., 2013; Xu and Oey, 2014).

A computationally inexpensive ensemble optimal interpolation (EnOI) approach is adopted in this study. EnOI runs only a single model member every time and has no risk of ensemble collapse (Pan et al., 2014). Its analysis formula is identical to that of the local ensemble transform Kalman filter (LETKF, refer to Miyazawa et al., 2012 for details), except that its background error covariance is advanced from a prescribed 100 static ensemble members instead of a flow-dependent ensemble. In general, EnOI has many attractive characteristics such as multivariate assimilation and inhomogeneous and anisotropic covariance. In addition, the static ensembles for EnOI can be time-dependent (e.g., Oke et al., 2005, 2013; Fu et al., 2008) or seasonally varying. Consequently, EnOI has been used in many operational ocean forecast systems such as BODAS (bluelink ocean data assimilation system) at the Bureau of Meteorology in Australia (Oke et al., 2013).

A new generation of climate forecast system at the Beijing Climate Center is under development (Beijing Climate Center Climate System Model, BCC_CSM1.1) (e.g., Wu et al., 2010; Wu et al., 2014). BCC_CSM1.1 is a fully coupled climate system consisting of atmosphere, land, ocean, and sea ice components. The primary objective in regard to developing BCC_CSM1.1 is to generate a high-quality reanalysis dataset and improve predictions from sub-seasonal, seasonal, and up to decadal time scales. The development of a data assimilation system is crucial for this objective. One purpose of this study is to introduce the new ocean data assimilation system that is going to be adopted in the BCC_CSM1.1.

Observing system simulation experiments (OSSEs) are widely used to examine the impacts of current observations and the performance of data assimilation systems (Masutani et al., 2010). OSSEs are powerful tools for providing guidelines for the development of forecast model and assimilation systems. They usually incorporate the following parts (Masutani et al., 2010): (1) a nature run which provides the synthetic “truth”; (2) a generation of simulated observations including estimated errors; (3) a control run without assimilating observations; (4) a perturbations run assimilating observations; (5) a comparison on the forecast skills between the control run and the perturbation run.

The other purpose of this study is to investigate the impact of data assimilation of various available observations on seasonal forecasts using the OSSEs framework. For this purpose, the individual and combined contributions of sea surface satellite data to forecasting, such as sea surface temperature (SST), sea surface height (SSH), sea surface salinity (SSS), and temperature and salinity (T/S) profiles were evaluated. Model generated SST, SSH, and SSS were taken as pseudo-observations of satellites, and T/S profiles close to locations of Argo floats were chosen to represent pseudo-observations of Argo. The satellite sea surface data and Argo float data are major observational data sources nowadays, with global coverage and widespread availability in most of the ocean observing network. The satellite SST observations have been widely used in ocean assimilation applications since SST is a key geophysical variable in air–sea exchanges of heat (e.g., Tang et al., 2004). The SSS plays an important role in surface mixed layer dynamics, water mass formation, and global ocean circulation (Vernieres et al., 2014). The satellite observations of SSS have been available since the first satellite was launched by the European Space Agency to monitor SSS (Boutin et al., 2016). Global SSH data from TOPEX/Poseidon altimeters have been available since October 1992. The dynamic topography depicts the surface geostrophic flow field. Furthermore, large-scale variability of SSH has close connections with climate signals. For example, assimilation of SSH contributes to better understanding of the tropical Pacific variability (Carton and Giese, 2008) and El Niño forecasting (Ji et al., 2000). One major concern in assimilating SSH is how to project the surface information downward to subsurface quantities (Fu et al., 2011). At present, SSH data are assimilated into ocean models either by developing a statistical relationship between SSH and subsurface temperature/salinity (Behringer et al., 1998; Yan et al., 2004) or by the inherent multivariate relation derived from the ensembles by using some ensemble-based data assimilation methods (Oke et al., 2008; Zheng and Zhu, 2015). T/S profiles improve the representation of seawater density, which dictates water mass. In addition, T profiles have a direct influence on ocean heat content.

EnOI is implemented in a global ocean model (about 110 km in the horizontal) based on MOM4.0, which is the ocean model used in BCC_CSM1.1. An idealized twin experiment was carried out to test the assimilation and prediction system in a situation where the “truth” was known. The “observed” SST, SSH, SSS, and T/S were derived from free mode simulations and considered as the “truth.” This paper is organized as follows. Section 2 presents a brief introduction of the EnOI data assimilation system and experimental setup. The assimilation and forecast results of all experiments are presented in Section 3. The discussion and conclusions are given in Section 4.

2. The EnOI assimilation system and twin experiment setup

2.1 EnOI

EnOI is a simplified form of EnKF, which uses a stationary historical ensemble of model states to represent the background covariance matrix instead of time-dependent ensembles for EnKF. Consequently, it is more computationally efficient than EnKF, but is still multivariate and three-dimensional. In this study, EnOI is derived based on LETKF, an advanced version of EnKF (Miyoshi et al., 2010).

Here, the calculation equations of EnOI are given below (Oke et al., 2013).

where $ {{{\varphi }}^{\rm{a}}} $ is an m-dimensional vector representing the model analysis, $ {{{\varphi }}^{\rm{f}}} $ is an m-dimensional vector representing the model forecast, and ρ is a scaling factor used to represent the instantaneous forecast error variance, which is usually less than the historical error variance over a long time period. ρ is in the range between 0 and 1, and it was set to 0.5 here by tuning the assimilation results. A is the historical ensemble composed of model states, and $ {{A}}' $ is the centered historical ensemble (i.e., ${{A}}'={{A}}-\overline {{{A}}}$). $\overline {{A}}=\dfrac{1}{N}{\displaystyle\sum \limits_{i=1}^{N}}{{{{A}}}_{i}}$, and N represents the number of the historical ensembles. A is an m × N matrix. w^a is an N-dimensional vector calculated from the observational data, model forecast, and historical ensemble model simulations; it can be computed as follows:

where d represents the measurements, H is the measurement operator that interpolates the model space into the observational space, and R is the measurement error covariance.

Localized use of observation data is important in the method. The primary benefit of localization is to increase the rank of the forecast covariance, thus resulting in analysis fields that fit well with the observations (Oke et al., 2007). Localization is implemented explicitly in consideration of observational data from a region surrounding the target model grids. Two localized scale parameters are defined following Miyoshi et al. (2010) and Miyazawa et al. (2012):

where $ {\sigma _{{\rm{obs}}}}$ (the number of surrounding grids) and $ {\sigma _{{\rm{obsv}}}}$ (meters) are the horizontal and vertical localization scales, respectively. The localization scale is chosen to correspond to the distance at which the Gaussian function drops to e^–0.5 (Miyoshi et al., 2010). Observational data far from the target grid with horizontal distances larger than $ Dis{t_{{\rm{zero}}}}$ or vertical distances larger than $ Dis{t_{{\rm{zerov}}}}$ are not used. A factor, exp $\left\{0.5 \times \left[ {{{\left( {\dfrac{{Dist}}{{Dis{t_{{\rm{obs}}}}}}} \right)}^2}{\rm{ + }}{{\left( {\dfrac{{Dist}}{{Dis{t_{{\rm{obsv}}}}}}} \right)}^2}} \right]\right\}$, is multiplied to enhance observational errors of data far from the target grid (Miyazawa et al., 2012). The resulting localization scales are approximately 110 km and 2 000 m in the horizontal and vertical, respectively.

2.2 The global ocean model

The EnOI algorithm are implemented in MOM4, which was originally developed at the Geophysical Fluid Dynamics Laboratory (Griffies et al., 2003). The model covers the global ocean with a horizontal resolution of 1° and at 50 vertical levels. In the meridional direction, the resolution increases to (1/3)° within 10° of the equator, and it smoothly reduces down to 1° poleward of 30°. To avoid a singularity at the North Pole, tripolar grids are adopted (Griffies et al., 2005). The physical parameterization schemes used in the simulation include the K-profile parameterization vertical mixing scheme, the isopycnal tracer mixing and diffusion, and the Laplace horizontal friction scheme, etc., the same as described in Griffies et al. (2005).

The model is driven by wind stress and heat fluxes estimated from 6-hourly atmospheric variables obtained from the National Centers for Environmental Prediction/National Center for Atmospheric Research Reanalysis I dataset (NCEP/NCAR, http://www.esrl.noaa.gov/psd/). The climatological river runoff (http://www.cgd.ucar.edu/cas/catalog/) is specified at the model coastlines. The surface temperature and salinity are relaxed to World Ocean Atlas (WOA09) monthly climatology (http://coastwatch.pfeg.noaa.gov/erddap/griddap/nodcWoa09mon1t.html), with restoring time scales of 90 d and 120 d, respectively. Tidal forcing is not included. Sea ice is simulated with the Sea Ice Simulator (Griffies et al., 2011).

In this study, the model was first spun up from 1948 to 2000, and a statistically quasi-equilibrium ocean field was established. This run was then continued from January 1, 1990 through 2009.

2.3 Twin experiment setup

The performance of EnOI assimilation system of BCC_CSM under initial condition errors was tested by identical twin OSSEs, in which model errors due to parameterization and numerical difference implementations were neglected (Zhang et al., 2015). An experiment (denoted as TRU), which was allowed to freely run from January 1, 2005 to December 31, 2006, was designed to produce pseudo-observations and make comparisons with the assimilative analysis and model predictions. Another free model run (denoted as CTR), which was the same as TRU but initialized on a start date of June 1, 1990, was used to create large biases for initial conditions. In the identical twin OSSEs, the difference between experiments results is small when integrating at different periods. A time-invariant ensemble was used to approximate the system’s background error covariance matrix in the EnOI. The output from a historical run of BCC_CSM was used to construct the ensemble of seasonal anomalies. One hundred ensemble members for estimates of the background error covariance were sampled from the free-run simulation at time intervals of 25 days from January 1, 1995 to December 31, 2009. The relationships (e.g., T-S) between different control variables were derived from the selected 100 ensembles.

The sea surface pseudo-observations including SST, SSH, and SSS were selected every three points in the model grids meridionally and zonally. The pseudo T/S profiles were selected at model grids that were as close to the locations of Argo floats on June 1, 2005 as possible (Fig. 1). This time was chosen because it was the median time of the assimilation period (January 1, 2005 to December 31, 2005). Considering the slow drifts of most Argo floats, the locations of pseudo T/S profiles did not change with time because of the relatively low model resolution of about 1°×1°. The vertical levels are set as the same as the model levels. The altimeter SSH errors generally vary from 1 cm to 4 cm (Chambers et al., 2003), and thus, the pseudo SSH error was specified as 3 cm. The SST error was set to be 0.3°C according to Guan and Kawamura (2004). The SSS error was set to be 0.1 in consideration of the rapid development in inversion algorithms for satellite salinity (Peng et al., 2016). Similarly, for T/S profiles, the temperature and salinity errors were prescribed to be the same as the SST error and the SSS error, respectively.

Figure  1.  Locations of Argo floats on June 1, 2005, which were used to produce pseudo-observations of T/S profiles.

DownLoad: Full-Size Img PowerPoint

Assimilation experiments, E01–E07 initialized as CTR, assimilated pseudo-observations from January 1 to December 31, 2005 (Table 1). E01 assimilated SST only, E02 assimilated SSH only, E03 assimilated both SST and SSH, E04 assimilated SSS, E05 assimilated all the SST, SSH, and SSS data, E06 assimilated T/S, and E07 assimilated all the variables. Then, seven 12-month test forecasts are conducted starting from January 1, 2006 corresponding to the seven analyses. By comparing the model states from CTR, E01–E07 against the known “true” states, it is convenient to investigate the performance of the assimilation system and forecast skills.

Table  1.  Experimental design

	SST	SSH	SSS	T/S profiles	Assimilation /month	Free run /month
CTR	–	–	–	–	1–12	13–24
E01	Yes	–	–	–	1–12	13–24
E02	–	Yes	–	–	1–12	13–24
E03	Yes	Yes	–	–	1–12	13–24
E04	–	–	Yes	–	1–12	13–24
E05	Yes	Yes	Yes	–	1–12	13–24
E06	–	–	–	Yes	1–12	13–24
E07	Yes	Yes	Yes	Yes	1–12	13–24
Note: – represents no data.

 | Show Table

DownLoad: CSV

3. Results

3.1 Assimilation performance measures

To evaluate the performance of data assimilation experiments, the domain-averaged root mean square error (RMSE) of SST, SSH, SSS, temperature, and salinity in the upper ocean (0–500 m) and the deep ocean (500–1 500 m) with respect to the pseudo observations of TRU experiment from Months 1 to 12 are computed (Fig. 2). All assimilative experiments generally showed improvements over CTR, but the improvements varied among different experiments. The SST RMSEs of E02 and E04 were comparable with that of CTR, while the other assimilation experiments approximately reduced the RMSEs by half (Fig. 2a). These results indicate that assimilation of SSH and SSS alone do not contribute to the SST analysis in the system. In EnOI, the projection of SSH and SST into subsurface is achieved through the multivariate relationship derived from the static ensembles. Assimilating surface observations of one type will update all control variables in the assimilation step. The SSS RMSEs of E04, E05, and E07 were reduced by about 80% compared to CTR (Fig. 2b). The RMSEs of E03 and E06 were reduced by about 30%. E01 and E02 only slightly improved SSS estimates. So experiments with SSS assimilation improve the SSS the most, while assimilating T/S profiles alone (E06) or SSH and SST (E03) can only improve SSS to a limited extent. For SSH, E02, E03, E05, and E07, the RMSE was reduced by about 82%, thus indicating the importance of the SSH assimilation (Fig. 2c). The T/S profile assimilation (E06) reduced the SSH RMSE by about half. SST (E01) and SSS (E04) only improved the SSH by about 20%. For the analysis of temperature and salinity at depth, all experiments showed improvements (Figs 2d–g). E07 had the smallest RMSEs among all experiments. The RMSEs obtained when assimilating SST alone (E01) were larger than those obtained when assimilating T/S profiles (E06). Similarly, the RMSEs obtained when assimilating SSS alone (E04) were larger than those obtained when assimilating T/S profiles (E06). These results demonstrate the importance of the assimilation of T/S profiles in the global data assimilation system.

Figure  2.  Time series of domain-averaged RMSEs for SST (a), SSS (b), SSH (c), temperature in the upper 500 m (d), temperature from 500 m to 1 500 m (e), salinity in the upper 500 m (f), and salinity from 500 m to 1 500 m (g) from all the experiments (CTR and E01–E07) over the assimilation period (Months 1–12).

DownLoad: Full-Size Img PowerPoint

3.2 Predictions

The impacts of data assimilation on seasonal forecasts were investigated by conducting a 12-month forecast initialized from restart files produced by CTR, E01–E07. The forcing was identical for all cases.

The time series of spatial RMSEs for temperature and salinity among all forecast experiments and TRU are shown in Fig. 3. All RMSEs from the forecasts initialized from the assimilated runs were smaller than that from CTR. The E07 forecast, initialized from joint assimilation of SST, SSS, SSH, and T/S profiles, had the smallest RMSEs for temperature and salinity compared to the others. Figure 3 also shows the “persistence” curves (black lines) based on TRU, i.e., the temperature, salinity, and SSH from TRU (January 1 to December 31, 2005) were assumed as “repeat” for the subsequent 12 months. The model forecasts from E03, E05, and E07 beat the persistence in the upper ocean (Figs 3a–d and f), while the other experiments showed some deficiencies, such as E01 for the SSS forecast (Fig. 3b) and E04 for the SSH forecast (Fig. 3c). In contrast, in deep water (500–1 500 m), the persistence beat the model forecasts because of the large bias from the initialization starting on June 1, 2005 in the CTR run and all assimilation runs (Figs 3e and g). These results demonstrate that the deep ocean bias cannot be completely corrected after one-year assimilations, though improvements are possible.

Figure  3.  Time series of domain-averaged RMSEs for SST (a), SSS (b), SSH (c), temperature in the upper 500 m (d), temperature from 500 m to 1 500 m (e), salinity in the upper 500 m (f), and salinity from 500 m to 1 500 m (g) from all the experiments (CTR and E01–E07) over the prediction period (Months 13–24).

DownLoad: Full-Size Img PowerPoint

In addition to the aforementioned temporal variability, the spatial variability of ocean predictions was evaluated in different experiments as well. The SST and surface currents were compared at first. Ocean heat content (OHC) and mixed layer depth (MLD) were used to examine the subsurface predictions. Besides, the values of the Niño 3.4 index were compared, too, because it is an important climate signal in the tropical Pacific.

3.2.1   SST

Figure 4 shows the spatial distribution of the RMSE for SST from CTR and the differences with respect to E01–E07 over the prediction period (Months 13–24). Reductions of the RMSE were generally found in all experiments in regions where the RMSE was relatively high in CTR, such as the tropical eastern Pacific, the subarctic, and the Southern Ocean. The domain-averaged RMSE from E07, which was about 0.11°C, was the smallest. Interestingly, the experiment initialized from the assimilation of T/S profiles (E06) produced a comparable RMSE to the one from the experiment initialized from the joint assimilation of all surface data (E05). The values of the domain-averaged RMSEs from these two experiments were much smaller than those of E01, E02, and E04. These results indicate the importance of multivariate assimilation on initial conditions, and subsequently the forecasting.

Figure  4.  Global distribution of the RMSE for SST from CTR, and the bias of RMSE between assimilation run (E01–E07) and CTR from Months 13 to 24. Negative values mean that the RMSE from the assimilation runs is smaller than that from CTR. The black box indicates an area with enhanced RMSE. The domain-averaged RMSE is shown on the top left of each panel in degrees Celsius.

DownLoad: Full-Size Img PowerPoint

Noticeably, in the Southern Ocean south of Africa (48°–60°S, 0°–60°E, the black box in Fig. 4), the RMSEs for SST were much larger in E01, E02, and E03 than that in CTR. To explore the reasons for the high RMSEs, the time evolution of vertical profiles of temperature averaged over the high RMSE region is examined. Figure 5 shows the vertical profiles of temperature and the corresponding RMSEs in January and September from all experiments. In January, SST estimates in E01, E03, E05, and E07 were much better than that in CTR (Figs 5a and c). Conversely, in the subsurface layer (50–150 m), the values of RMSE in E01 and E03 were about 0.3°C larger than that in the CTR. In September, a thick mixed layer about 100 m deep developed after the austral winter (Figs 5b and d). A pre-existing subsurface bias of temperature in E01, E03, and E05 emerged near the surface due to the winter mixing. Since temperatures in both the surface and subsurface in E02 and E04 were more biased than those in CTR, the values of RMSE remained relative high. On the other hand, the temperatures in deep ocean waters (below about 200 m) were improved in all experiments compared to the CTR. Thus, the high RMSE found in the black box of Fig. 4 mainly developed as a result of the large temperature bias in the subsurface. The subsurface bias probably came from inaccurate estimates of the background error covariance in the multi-fronts Antarctic circumpolar region during the surface data assimilation.

Figure  5.  Vertical profiles of the box-averaged temperature (°C) in all experiments (a, b) and the RMSE for temperature (c, d) in January and September (in 48°–60°S, 0°–60°E) derived from eight experiments along with the TRU experiment.

DownLoad: Full-Size Img PowerPoint

3.2.2   Ocean surface currents

Large-scale ocean circulation is primarily geostrophic a few degrees away from the equator. Because of the availability of long-term satellite altimeter data, geostrophic parts of any model generated currents can be easily evaluated by using altimetry data. It is hence interesting to assess all forecasting experiments in this regard. Figure 6 compares the RMSEs of the predicted SSH from all experiments. Compared to CTR, significant improvements of SSH were found in E02, E03, E05, E06, and E07. Similar to the SST RMSE reduction, the large reduction primarily occurred in regions where the RMSE was large in CTR. However, E01 and E04 showed almost no reduction, thus indicating that assimilation of SST or SSS alone cannot largely improve SSH forecasting.

Figure  6.  Global distribution of the RMSE for SSH from CTR, and the bias of RMSE between assimilation run (E01–E07) and CTR from Months 13 to 24. Negative values mean the assimilation run is better than that from CTR.

DownLoad: Full-Size Img PowerPoint

The RMSEs of the predicted surface current speed from all experiments are compared in Fig. 7. The largest reduction of the overall RMSE was observed in E07. E02, E03, and E05 showed clear improvements near the tropical Pacific and the western boundary current regions such as the Gulf Stream, the Kuroshio extension. Comparisons of E01–E05 revealed that SSH assimilation was the dominant factor accounting for the forecast improvements. In contrast, major improvements in E06 were observed in the tropical regions, and these were more prominent than those in E05. These results indicate that the ocean surface currents in the tropics are better predicted when initial conditions are produced by assimilating T/S profiles, while surface data assimilation becomes more important at higher latitudes, particularly near the western boundary currents. This can be attributed to the dominant effects of SSH assimilation on geostrophic parts of surface currents away from the equator.

Figure  7.  Global distribution of the RMSE for surface current speed from CTR, and the bias of RMSE between assimilation run (E01–E07) and CTR fromMonths 13 to 24. Negative values mean that the RMSE from the assimilation runs is smaller than that from CTR.The domain-averaged RMSE is shown on the top left of each panel in m/s.

DownLoad: Full-Size Img PowerPoint

3.2.3   OHC

Ocean heat content is an important variable in climate studies, and it reflects the internal energy that the ocean has. To assess the standalone and joint effects of assimilation of surface data and T/S profiles on ocean predictions, the upper 700 m OHC estimates from all experiments are analyzed. Figure 8 shows the global distribution of the time-averaged upper 700 m OHC per unit area from all forecast experiments relative to that from the “truth”. E07 had the smallest RMSE for OHC compared to TRU (Fig. 8h). The RMSE in E06 was about 1.2×10⁸ J/m² larger than that in E05, and this was primarily caused by the large bias in the subpolar regions (Figs 8f and g), where the T/S profiles were relatively sparse compared to the gridded satellite surface data. In the lower latitudes, the difference between E06 and TRU was smaller than that in E05. Interestingly, none of the standalone assimilation of surface data experiments (E01, E02, or E04) significantly improved the OHC estimates (Figs 8b, c and e). The joint assimilation of SST and SSH had already reduced the deficiency of OHC predictions to a large extent (Fig. 8d). When SSS was assimilated, the reduction was more significant (Fig. 8f). Thus, both surface variables and T/S profiles are important for OHC predictions.

DownLoad: Full-Size Img PowerPoint

Figure  8.  Global distribution of the time-averaged upper 700 m ocean heat content (OHC) per unit area from the “truth” experiment from Months 13 to 24 (i) and the differences of the mean values between the forecast experiments and the “truth” (a–h). The domain-averaged RMSE of OHC per unit area is shown on the top left of each panel in J/m².

DownLoad: Full-Size Img PowerPoint

3.2.4   MLD

Mixed layer depth is one of the most important quantities in the upper ocean. Here, MLD is defined following Breugem et al. (2008) as the depth (z) at which the potential density is $ {\sigma }_{z}=\sigma ({S}_{10\;{\rm{m}}},{\theta }_{10\;{\rm{m}}}\!-\!0.2) $, where θ_{10 m} and S_{10 m} are the potential temperature and salinity at a depth of 10 m, respectively, and σ is the potential density. The MLD in the tropical Pacific is examined since its variability is closely related to the El Niño Southern Oscillation (ENSO). Figure 9 shows the spatial distribution of differences of the time-averaged MLD (over Months 13–24) from all experiments relative to TRU in the tropical Pacific. A weak Central Pacific (CP) El Niño appeared in TRU. Compared to CTR, improvements in MLD were seen in all seven experiments, but the improvements were most significant in E06 and E07 (Figs 9g and h). Noticeably, the differences of MLD from E02 were much smaller than those from E01, which suggests that assimilation of SSH instead of SST is important for improving MLD forecasting. E06 performed better than E05 because the assimilated T/S profiles had a direct influence on the MLD of initial conditions used for forecasting (Figs 9f and g). Since the surface variables are influenced by intense air–sea interactions, there are thus more uncertainties when calculating the background error covariance. Consequently, the deficiency of initial conditions resulting from the assimilation of surface data is then inherited in the forecast.

Figure  9.  The time-averaged mixed layer depth (MLD) in the tropical Pacific from TRU (i) and the differences between all prediction experiments and TRU from Months 13 to 24 (a–h). The domain-averaged RMSE of MLD is shown on the top left of each panel in m.

DownLoad: Full-Size Img PowerPoint

3.2.5   Niño 3.4

Figure 10 compares the Niño 3.4 index from all experiments. Over the assimilation period, the values of the index from all experiments including CTR gradually approached towards the “truth”, though relatively large deficiencies still existed in some experiments such as CTR and E02 (Fig. 10a). It is interesting to note that E04, which involved assimilating SSS, significantly improved the Niño 3.4 data. Over the prediction period, the values of the index started to diversify after four months (Fig. 10b). CTR and E04 did not produce a large positive Niño 3.4 index at the end of the prediction period, thus suggesting that standalone assimilation of SSS cannot well capture an El Niño event one year in advance. Conversely, E02 produced a too strong Niño 3.4 index, thus implying that standalone assimilation of SSH might overestimate an incoming El Niño event. E03, E05, and E07 produced the best estimates of Niño 3.4 index. These findings tell us that the joint assimilation of surface variables and T/S profiles are crucial for ENSO predictions. The newly developed system can well predict an El Niño event one year ahead.

Figure  10.  Comparison of the estimated Niño 3.4 index from all experiments from Months 1 to 12 (a) and from Months 13 to 24 (b).

DownLoad: Full-Size Img PowerPoint

Equatorial long wave propagations are critical to ENSO development. According to Boulanger and Menkes (1995), the coefficient amplitude of the SSH and surface zonal current represents the projection of the ocean variations onto these equatorial long waves (Kelvin waves and Rossby waves). In this study, the prediction of wave propagation is examined at a given longitude and time. Figure 11 shows the coefficients for the long waves, which were computed from SSH and surface current anomalies of TRU, and the differences from CTR and E01–07 relative to TRU. Eastward propagation of Kelvin waves and westward propagation of Rossby waves clearly appeared (Fig. 11i), and subsequently contributed to the occurrence of the CP El Niño (Fig. 9i). E07 produced the best estimates of the wave propagation (Fig. 11h). The second best estimate came from E05 (Fig. 11f). Comparing E03 and E05, the findings suggest that assimilation of SSS can improve the prediction of wave propagation. Besides, the prediction initialized from joint assimilation of surface data (E05) slightly outperformed that from assimilation of T/S profiles (E06) in terms of ENSO predictions.

Figure  11.  Coefficients of the equatorial waves for TRU and the difference between all experiments including CTR (a), E01–E07 (b–h), and TRU (i). The coefficients are nondimensional and computed from the surface zonal current and the sea surface height anomalies.

DownLoad: Full-Size Img PowerPoint

4. Discussion and conclusions

In this study, EnOI with a global ocean model (MOM4.0) was applied to estimate three-dimensional global ocean states when assimilating various variables, e.g., SST, SSH, SSS, and T/S profiles, in an idealized twin experiment framework. Tests were conducted step-by-step to explore the sensitivity of estimates to each variable. The results were compared against quantities from the TRU experiment to assess the analysis and forecasting skills. The major findings are as follows.

(1) Data assimilation generally improved all investigated quantities; assimilation of all the variables together gave more improved results than assimilating them separately.

(2) A 12-month test forecast showed that initializations from E07 produced significantly improved forecasts compared to the others.

(3) The SST forecasts initialized from joint assimilation of surface data (E03 and E05) produced comparable global averaged RMSEs to that from assimilation of T/S profiles (E06), but the assimilation of T/S profiles should not be overlooked because subsurface deficiencies can develop into the surface during forecasts, particularly for highly nonlinear flow regions.

(4) The ocean surface currents in the tropics were better predicted when initial conditions were produced by assimilating T/S profiles (E06), while surface data assimilation (E05) became more important at higher latitudes, particularly near the western boundary currents.

(5) The development of a CP El Niño was well predicted in E05 and E07, thus indicating that it is important to jointly assimilate SST, SSH, and SSS for ENSO predictions.

To apply the data assimilation method in actual operations, the forecast run number can be decreased to 1 in LETKF. This time evolving run is combined with 100 static ensemble members to build the background error covariance for data assimilation. The results showd that increasing the number of time evolving ensemble members from 1 to 5 was not effective for obtaining more accurate ensemble covariance using the LETKF assimilation scheme. If the number increases to 20, the covariance is more accurate. This has been demonstrated in many previous studies such as Miyazawa et al. (2012), Xu et al. (2013), and Xu and Oey (2014). However, the computer costs of more than 10 ensemble numbers limits the ability of this approach in operational applications.

The inclusion of large errors in the initial conditions was aimed at testing the ability of the ocean assimilation system to correct the errors. A series of experiments similar to the above experiments but initialized from January 5, 1990 are conducted to reduce initial errors. Note that TRU was initialized starting from January 1, 1990. The results showed improvements of the analysis and forecasts as well, though not as significant as those from the experiments listed in Table 1.

The newly developed system was tested in a twin experiment framework. This approach allowed for extensive tests of system accuracy, and such an approach has been widely used in data assimilation studies (e.g., Counillon et al., 2014; Zhou et al., 2016). Even though the results were encouraging, the next plan is to conduct a comprehensive test in a realistic framework before the system is put into operation.

Finally, coupled assimilation with ice properties and atmospheric fluxes has been shown to be advantageous in climate system analyses (e.g., Zhang et al., 2009; Zheng and Zhu, 2010). Thus, the assimilation of those variables will be conducted in future studies.