First outcomes from the CNR-ISAC monthly forecasting system

Abstract. A monthly probabilistic forecasting system is experimentally operated at the ISAC institute of the National Council of Research of Italy. The forecasting system is based on GLOBO, an atmospheric general circulation model developed at the same institute. The model is presently run on a monthly basis to produce an ensemble of 32 forecasts initialized with GFS-NCEP perturbed analyses. Reforecasts, initialized with ECMWF ERA-Interim reanalyses of the 1989–2009 period, are also produced to determine modelled climatology of the month to forecast. The modelled monthly climatology is then used to calibrate the ensemble forecast of daily precipitation, geopotential height and temperature on standard pressure levels. In this work, we present the forecasting system and a preliminary evaluation of the model systematic and forecast errors in terms of non-probabilistic scores of the 500-hPa geopotential height. Results show that the proposed forecasting system outperforms the climatology in the first two weeks of integrations. The adopted calibration based on weighted bias correction is found to reduce the systematic and the forecast errors.


Introduction
Monthly and subseasonal dynamical forecasts are nowadays issued by various meteorological centres.For instance, the European Centre for Medium-Range Weather Forecasts (ECMWF) produces 32-day forecasts (Vitart et al., 2008), the Australian Bureau of Meteorology emits 6-week forecasts (Hudson et al., 2011), and the National Centers for Environment Prediction (NCEP) daily operates the Climate Forecasting System version 2 (Saha et al., 2012) issuing monthly means deduced from seasonal forecasts.
To obtain a successful extended-range forecast, the adopted forecasting system not only must rely on the atmospheric initial conditions but, also, it must be capable to catch the evolution of those slowly varying "boundary" processes whose possible initial anomalous state determine the future atmospheric anomalies (Palmer and Hagedorn, 2006).To accomplish this goal, because of the intrinsic limits of atmospheric predictability, an extended-range forecasting system must be probabilistic: by generating an ensemble of deterministic integrations, the forecasting system samples initial uncertainties to produce an evolution of the probability dis-tribution of the atmospheric state.If the predicted distribution correctly describes an anomalous atmospheric state with respect to the reference climatological distribution, the forecast can be considered skilful.Initial conditions of land and ocean have been found to positively influence predictability on the monthly time scale, mainly over extratropical and tropical areas, respectively (Chen et al., 2010).Atmospheric initial conditions influence monthly predictability especially when they already feature the signal of a persisting major mode of atmospheric variability (Reichler and Roads, 2003), such as the Madden-Julian Oscillation, one of the most relevant sources of predictability on the monthly time scale over the Northern extratropics (e.g.Vitart and Molteni, 2010).
At the ISAC Institute of the Italian National Research Council (CNR), a new grid-point atmospheric general circulation model, named GLOBO, has been recently developed (Malguzzi et al., 2011).In the framework of a cooperation with the National Italian Civil Protection Agency, GLOBO is currently operated to produce, once a day, medium-range (up to 6 days) deterministic forecasts and, once a month, monthly probabilistic forecasts of atmospheric anomalies.In this paper, we present the monthly forecasting system (Sects.2 and 3) with the aim of evaluating (Sect.4) the results obtained from an experience of 13 months of forecasting.

Numerical setup and forecasting strategy
The GLOBO model integrates the atmospheric equations on a regular latitude/longitude grid that, for the monthly forecast, is set up with a horizontal grid spacing of 0.75 × 1.0 • and 50 vertical hybrid levels.A full description of the dynamics and physical parameterizations of the model is given in Malguzzi et al. (2011).
The sea surface temperature (SST) evolution is modelled with a climatological mixed layer ocean with depth H.The SST is relaxed to a prescribed field, computed as a sum of the climatological temperature T CLIM and of the initial observed anomaly T AN (the latter is assumed to decay slowly in time in the amount of 15 % per month), as follows: where F NET denotes the sum of turbulent and radiative surface fluxes, ρC the thermal capacity of sea water, and where the constant γ is set to (2.3 days) −1 .The sea ice fraction is computed starting from the observed initial state and applying an observed climatological tendency.Both T CLIM and climatological sea ice cover have been computed from a 20-yr dataset extracted from the ECMWF ERA-Interim reanalysis (Simmons et al., 2007;Dee et al., 2011).The monthly forecasting system proposed here is based on an ensemble technique aimed at forecasting the probability distribution of atmospheric parameters in the extended range.Forecast anomaly of a single variable is computed by removing from the forecast ensemble mean the modelled climatology of the month obtained through reforecast simulations initialized with ERA-Interim reanalyses (see Sect. 3).The initial conditions for the ensemble forecast are defined using the NCEP Global Forecasting System (GFS) analyses, which are available in real time.Specifically, for a given initialization date of the monthly forecast, meteorological fields from fourteen 00:00 UTC and eighteen 12:00 UTC perturbed analysis of the GFS system are used as initial conditions to produce an ensemble of 32 members.
Over land, the temperature and water content of the deepest, climatological layer are specified from the corresponding climatological data computed from ERA-Interim reanalyses, both for forecasts and reforecasts.This is done to reduce the impact on the surface fields due to the use of different datasets for forecasting and calibration (NCEP and ERA-Interim analyses, respectively).
The model is affected by systematic errors or bias that can be a priori determined and removed from forecast anomalies (calibration).The next section describes how we treat model bias.

