The development of the Istituto Nazionale di Geofisica e Vulcanologia (INGV)–Centro Euro-Mediterraneo per i Cambiamenti Climatici (CMCC) Seasonal Prediction System (SPS) is documented. In this SPS the ocean initial-conditions estimation includes a reduced-order optimal interpolation procedure for the assimilation of temperature and salinity profiles at the global scale. Nine-member ensemble forecasts have been produced for the period 1991–2003 for two starting dates per year in order to assess the impact of the subsurface assimilation in the ocean for initialization. Comparing the results with control simulations (i.e., without assimilation of subsurface profiles during ocean initialization), it is shown that the improved ocean initialization increases the skill in the prediction of tropical Pacific sea surface temperatures of the system for boreal winter forecasts. Considering the forecast of the 1997/98 El Niño, the data assimilation in the ocean initial conditions leads to a considerable improvement in the representation of its onset and development. The results presented in this paper indicate a better prediction of global-scale surface climate anomalies for the forecasts started in November, probably because of the improvement in the tropical Pacific. For boreal winter, significant increases in the capability of the system to discriminate above-normal and below-normal temperature anomalies are shown in both the tropics and extratropics.
The scientific basis for seasonal predictions lies in the interaction of the atmosphere with slowly varying components of the climate system such as the ocean (e.g., Navarra 2002; Shukla and Kinter 2006). Early studies showed that El Niño can be predicted seasons in advance using numerical models of the coupled ocean–atmosphere covering the tropical Pacific (Cane et al. 1986; Zebiak and Cane 1987). Since then, there have been many developments in forecasting sea surface temperature (SST) anomalies in the tropical Pacific (e.g., Latif et al. 1998; Palmer 2006; Balmaseda et al. 2007).
The ability of models to predict ENSO is critically important, as the most significant climate variability on the interannual time scale is related to this phenomenon (Ji and Leetmaa 1997; Trenberth et al. 1998; Wallace et al. 1998) and as the SSTs in the tropical Pacific have a global impact on atmospheric circulation (e.g., Shukla and Wallace 1983; Trenberth et al. 1998). However, models skill in predicting tropical Pacific SSTs is still limited. For instance, most of the seasonal prediction systems underestimate or do not predict the onset of the exceptional 1997/98 El Niño (e.g., McPhaden 1999; Vitart et al. 2003).
The importance of the oceanic subsurface memory—as expressed by slow variations in the equatorial Pacific upper-ocean heat content—for the evolution of ENSO has been shown in many observational and modeling studies (e.g., Chen et al. 1995; McPhaden et al. 1998; Latif et al. 1998; McPhaden 1999; Navarra et al. 2008). Subsurface data assimilation can contribute in obtaining skillful seasonal forecasts and beneficial effects on predictability have been also reported increasing the space–time coverage of the observational network (Rosati et al. 1997; Alves et al. 2004; Ji and Leetmaa 1997; Wang et al. 2002; Balmaseda et al. 2007; Vidard et al. 2007). These studies have often showed reductions of the model errors. However, the reduced model errors do not always correspond to significant increases in prediction skill and the results are substantially depending on the model, the geographical region, the year, and the season under consideration (e.g., Ji and Leetmaa 1997; Vidard et al. 2007; Balmaseda et al. 2007).
To be useful for decision making, seasonal climate predictions need to be probabilistic and the capability of probability forecasts to provide valuable information needs to be assessed (e.g., Richardson 2006). Specifically, it would be desirable to evaluate how well a set of probability forecasts is able to discriminate among the occurrence of mutually exclusive and collectively exhaustive climate events, with the simplest possible situation represented by dichotomous yes–no cases (e.g., temperature above normal or not). Considering predictions of dichotomous events, the joint distribution of the observed predictands and of the respective probability forecasts can be conveniently analyzed through the likelihood-base-rate factorization (Wilks 2006). This kind of analysis can be used for a direct quantification of the forecasts ability to discriminate among the occurrence of one or the other of a pair of dichotomous events.
This work documents the development of the Istituto Nazionale di Geofisica e Vulcanologia (INGV)–Centro Euro-Mediterraneo per i Cambiamenti Climatici (CMCC) Seasonal Prediction System (SPS), which includes an assimilation of in situ vertical profile observations in the oceanic model in order to produce initial conditions (ICs). The ocean-data assimilation system has been developed at CMCC-INGV (Di Pietro and Masina 2009; Bellucci et al. 2007) and it has been used in order to assimilate observed profiles of temperature and salinity through the water column. The main focus of the paper is on the assessment of the impact of the assimilated initial conditions on the forecast skill.
The paper is organized as follows. Section 2 describes the seasonal prediction system, the experiments performed, and the data used for validation. After a description of systematic errors and skill performance of the latest release of the system, sections 3 and 4 contain the comparison with control forecasts (i.e., without ocean assimilation) to analyze the effects of the improved ocean ICs on SST bias and on prediction skill, respectively. The effects of the assimilation on the skill of global-scale probability forecasts of dichotomous predictands are addressed in section 5. Section 6 contains the discussion of the main results and a summary of the conclusions of the study.
2. The Seasonal Prediction System
The SPS documented in the present study represents the evolution of the system described in Gualdi et al. (2004) and developed in the framework of the European Union (EU) project Development of a European Multimodel Ensemble System for Seasonal-to-Interannual Prediction (DEMETER; Palmer et al. 2004). The ICs for the ocean–atmosphere system are prepared separately for the atmosphere and for the ocean. For each start date the atmospheric ICs are obtained from prescribed SST simulations. Differently, the ocean component is obtained from the ocean data assimilation or simply from a flux-forced ocean simulation for the control forecasts. Figure 1 summarizes the hindcasts generation strategy of the latest versions of our seasonal prediction system. The details of the setup and integrations performed in this study are described in section 2c.
a. The coupled model
The coupled model included in the system is the ocean–atmosphere coupled general circulation model (CGCM) Scale Interaction Experiment-Frontier (SINTEX-F; Gualdi et al. 2003a,b; Luo et al. 2005). SINTEX-F is an evolution of SINTEX (Gualdi et al. 2003b), where both oceanic and atmospheric components have been improved. The model components are Ocean Parallelise (OPA) 8.2 for the ocean and ECHAM4 for the atmosphere.
OPA8.2 (Madec et al. 1998) is used in the ORCA2 global implementation (see details at http://www.nemo-ocean.eu/). It is a finite-difference oceanic GCM and solves the primitive equations with a nonlinear equation of state on an Arakawa C grid. The horizontal mesh is orthogonal and curvilinear on the sphere, and its spatial resolution is roughly equivalent to a geographical mesh of 2° × 2° with a meridional resolution of 0.5° near the equator. A total of 31 vertical levels are used with 10 levels in the top 100 m. In the configuration used there is no interactive model for the dynamics of sea ice, whose area coverage is relaxed toward observed monthly climatology.
ECHAM4 (Roeckner et al. 1996) is the fourth generation of the ECHAM atmospheric general circulation model developed at the Max-Planck-Institut Fur Meteorologie in Hamburg, Germany. The model equations are solved on 19 hybrid vertical levels (top at 10 hPa) by using the spectral transform method. In these simulations, ECHAM4 is used with a triangular truncation T106, corresponding to an associated Gaussian grid of approximately 1.1° × 1.1°. As shown in Gualdi et al. (2004), this relatively high resolution improves considerably the prediction skill compared to a coarser atmosphere (T42) and gives a better representation of the delayed oscillator mechanism (Navarra et al. 2008). An exhaustive description of the dynamical and physical structure, and of the simulated climatology of ECHAM4, is given by Roeckner et al. (1996).
Atmospheric and oceanic components are coupled through the Ocean Atmosphere Sea Ice Soil coupler (OASIS2.4; Valcke et al. 2000). No flux adjustment or restoring was used in the simulations. Air–sea fluxes and SST between atmosphere and ocean were exchanged every 2 h. The features of the SINTEX-F climatology and variability have been widely described in the past (e.g., Gualdi et al. 2003a,b).
b. Ocean assimilation
The latest version of the CMCC-INGV Global Ocean Data Assimilation System (CIGODAS; Di Pietro and Masina 2009; Bellucci et al. 2007) has been used in order to assimilate observed profiles of temperature and salinity through the water column of the global configuration of the OPA8.2 ocean model. The assimilation scheme used in CIGODAS is based on the System for Ocean Forecasting and Analysis (SOFA; De Mey and Benkiran 2002), which is a reduced-order multivariate optimal interpolation scheme. As described in Bellucci et al. (2007), CIGODAS considerably corrects the subsurface thermal structure of the oceanic model. In particular, tropical Pacific and western boundary currents regions show a beneficial impact from the assimilation. Details of the CIGODAS and of the effects on the ocean model simulated climatology and variability are found in Di Pietro and Masina (2009) and Bellucci et al. (2007).
The temperature and salinity profiles used for this study are taken from the EN3 package [an assembling of the World Ocean Database 2005 (WOD05), the Global Temperature–Salinity Profile Program (GTSPP), and Argo databases, as summarized in Table 1; more information is available online at http://hadobs.metoffice.com/en3/]. Only the profiles that passed all the quality checks described in Ingleby and Huddleston (2007) have been retained for assimilation. The temperature profiles spatial coverage over different latitudes and regions of the globe are reported in Table 2 for the period 1991–2003.
c. Experiments and data
We performed an experiment using 5-month seasonal forecasts for the period 1991–2003. To consider the possible impact of the seasonal cycle on the forecasts, the simulations are started from two different dates of the year: 1 May and 1 November. Two sets of nine-member ensemble forecasts have been produced taking the same ICs for all the coupled model components but the ocean. In the first set the ocean initial states were estimated through the use of the data assimilation system described in section 2b (hereafter DAS) while in the second no observed in situ data were assimilated (hereafter NODAS). For both DAS and NODAS, the ocean model was forced starting from 1955 with momentum, heat, and freshwater flux data from the 40-yr European Centre for Medium-Range Weather Forecasts (ECMWF) Re-Analysis (ERA-40; Uppala et al. 2005) before 2002 (ERA-40 only covers up to August 2002) and from the ECMWF operational analysis after 2002. Furthermore, in order to keep the simulated SST close to observations, the model field was damped with a time scale of 7 days toward the Reynolds SSTs (Reynolds and Smith 1994) from 1982 onward and the ERA-40 SSTs before. In practice, the NODAS experiment ICs are simply forced from atmospheric fluxes and relaxation to surface SST, whereas the DAS experiment further included assimilation of in situ ocean profiles.
The scientific basis for seasonal predictions lies in the interaction of the atmosphere with slowly varying components of the climate system. As such, seasonal climate predictions are believed to be first an initial value problem for the slow ocean component (Palmer 2006; Shukla and Kinter 2006), while the solution for the atmosphere can be conceived as a boundary value problem (Navarra 2002; Palmer 2006). Consistently, the atmospheric initial conditions were obtained through an Atmospheric Model Intercomparison Project (AMIP)-type simulation (i.e., by prescribing observed SST boundary forcing to the atmospheric model). It is argued that, compared to the method based on atmospheric data assimilation (as widely developed for numerical weather predictions), with the AMIP-type approach atmosphere and oceanic conditions (i.e., SST) are more in balance and this could minimize the initial coupling shock (e.g., Tribbia and Troccoli 2008). The AMIP-type run was performed by using the observed SSTs from the Met Office Hadley Centre’s Global sea ice and SST dataset (HadISST1.1; Rayner et al. 2003) for the period 1985–2003. To represent the uncertainties in the initial state of the system, an ensemble of nine atmospheric ICs has been produced by taking lagged days as initial states. For each starting date we consider the reference date (1 May or 1 November) but also the 4 days before and after (Fig. 1).
In summary, for each starting date an ensemble of nine atmospheric initial states were created. Starting from these ICs, for both DAS and NODAS (i.e., with and without oceanic assimilation of in situ data), the coupled model has been integrated for five months, producing two sets of nine-member ensemble forecasts covering the period 1991–2003.
1) Observed and reanalysis datasets
The predictive skill of the model is assessed comparing the forecasts with analyses and observational products. The ERA-Interim reanalysis (Uppala et al. 2005; Berrisford et al. 2009) is used for verification of the forecasts. For precipitation, we use the Climate Prediction Center (CPC) Merged Analysis of Precipitation (CMAP; Xie and Arkin 1997) dataset. The model and the observed (reanalysis) anomalies are defined as the deviations from the respective climatology for the period 1991–2003.
3. Improved ocean IC effect on SST bias
An important source of inaccuracy of the forecasts performed using CGCMs is represented by the model systematic errors (e.g: Gualdi et al. 2004). This problem is especially important for the prediction of the SST anomalies. The magnitude of the model SST bias, in fact, can be as large as the amplitude of the observed SST anomalies that should be predicted. In this section, we present a description and discussion of the effect of oceanic subsurface assimilation on the bias of the SPS forecasts. Before analyzing the details of the impact of the improved ocean IC, we discuss briefly the main features of the bias in the DAS experiment, with the main focus in the tropical Pacific.
a. Bias of the system
In Fig. 2 the DAS SST systematic error for both May (Fig. 2a) and November (Fig. 2b) start dates are shown. The biases are defined as the difference between the forecast ensemble means and ERA-Interim SST climatologies of the period 1991–2003. Months 2 to 4 of the forecast period are considered in the average, which means that for the 1 May start dates we used the average of the monthly means for June, July, and August while monthly means for December, January, and February are used for the 1 November forecasts. The results shown in Fig. 2 indicate that the systematic error of the model in predicting the SST field is moderate in most of the tropics. Over a large portion of the tropical belt, in fact, the model exhibits an averaged cold bias less than 1°C. The error is remarkably small in the tropical Indian Ocean. In the southeastern tropical Pacific and in the upwelling regions off the American coast and in the Gulf of Guinea the model is too warm and the averaged error is greater than 1°C. The equatorial cold tongue is too pronounced and penetrates too far into the west Pacific, producing SST patterns too symmetrical around the equator. It follows a tendency to produce a double intertropical convergence zone (ITCZ) in the tropical Pacific, consistent with what is described in Gualdi et al. (2004).
Figure 2 shows a seasonal dependency of the systematic errors on the date of the ICs. For example, the SST warm bias in the tropical southeastern Pacific and Atlantic Oceans and Southern Hemisphere midlatitudes appears to be more pronounced for the forecasts with start date in November (Fig. 2b), whereas a warm bias in the tropical northeastern Pacific and Atlantic is found in the forecasts with start dates in May (Fig. 2a). The error in the equatorial Pacific cold tongue, on the other hand, is more evident for the forecasts starting in November.
b. Sensitivity to ocean assimilation
Figure 2 compares the systematic error in the DAS SST (Figs. 2a,b), averaged over the 2–4 forecast months, with the one suffered by NODAS (Figs. 2c,d). Figures 2e,f report the SST bias difference between DAS and NODAS, with the shading evidencing the areas of significant (10% level, bootstrap method) systematic error increase (light) and decrease (dark) in DAS. The assimilated ocean IC estimate in DAS leads to a reduced mean bias over the tropical belt. For both the May (Figs. 2a,c,e) and November (Figs. 2b,d,f) start date, the warm bias in the upwelling regions off the American Coasts and in the Gulf of Guinea is significantly reduced over large areas. Similarly, the systematic error in the subtropical south-central Pacific is reduced. The SST bias in the regions influenced by the western boundary currents in the extratropical Pacific and Atlantic appears also to be affected by the subsurface profile assimilation. During boreal summer, the temperature bias over the regions influenced by the Kuroshio is considerably reduced over large areas. Similarly, the forecasts started in November display a reduction of the SST systematic error over the Gulf Current regions in the Atlantic sector poleward of 45°N. In contrast, in the Atlantic southward of 45°N, the SST bias in some regions bounded by the Gulf Current appears to increase in DAS.
As the tropical Pacific cold tongue region was shown to be the region with the strongest bias for the forecasts started in both May and November, we focused on the mean systematic errors averaged over the Niño-3.4 region (5°S–5°N, 170°–120°W) for DAS (dashed) and for NODAS (dashed–dotted) experiments (Fig. 3). The forecasts started in May have relatively moderate drifts, reaching ∼1 K after 5 months in both NODAS and DAS (Fig. 3). Stronger drifts, exceeding 1.5 K after 5 months, are found for the forecasts starting in November. The cold bias of ∼0.5 K already present in the first month (Fig. 3a) indicates rapid adjustments going on due to quite a prominent “coupling shock.” This initialization-related adjustment appears to be more effective for the November start date.
The DAS forecasts have a smaller drift than NODAS during the first part of the predictions. In particular, the bias of the first month forecasts is reduced by ∼0.25 K (≃40%) in DAS for the November start date (Fig. 3a). However, the bias tends to converge to similar values in the latter stages of the predictions (from month 4 onward), as the forecasts tend toward the modeled climatology (e.g., Alves et al. 2004; Jin and Kinter 2008). Our results indicate that this tendency is reduced for the first 3 months of the forecasts by the use of consistent subsurface temperature and salinity information in order to initialize the ocean component. Please note that the uncertainty associated with analyzed SST data (based on satellites and in situ measures) is on the order of 0.2 K (e.g., Rayner et al. 2006) and this could limit the significance of the systematic error differences reported above.
4. Improved ocean IC effect on predictability
a. Skill of the system
Figures 4a,b show the point-by-point correlations between predicted (DAS) and observed (ERA-Interim) surface air temperature anomalies (hereafter TAIRA). Time correlations are computed retaining for each year all the monthly means from the 1-month lead-time seasonal predictions (forecast months from 2 to 4: June, July, and August for 1 May start dates and December, January, and February for 1 November start dates). Higher correlation values are found for both the May start date (Fig. 4a) and the November start date (Fig. 4b) over the tropical Pacific. Positive significant correlations are also found over most of the tropical Indian and Atlantic Oceans with relatively higher values and significant area coverage for the forecasts started in November. The high correlations over the tropical Pacific display the skill of our coupled model in predicting ENSO (see the Niño-3.4 index in Fig. 5). As summarized in Table 3 the SPS performs particularly well in predicting the Niño-3.4 index both for the seasons with 1-month lead time (months from 2 to 4, hereafter lead-1 seasons) and for the seasons with 2-month lead time (months from 3 to 5, hereafter lead-2 seasons). The correlation between the predicted and observed monthly Niño-3.4 index always exceeds 0.9 for both May (0.94 and 0.91 for lead-1 and lead-2 seasons, respectively) and November (0.97 at lead 1 and 0.95 at lead 2) start dates. Correspondingly, the root-mean-square error (RMSE) is moderate with values below 0.4 for November. In the cases with a start date in May RMSE is 0.44 and 0.53 for lead-1 and lead-2 seasons, respectively (Table 3).
From the tropical Pacific, the positive significant correlations tend to irradiate toward the whole tropical belt and toward the extratropics (Figs. 4a,b). Some positive correlations are also found in land regions strongly influenced by ENSO teleconnections (Shukla and Wallace 1983; Trenberth et al. 1998). During boreal winter, such regions are identified over central and southern Africa, the Amazon basin, and southeastern South America, the midlatitude western Pacific coasts of Asia, and northwestern as well as northeastern North America. Boreal summer positive correlations are found in northern Australia and the Indonesian Archipelago, Gulf of Mexico, and Central America and over southeastern Asia. Remarkably, significant correlations between the forecasts and observed temperatures are found in the Euro–Mediterranean region as well as the Middle East during the boreal summer, suggesting some predictability over these regions for the forecasts started in May. However, most areas evidencing significant correlations are found over the oceans (Figs. 4a,b). This indicates a lower predictability over lands and it may follow at least in part from the absence, in the experiments we performed, of any kind of assimilation in order to suitably initialize the long persisting land variables (e.g., Koster et al. 2004, 2006; Ferranti and Viterbo 2006; Alessandri and Navarra 2008).
Although the Niño-3.4 index shown in Fig. 5 is frequently used to characterize the state of the ENSO and to quantify the quality of simulations and predictions of the oscillation, it describes the averaged SST over only a small portion of the Pacific Ocean. The ocean anomalies associated with ENSO, on the other hand, affect the whole tropical basin with the development of widescale anomaly patterns. It is therefore of interest to check the skill of the forecasting system to predict the evolution of the SST pattern anomalies over the entire tropical Pacific. To this aim, spatial anomaly correlation coefficients (ACCs) and spatial RMSEs for the forecast SST anomalies in the tropical Pacific (defined as 25°S–25°N, 140°E–80°W) have been computed. Figure 6 reports the RMSE (Figs. 6a,b) and the ACCs (Figs. 6c,d) for the forecasted SST anomalies in the tropical Pacific together with the results obtained for persistence forecasts (dashed line). Thick solid lines (and filled circles) are the ensemble means while thin lines stand for each ensemble member of the DAS experiment. Filled triangles also report the results for the NODAS ensemble mean forecasts—see section 4b for a comparison between DAS and NODAS. The ACCs and the RMSEs have been computed relative to the ERA-Interim SST and persistence forecasts are made by continuing the ERA-Interim monthly anomaly from the month prior to the start date of the model forecasts. For example, SST persistence forecasts for the period May–September 1991 have been made by continuing the SSTA found for the observed April 1991. The results in Fig. 6 indicate that the ensemble mean forecasts have better skill than any ensemble member. This is in agreement with the results found with other coupled model forecast systems, and can be explained, at least in part, by the fact that the ensemble average reduces the internal dynamics noise present in the individual forecasts (Kirtman and Shukla 2002), thus potentially increasing the correlations and reducing the RMSE.
Overall, DAS displays a good skill in reproducing the observed SST over tropical Pacific, with the ensemble mean forecasts (hereafter, simply “forecasts”) usually performing better than persistence forecasts (dashed lines), especially on lead times greater than 1 month. However, both the anomaly correlations and the RMSEs show some seasonal dependence, with some evidence of the so-called spring predictability barrier in our system. A higher capability to beat the persistence forecasts, in fact, is found for the forecasts that start in May compared to the forecasts starting in November. These results are consistent with several previous studies (e.g., Gualdi et al. 2004; Schneider et al. 2003).
b. Sensitivity to ocean assimilation
Figure 4 compares the correlations between observed TAIRA and DAS (Figs. 4a,b) and NODAS (Figs. 4c,d) 2–4-month forecasts, respectively. Figures 4e,f report the DAS minus NODAS difference in correlations, with the shaded areas evidencing grid points with significant (10% level, bootstrap method) increase (red) and decrease (blue) in DAS. To better illustrate the skill difference as shown between Figs. 4a,b versus Figs. 4c,d, in Table 4 we show the fractional areas with significant correlations as well as with correlations exceeding given thresholds for both DAS and NODAS. Overall, the global area fraction with significant correlations between November forecasts and observations is 0.60 in DAS and only 0.55 in NODAS (Table 4). For the forecasts started in May the difference in global fractional areas reduces to 0.03 (0.61 in DAS versus 0.58 in NODAS). Considering the areas with correlations above 0.8, DAS still does better than NODAS with 0.08 versus 0.05 of the global area exceeding this threshold in the forecasts started in November (Table 4). In contrast, the boreal summer forecasts do not display any difference in fractional area with correlations above 0.8 (0.05 in both DAS and NODAS). The areas evidencing significant correlation difference between DAS and NODAS (Figs. 4e,f) are mostly placed over subtropical and midlatitude oceans. For the forecasts started in November (Fig. 4f), DAS displays areas with significant increase over the equatorial central Pacific, subtropical central South Pacific, subtropical Indian and Atlantic Ocean, as well as the Northern Hemisphere western boundary current regions. DAS also shows evidence of some significant correlation increases over the continents, and in particular, the coastal areas adjacent to the North Pacific and Atlantic Oceans. The May start date forecasts have fewer grid points that display significant correlation improvements in DAS (Fig. 4e), which appears to be outperformed by NODAS over large areas surrounding the Kuroshio. On the other hand, DAS is significantly better than NODAS over the northeastern Pacific.
1) Tropical Pacific
The prediction of SST over the tropical Pacific appears to be affected by the assimilation of temperature and salinity with stronger improvement for the November start date. From Figs. 6b,d it is clear that DAS (filled circles) improves noticeably compared to NODAS (filled triangles), particularly in terms of ACCs from the third month on. Less clear is the impact for the May start date (Figs. 6a,c) where both the ACC and RMSE values appear to be only slightly affected. While DAS improves to some extent in the forecast months 1 and 4, the results are uncertain for month 5 and the skill is almost identical for months 2 and 3.
The significance of the above results can be more clearly evaluated by means of scatterplots of the ACCs and RMSEs computed on the 1-month lead-time seasonal mean predictions (averages of forecast months from 2 to 4: June, July, and August for 1 May and December, January, and February for 1 November start dates). Figure 7 compares the forecasts performed with assimilated ICs (DAS) with the control (NODAS) in the tropical Pacific. ACCs (top panels) and RMSEs (bottom panels) for each forecast year (diamonds) as well as the average of the values over all the 13 forecast years (crosses) are displayed. In some years, an increase of the ACC for DAS in both May and November is visible (Fig. 7). Using a Monte Carlo bootstrap procedure, we checked the significance of the difference in the 13-yr means (Table 5). We found that the 5% level of significance is verified only for the November case. The ACCs for the lead-2 seasons (averages of forecast months from 3 to 5) display similar results (Table 5), evidencing a significant improvement in DAS for the November start dates. The impact of the assimilation of temperature and salinity profiles on the SST RMSE appears to be smaller than that for ACCs with the 13-yr averages, which do not pass the significance test for the difference at the 5% level. For completeness, we checked the results also over the Niño-3.4 region (Table 6). Compared to the tropical Pacific the ACCs decrease considerably, reflecting the fact that this smaller region has a strong SST signal and may develop anomaly patterns characterized by steep gradients. Nevertheless, for the lead-1 season, the effect of subsurface assimilation is similar to what evidenced by the comparison over the whole tropical Pacific basin. In contrast, for the lead-2 season, our analysis does not show evidence of any significant difference at the 5% level between DAS and NODAS over the Niño-3.4 region (Table 6).
A synthesis of the time–space SST variability over the tropical Pacific is reported in Table 7. For each forecast month, the standard deviations for the ensemble mean SST anomalies in DAS (first column), NODAS (second column), and for the observations (third column) are computed retaining both the space and interannual time variability. Compared to NODAS, in both May and November start dates the DAS ensemble mean prediction variability is closer to the observed value, in particular for the forecast months from 3 to 5. This is due to the fact that DAS do not display the marked progressive weakening of the predicted ensemble mean anomalous signal that characterizes NODAS. This result indicates an increased signal-to-noise ratio of the forecasts performed in DAS from months 3 to 5, driven by the improved subsurface initialization of the ocean. As previously discussed, this appears to produce an increased skill only for the forecasts started in November; for May the enhanced signal appears not to correspond to a better fit to the observed anomalies.
2) 1997/98 El Niño
The Niño-3.4 index for all start dates in the forecasts performed with NODAS (black) and DAS (red) is reported as a box plot time series in Fig. 8. The distribution of predicted monthly mean anomalies is represented by boxes (25th–75th percentiles) and the median is represented by the inside box mark. The shaded band includes the interannual standard deviation in the observations while the dashed lines refer to the forecasts. DAS appears to better represent the observed (green filled circles) anomalies in the Niño-3.4 index compared to NODAS (correlation coefficient 0.94 versus 0.91; see also Table 8 for the comparison of the skill in the lead-1 and lead-2 seasons for each start date). The evolution of the two major El Niño events (1991/92 and 1997/98) appear to be better represented in DAS. Particularly, DAS considerably improves the onset amplitude of 1997/98 El Niño. The strength and length of the 1991/92 El Niño, staying well above one standard deviation until May 1992 in the observations, is also better represented in DAS. However, DAS and NODAS appear to anticipate the onset of the 1991/92 El Niño, both providing a false alarm in the forecast started May 1991. The considered forecasting period also includes the 1994/95 and 2002/03 moderate-intensity El Niño events. They appear to be quite well predicted by both DAS and NODAS with the 1994 El Niño onset better captured in NODAS. On the other hand the 2002 onset appears to be slightly better represented in DAS.
The observed 1997/98 El Niño onset amplitude is very well reproduced by DAS, with an observed averaged anomaly of ∼2.3 K in September 1997. In contrast, NODAS reaches an anomaly of only 1.5 K at that time. The Hovmöller diagram of the evolution of the heat content anomaly (shaded) and zonal wind stress anomaly (contours) for NODAS (Fig. 9a), DAS (Fig. 9c), and the ocean analysis (Fig. 9b) is reported in Fig. 9. Note that for the model, ensemble-mean anomalies are reported, thus averaging out the “nonsignal variability” intrinsically working at the shorter time and smaller space scales (Kirtman and Shukla 2002). As shown in Fig. 9, DAS well represents the initial positive heat content anomaly along the equatorial Pacific and off the coast of Peru. Similarly, the Kelvin wave train traveling eastward through the tropical Pacific as well as the associated heat content anomaly propagation appears to be captured correctly by DAS. In contrast, the NODAS heat content anomaly is much weaker than observed and the Kelvin wave signal appears to be less clear and somewhat delayed compared with observations. This results in a reduced development and west–east propagation of the heat anomaly and in the consequent weak simulation of the onset of the 1997/98 El Niño in NODAS. It is important to note here that this weak initial development of the 1997/98 El Niño has been pointed out to be a major problem for most of the dynamical as well as statistical ENSO forecast models (e.g., McPhaden 1999; Vitart et al. 2003).
Interestingly, the second Kelvin wave train does not appear to be captured by DAS, probably because of uninitialized intraseasonal wind bursts (McPhaden 1999). Nevertheless the initial heat content anomaly to the east of the date line is correctly propagated eastward in this experiment. The resulting heat content prediction is still close to that observed at the end of the DAS forecast.
5. Probabilistic forecasts of dichotomous predictands
In the previous section we have evaluated the impact of the improvement of the oceanic ICs estimation on the prediction of the tropical Pacific SST, which represents the main source of climate predictability at the seasonal time scale (e.g., Ji and Leetmaa 1997; Trenberth et al. 1998; Wallace et al. 1998). In the following we will give an estimate of the associated global-scale impact on the performance of probabilistic forecasts of dichotomous observed predictands. In particular, we will concentrate on the forecasting of below-normal (i.e., below lower tercile of the sample distribution) and above-normal (i.e., above upper tercile of the sample distribution) TAIRA and precipitation anomalies (hereafter PRECA).
The joint distribution of forecasts and observed dichotomous predictands can be conveniently analyzed and displayed graphically through the likelihood base-rate factorization (Wilks 2006):
Here the conditional distributions p(yi|oj) express the likelihoods that each of the I allowable discrete forecast values yi would have been issued in advance of each of the observed dichotomous events oj (occurrence j = 1; no occurrence j = 0). Together with the associated sample climatological probabilities p(oj), it completely represents the information of the full joint distribution. Specifically, the conditional likelihood distributions, p(yi|oj), are directly indicative of how well a set of forecasts are able to discriminate among the events oj. Graphically, this can be appreciated through diagrams consisting of superimposed plots of the two likelihood distributions as a function of the forecast probability yi (hereafter discrimination diagrams). As pointed out by Wilks (2006), the previously mentioned characteristics of the two likelihood distributions could be used effectively to recalibrate the probability forecasts by calculating posterior probabilities for the two events given each of the possible forecast probabilities.
Differently from the previous sections, which were mostly focused on climatology and predictability of seasonal means, here we focus on the probabilities of monthly means to exceed tercile thresholds and on the capability of the forecasts to discriminate the corresponding observed outcomes. To this aim we retain all the forecast months from 2 to 5 in the analysis that follows. Figure 10 shows a comparison between DAS and NODAS of the discrimination diagrams evaluated over the tropics (25°S–25°N, 0°–360°) and considering all the monthly means from months 2 to 5 of the forecasts started in November (Figs. 10a,b; December, January, February, and March monthly means considered) and in May (Figs. 10c,d; June, July, August, and September monthly means considered). Figures 10a,c relate to the forecast probabilities of the event of TAIRA being below the lower tercile of the climatological sample distribution (ET−). In contrast, Figs. 10b,d refer to the forecast probabilities to exceed the upper tercile (ET+). In Fig. 10 the dashed lines represent the likelihood distribution given the no occurrence of the event [p(yi|o0)], while the solid lines are the likelihood distributions verified o1 [p(yi|o1)]. Both NODAS in red and DAS in blue are displayed in the same diagram. In the forecasts started in November (Figs. 10a,b), DAS displays larger p(yi|o0) values (dashed lines) for the smaller forecast probabilities for both ET− and ET+. Similarly, conditional probabilities given o1 (solid lines) are higher in DAS compared to NODAS for the larger probability forecast outcomes. This determines an increase in DAS of the separation between the two respective likelihood distributions, indicating an improved ability of the forecasts to discriminate warm and cold events.
The separation of the two likelihood distributions is also plotted in the same figure (bottom horizontal bars) in the form of the discrimination distance d, a scalar attribute defined as the difference between the means of the two likelihood distributions μ following Wilks (2006):
In Fig. 10, asterisks are placed in correspondence of the discrimination values to indicate that they are significantly (at the 5% level, using a Monte Carlo bootstrap method) higher than the other experiment. DAS increases significantly d compared to NODAS for both and (Figs. 10a,b). This result indicates that, for the forecast started in November, the assimilation of temperature and salinity profiles improves the ability of the model to discriminate between the occurrence of below-normal or above-normal events over the tropics. Differently from the November start date case, the forecasts started in May (Figs. 10c,d) do not display any significant (5% level) difference in d. In fact, considering both and , the two likelihood distributions, compared between DAS and NODAS, are, respectively, very close to each other, resulting in discrimination distances almost identical.
The panel insets in the discrimination diagrams in Fig. 10 also report the refinement distributions, p(yi). The dispersion of p(yi) reflects the overall confidence of the forecasts, so that forecasts that deviate rarely from their average value exhibit little confidence (Wilks 2006). Interestingly, all the refinement distributions in Fig. 10 evidence an increased confidence in DAS. This is verified even when the discrimination distance is not affected considerably, as it is the case for the May start date. This result indicates that the addition of subsurface information to the ocean IC tends to increase the signal-to-noise ratio of the predictions and consequently drives the predictions to fall preferentially outside normal conditions. However, consistently with the results in section 4, this produces increased discriminations only for the November start date, while for May the enhanced signal does not appear to coincide with better correspondence to observations. Table 9 summarizes, for both November and May, the discrimination distances in the tropics, northern extratropics, and southern extratropics for and , where T and P indicate temperature and precipitation, respectively. Considering precipitation over the tropics, DAS shows an enhanced (5% significance level) and, interestingly, also the forecasts started in May display an improvement due to the initialization of the subsurface ocean. In fact increases significantly compared to NODAS.
From Table 9, it is shown that discrimination distances outside the tropics tend to decrease to values close to or below 0.1 in both DAS and NODAS. Nevertheless, extratropics discrimination appears to be affected by the assimilation of temperature and salinity as well. Considering northern extratropics TAIRA during boreal winter, and are improved considerably in DAS compared to NODAS (5% significance level verified). Differently, the PRECA field discrimination is affected during boreal summer in the northern extratropics, with increasing significantly in DAS. Similarly to the Northern Hemisphere case, the comparison of the d values for the southern extratropics display improvements for TAIRA during boreal winter. In this case is increased in DAS compared to NODAS. Interestingly, displays a higher value in NODAS (not significant however). This is probably due to the fact that only a very small number of profiles is available in the southern extratropics during austral winter and spring for initialization purposes (see Table 2).
Table 10 reports the comparison of the discrimination distances between DAS and NODAS for the three tropical ocean basin sectors. For the tropical Pacific, results emphasize the analysis performed in section 4. Compared to NODAS, the DAS discrimination distances over this ocean basin increase significantly (5% level) in the November start date forecasts for all the cases considered . Differently, for the May start dates the discriminations appear to be little affected, with no significant differences between DAS and NODAS over the tropical Pacific. On the other hand, the assimilation of subsurface observations leads to the enhancement of the May prediction performance in the Indian Ocean. A significantly increased discrimination has been verified for , , and (Table 10).
The tropical Atlantic behaves differently from the other tropical oceans. Considering the May start date, in this basin the assimilation of subsurface temperature and salinity leads to a clear worsening of the predictions in terms of discrimination distance. In fact, the discriminations for and increase significantly in NODAS compared with DAS (Table 10). Figure 11a shows a clear reduction of the overlapping between p(yi|o0) and p(yi|o1) in NODAS compared with DAS and this ends up in a significantly (5% level) higher value. Figure 11d documents the increased NODAS discrimination for . Noteworthy in this case, the refinement distribution (inset histogram) shows much less confidence than for the temperature cases, and correspondingly the discrimination values appear to be reduced. Nevertheless, there is a considerable improvement of NODAS compared to DAS. The forecasts over the tropical Atlantic that started in November (Fig. 12) display improved temperature discriminations in DAS. In fact, and are significantly better discriminated by DAS (Figs. 12a,b). To understand the opposite behavior of boreal winter and boreal summer forecasts over the tropical Atlantic, we compared the subsurface thermal climatology of the ICs taken from the ocean analysis with a long free simulation of the coupled model (i.e., radiative boundary conditions of the 1991–2003 period were used). We found that the equatorial Atlantic subsurface thermal structure is very badly represented in the coupled model during boreal spring and early summer (not shown). During these seasons the coupled model shows an opposite slope of the thermocline with respect to observations. This result suggests that the subsurface correction, due to the data assimilation during spring, drives the coupled model too far from and not “in balance” with the state it would have and thus leads to a negative impact on the forecasts started in May.
Table 11 compares the DAS and NODAS discriminations over the Pacific–North American (PNA; 40°–65°N, 150°E–60°W) and the Euro–Atlantic (35°–65°N, 80°W–40°E) regions. For the PNA region, the results indicate increased discriminations in the forecasts started in November for the TAIRA prediction. Differently, the performances for the precipitation are in general very little affected. Even if discriminations appear to decrease considerably compared to tropical regions, the refinement distributions show good confidences for both November (Figs. 13a,b) and May (Figs. 13c,d) start dates. In this context DAS exhibits a considerably higher sharpness of the forecast distribution. This result appears to drive the increased discrimination in DAS, which is verified for significance (5% level) for and .
Figure 14 is the same as Fig. 13 but considering the Euro–Atlantic region. Compared to Fig. 13, it is noted that discrimination performance over this area is considerably lower than for PNA (see also Table 11). However, for the forecasts starting in November, DAS improves the TAIRA results if compared with NODAS. In particular, for the ET+ case there is a significant (5% level) increase (2.5 times) of the discrimination from the very low value of 0.011 to 0.026, leading to a noticeable increase in the performance. This is an important result, as it leads to increased predictability in the region. However, increased discrimination is achieved by reducing the overlapping in the “not conventional” direction; in fact, DAS displays lower p(yi|o0) values for the smaller forecasts probabilities and the p(yi|o1) are smaller for the larger forecasts probabilities. From a physical point of view, this is of course deplorable, as it means that the model tends to represent the opposite of reality in this region. Nevertheless, as soon as it is systematic, it is useful for predictions.
The forecasts performed using the assimilation of subsurface temperature and salinity profiles during ocean initialization (DAS experiment) have evidenced an enhanced signal-to-noise ratio of the predicted surface temperature anomalies. In particular, in the tropical Pacific the magnitude of the time–space variability of the ensemble-mean SST anomalies predicted by DAS appears to be very well simulated when compared with observations. It is an improvement over NODAS (i.e., without subsurface assimilation during initialization), which considerably underestimates the predicted anomalous SST in the tropical Pacific from the third forecast month onward. In terms of ACCs and RMSEs, subsurface initialization improves the November start date forecasts over the tropical Pacific. In fact, the averaged (over the 13 forecast years) ACCs computed on the seasonal mean anomalies increase significantly (5% level) in the tropical Pacific for the forecasts started in November (considering the 1-month lead-time predictions it increases by 8% with respect to NODAS). The impact on the SST RMSE is smaller (5% reduction for the 1-month lead-time predictions) and the significance in this case can only be verified at the 10% level. In contrast, for the forecasts started in May, ACC and RMSE values are only slightly affected by the subsurface oceanic initialization and the averages over all forecast years are not significantly modified.
NODAS shows the tendency to underpredict the ENSO anomalies. This has been reported as a major problem also for other dynamical as well as statistical ENSO forecast models (e.g., McPhaden 1999; Vitart et al. 2003). In particular, when no subsurface data are assimilated, our system appears to considerably underestimate the development of the 1997/98 El Niño. In contrast, the assimilation of temperature and salinity profiles in DAS leads to a considerable improvement of the representation of the initial positive heat content anomaly along the equatorial Pacific and off the coast of Peru. Similarly, the observed Kelvin wave train traveling eastward through the tropical Pacific as well as the associated heat anomaly propagation appears to be captured correctly by DAS. These results further highlight the importance of subsurface data assimilation for the improvement in forecasting the development and evolution of El Niño events.
The assimilated ocean initial conditions, and the resulting enhanced tropical Pacific prediction, have been shown to improve the probabilistic predictions of dichotomous events globally. Considering the tropical belt, DAS displays significantly (5% level) enhanced capability to discriminate both warm (above upper tercile of the sample distribution, ET+) and cold (below lower tercile, ET−) surface air temperature events in the forecasts started in November. For boreal winter, the discrimination distance of anomalous temperature forecasts increases by 7% for ET− and by 9% for ET+ over the tropics. In some cases, the enhancement in the discrimination distance has been shown to be significant (5% level) not only for temperature but also for precipitation. In particular, wet events appear to be significantly better discriminated during boreal winter (12% increase with respect to NODAS), while dry events improve in boreal summer (6% increase with respect to NODAS). In contrast to November, the forecasts started in May do not display a significant increase in the discrimination of anomalous temperature events considering the whole tropics. And this is even though the refinement distributions show increased confidence of the forecasts in a way similar to the November case. This result indicates that the addition of subsurface information to the ocean IC tends to increase the signal-to-noise ratio of the forecasts and consequently drives the predictions to fall preferentially outside normal conditions. However, this improves the skill in discrimination only for the November start dates, while for May the enhanced signal does not appear to coincide with significantly better correspondence to observations.
The tropical Atlantic is found to behave differently from the other tropical basins, with a marked and significant worsening in DAS of the discrimination distance for the prediction of cold and wet events during boreal summer. By contrast, the November forecasts improve the prediction of the anomalous temperature quite similarly to the other tropical basins. The results in the Atlantic are in agreement with previous studies evidencing the tendency to a degradation of the prediction skill over tropical Atlantic when subsurface data assimilation is used (e.g., Vidard et al. 2007; Balmaseda et al. 2007). Our analysis suggests that the subsurface correction due to the assimilation of observed profiles in spring leads to a change too drastic in the tropical Atlantic and drives the system too far from the state that the coupled model would have there—in this season the slope of the equatorial Atlantic thermocline in free coupled model long runs is opposite than observed—thus leading to a negative impact on the forecasts started in May.
The assimilation of subsurface data in the ocean IC impacts the extratropics as well in the forecasts started in November. For the boreal winter forecasts, the discrimination distances appear to increase significantly (5% level) for anomalous temperature events in both the northern midlatitudes (compared to NODAS, the discrimination of both warm and cold events increase by 24% and 38%, respectively) and the southern midlatitudes (cold events discrimination increases by 10%). The Pacific–North American sector displays the strongest improvement in the Northern Hemisphere midlatitudes. In this region, a significant increase in the discrimination for both warm (46% increase with respect to NODAS) and cold (39% increase with respect to NODAS) events for the forecasts started in November is evidenced. Remarkably, DAS increases the capability to discriminate warm dichotomous events during boreal winters even in the Euro–Atlantic region. However, this is achieved through smaller conditional likelihoods for the smaller forecast probabilities, given the nonoccurrence of the event. Concurrently, when the event occurs, smaller conditional likelihoods correspond to the higher forecast probabilities.
Dynamical weather and climate prediction is challenging: progress did not occur in the last 30 yr because of drastic breakthroughs but rather because of slow incremental progress and through a great deal of hard work (Shukla and Kinter 2006). This study has evidenced beneficial effects on the boreal winter forecasts derived from the subsurface initialization of the ocean model in our prediction system. However, the impact of the ocean data assimilation is small and mostly negligible for the forecasts started in May. Large SST biases in the central tropical Pacific and quite a prominent initial “coupling shock” have also been evidenced in our system. Thus, more efforts are needed on the initialization of the coupled model and on reducing SST biases as they are probably limiting the positive impact of subsurface assimilation in our system.
In this study we considered the climate predictability using our SPS at scales up to the first 5 months and mostly focusing on the 1-month lead-time seasonal (months 2–4) forecasts. The predictability of seasonal climate at longer lead times is also an important issue that requires a large effort and it will be considered in future works. Furthermore, given the relatively short prediction scales considered in this work, the atmospheric initial conditions may have an impact on the skill of the system. Recently, some dynamical forecasts of intraseasonal oscillations produced rather credible simulations of the Madden–Julian oscillation, with evidence of some prediction skill out to a lead time of about 2 weeks (Kim et al. 2007; Vitart et al. 2007). Whether or not this kind of atmosphere initialization can affect the prediction of ENSO and of climate at the seasonal time scale still needs to be evaluated.
Many thanks for the generous support from several CMCC staff members. We especially thank Alessio Bellucci, Fabrizio Massari, Loredana Amato, Alberto Troccoli, Antje Weisheimer, and Francisco D. Reyes for their help and for the inspiring discussions.
Corresponding author address: Andrea Alessandri, Centro Euro-Mediterraneo per i Cambiamenti Climatici, via Aldo Moro 44, 40127 Bologna, Italy. Email: firstname.lastname@example.org