Model calibration and systematic error
Let B Ψ (λ, ϑ, v, c) indicate the model bias for the generic atmospheric parameter Ψ (geopotential height and temperature on standard pressure levels and daily precipitation).The model bias is a function of longitude λ, latitude ϑ, validity day ν = 0, 1...31, and calendar day c = 1, 2...365.Of course, the bias is zero for ν = 0.
Our goal is to compute B Ψ from a large number of model simulations of past days (reforecasts).Each simulation is initialized on the ERA-Interim reanalysis corresponding to the 1st and the 15th of each month of the 21-yr period ranging from 1 January 1989 to 15 December 2009.We denote by C i (i = 1, 2...24) the calendar days corresponding to the 1st and the 15th of each month (for the sake of clarity, the problem of leap years is not considered).To increase the statistical ground, a two-member ensemble is performed for each C i by starting the model integration from both the 00:00 and 12:00 UTC analysis.This gives a total of 1008 (31-day long) reforecasts that can be averaged to provide model climatology as a function of validity day.
Let Ψ i,y,h (λ, ϑ, v) represent the value of Ψ(λ, ϑ) predicted after ν days of integration by the reforecast started at the calendar day C i of year y = 1, 2...21 and initial hour h = 1, 2. The model error of this particular reforecast is given by the following difference: where Ψ T denotes the verifying ERA-Interim reanalysis of this reforecast.By contracting the previous expression over the index y and h we obtain a first (hereafter referred to as "raw") evaluation of the model bias, which is valid at calendar day C i : Although we may expect that the model bias changes significantly from season to season, it is reasonable to assume that (possible) large differences between two consecutive B Ψ i are due to unpredictable atmospheric fluctuations that do not average out perfectly (given the limited amount of data available here).Therefore, in order to smooth the annual cycle of the model bias, we introduce a Gaussian weighted running average over the calendar days which also allows us to define the model bias for any arbitrary calendar day c as follows: where the normalizing factor W is given by the sum of the Gaussian weights, and where the half-width ∆ is a parameter that can be adjusted to obtain the best results from model calibration.In practice, the raw biases given by Eq. ( 2) with C i close to c give the most significant contribution to the bias at day c (the computation of the distance between two calendar days must take into account the periodic nature of this particular variable).Expressions (3) and (4) give us the basic elements to calibrate the monthly forecasts for any initial date.In the following, ∆ = 45 days is used for 500-hPa geopotential and 850-hPa temperature.
Figure 1 shows the mean model bias (average of Eq. ( 2) over i, y, h) of the 500-hPa geopotential height for the 1st week (average over ν from 1 to 7) and 4th week (ν from 22 to 28) of validity.During the first week, the model already develops the main features characterizing the long-term bias: a positive bias centred over the North Pole, a negative one centred over the South Pole, and a weaker negative bias spread over the intertropical belt and adjacent extratropical low-latitude areas.During the model integration, bias values intensify while some pattern differences arise over the southern portion of the Southern Hemisphere (SH) with, in particular, a sign reversal occurring over the area between Southern America and the Weddel Sea.Moreover, a negative bias over the low-latitude Northern Hemisphere (NH) tem-perate zones extends slightly northward.This implies that our model, in the extended range, reproduces a less intense mid-tropospheric westerly flow than that characterizing the real climate.Furthermore, the model tends to reduce the amplitude of the climatological large-scale wave across the Rocky Mountains.
The amplitude of the spatial pattern of model bias can be quantified by computing its root mean square (RMS) over a spatial domain.Typically, two domains are considered: the NH (latitudes from ϑ 1 = 20 • to ϑ 2 = 90 • ) and the SH (latitudes from ϑ 1 = −90 • to ϑ 2 = −20 • ) extratropics.The RMS is defined by the following expression: The evolution of the model bias as a function of the validity day is shown in Fig. 2 for the NH and SH extratropics.The dashed curves are the average of the RMS of the raw biases (averaged over i from 1 to 24 of R B Ψ i (v)).Over both domains, the raw RMS steadily grows up to about the day 20, where it reaches the asymptotic value of about 40 m.The solid curves represent the same quantity but for the weighted biases (averaged over i of R B Ψ i (v,C i )).Compared to raw RMS, weighted RMS is smaller for all validity days, especially on the SH extratropics.This shows that the weighting procedure described above effectively reduces the estimation of systematic error.This result is in part expected for merely algebraical reasons; in the next section we will show that calibration performed by removing weighted biases effectively improves the model ability to predict anomalies.

Preliminary evaluation of the monthly forecasting system
The preliminary verification of the monthly forecasts is based on the comparison of 500-hPa geopotential height daily fields with ERA-Interim verifying reanalysis.Both dataset are scaled on a 1 × 1 • regular latitude/longitude grid and time averaged 00:00 and 12:00 UTC instants to obtain a single daily field.
The performance of the probabilistic forecasting system is evaluated over 13 monthly outputs (August 2010-August 2011), which are the only months that can be presently verified.Such dataset is neither statistically significant nor allows us to evaluate the forecasting system in a probabilistic framework, which constitutes the proper method to be applied to probabilistic forecasting.However, an overview of the forecast performance can still be drawn by applying non-probabilistic scores to weekly averaged forecast ensemble means.Weekly averaging removes part of the shorterscales, more directly affected by the chaotic atmospheric nature, leaving the larger-scale signal, which a monthly forecasting system is designed for.
The top panels of Fig. 3 shows weekly time series of uncentred anomaly correlation (AC; Wilks, 2006), while the bottom panels of the same figure report the RMS forecast error for the two extratropical areas.All series are obtained by averaging over the 13 cases available.Blue (red) markers refer to calibrated (uncalibrated) forecasts.On both domains, AC gets lower than 0.6 during the second week of integration.The RMS forecast error reaches the value of the climatological RMS (black markers) by the third week.The calibration method presented in the previous section slightly improves the AC up to the third week, with the highest improvement occurring during the first week over the NH extratropics (this is particularly evident for the first days of forecast, which are not resolved in the figure).The RMS forecast error is improved throughout the whole forecast range in both extratropical zones.If calibration were to be performed by subtracting raw biases instead, worse AC would be obtained (particularly for the first forecast days -not shown).This is a clear indication that reforecasting 21 yr in the past is not enough to separate the model systematic and random errors for a given month, so that we must resort to a weighted average of model errors over nearby months to improve AC by calibration.

Summary and conclusions
In this paper, the monthly forecasting system run experimentally at the ISAC institute of the National Research Council of Italy has been presented.Ensemble forecasts of 500-hPa geopotential height, 850-hPa temperature, and precipitation anomalies have been produced on a monthly basis through the GLOBO model.Forecast anomalies have been computed based on reforecast and forecast simulations, the former providing the model climatology used to calibrate the forecast ensemble mean.
Through the comparison with ERA-Interim reanalyses, the large number of reforecast simulations has allowed to evaluate the yearly averaged systematic error of the 500-hPa geopotential height field as simulated by GLOBO.The comparison shows that the main structure of the bias pattern is already developed during the first week of integration, and that the model smooths to some extent the main climatological features of the geopotential height.
The 13 monthly forecasts available so far have been analysed to gain some indications on the forecast error, that has been evaluated based on anomaly correlation (AC) and root mean square (RMS) of forecast error.These (typically deterministic) indices have been computed using weekly averages of modelled and observed atmospheric fields to directly compare statistical quantities (model ensemble means) with the single realization available, through reanalysis, of the real atmosphere.Results show that, even if the AC gets lower than 0.6 by the second week, the RMS forecast error outperforms the climatology in the first two weeks of integration.Moreover, the adopted calibration based on bias correction has been found to reduce the systematic and the forecast RMS errors.Therefore, preliminary results suggest that the application of the ISAC monthly forecasting system provides some predictive skill in the extended range, in terms of weekly means.Further developments of the GLOBO model, and of the forecasting system, are expected to impact the forecast performance in the next future.Possible evolutions rely on implementing a more sophisticated oceanic module (Rendina, 2011) and different initialization data.

Figure 1 .
Figure 1.500-hPa geopotential height mean bias (m) computed from reforecasts and corresponding ERA-Interim reanalyses for the 1st (a) and 4th (b) week of simulation.The zonally averaged bias is shown in side boxes.

Figure 2 .
Figure 2. RMS (m) spatially averaged over the Northern (black curves) and Southern (grey curves) Hemisphere extratropics.Dashed lines refer to raw biases (see text).

Figure 3 .
Figure 3. Anomaly correlation (top panels) and root mean squared (RMS) error (m, bottom panels) computed from 13 monthly forecast anomalies of 500-hPa geopotential height over the Northern (left panels) and Southern (right panels) Hemisphere extratropics.Blue markers are obtained from calibrated anomalies, red markers from uncalibrated (raw) anomalies.The black markers in the bottom panels represent the RMS computed from 1989-2009 ERA-Interim climatology.