Improved forecast of 2015/16 El Niño event in an experimental coupled seasonal ensemble forecasting system

To improve NOAA’s seasonal forecasting capabilities, a new coupled system within the Unified Forecast System (UFS) framework is being developed through a community-wide effort led by NOAA’s Environmental Modeling Center targeting the configuration of a future operational Seasonal Forecast System (SFS v1). An experimental version of this ensemble seasonal forecasting system is tested on forecasting the strong El Niño of 2015/16. The then-operational systems and NCEP real-time seasonal forecasts (CFSv2) underestimated its strength towards the end of 2015 and beginning of 2016. In addition to perturbing the atmospheric initial conditions, run-time stochastic physics-based perturbations are applied in both atmosphere and ocean components of this new coupled system to represent the model uncertainties. The UFS ensembles are initialized on June 1st, 2015 and run through a 9-month period. Compared to CFSv2, the forecast of Niño 3.4 SST and intra-seasonal zonal windstress for the 2015/16 El Niño in the UFS system are improved, as is the ensemble spread. A cold SST forecast error develops in the central equatorial Pacific, likely from excess evaporative cooling, shallower thermocline, and an excessively strong vertical current shear driven cooling. Near the eastern equatorial Pacific coast, on the other hand, warm surface and cool subsurface errors persist from initialization until the end of the forecast. The results suggest that further improvement in the seasonal forecast may be achieved by a combination of factors, including, but not limited to, improving the coupled system initialization, along with the atmospheric physics.


Introduction
The year 2015/16 witnessed one of the strongest El Niño events in recent decades following the 1997/98 episode, that persisted multiple years (Sanchez and Karnauskas 2021;Santoso et al. 2017) and was accompanied by a marine heat wave in the North Pacific (Bond et al. 2015), and drought conditions in the US West coast. Many dynamical and statistical seasonal forecasts of the time failed to predict not only the enormity of the 2015/16 El Niño but also its teleconnections on North American precipitation patterns (Dole et al. 2018;Jong et al. 2018;L'Heureux et al. 2017). Although the mean signal of the winter precipitation across North America following a canonical El Niño event includes dry conditions in the north and wet conditions in the south (Livezey et al. 1997;Ropelewski and Halpert 1989;Trenberth et al. 1998), the models failed to predict the dry winter precipitation anomaly in Southern California in 2016. However, factors going beyond El Niño teleconnections, such as strong SST anomalies in the North Pacific (Bond et al. 2015;Hartmann 2015;Hu et al. 2017;Seager et al. 2015), the phase of the Arctic Oscillation (Lim et al. 2018;Singh et al. 2018), and subseasonal variability unrelated to the El Niño signals (Zhang et al. 2018) have the potential to influence the precipitation over the continental US in 2015/16 too. The failure of many operational systems in early predictions of the enormous 2015/16 El Niño event underscores the need to improve the current prediction systems and include the 1 3 broad uncertainty associated with the internal variability of the system.
Studies show that the timing and location of very strong wind events-westerly wind bursts (WWB) in the westerncentral equatorial Pacific-were critical for initiating the 2015/16 event, and also defined the SST spatial anomaly of the event (Chen et al. 2017;Fedorov et al. 2015;Hong et al. 2017;Hu and Fedorov 2019;Ineson et al. 2018). The relevance of WWBs is further highlighted when the developing El Niño of 2014/15 was halted either due to lack of WWBs or following a strong easterly wind burst (Hu and Fedorov 2016;Ineson et al. 2018). Whether initiating an anomalous equatorial warming leading to an El Niño event or occurring during a growing El Niño thereby sustaining the growth or amplifying it, WWBs has been an integral part of El Niño predictions and predictability (Fedorov et al. 2015;Hu et al. 2014;Puy et al. 2019;Vecchi and Harrison 2000). Due to the semi-stochastic nature of WWBs, modulation of WWBs by equatorial Pacific SST is critical to El Niño dynamics (Eisenman et al. 2005;Gebbie et al. 2007). WWBs can be generated through non-linear interaction with the mean state of the SST (Gebbie et al. 2007;Hayashi and Watanabe 2019;Hu et al. 2014;Tziperman and Yu 2007), and/or coupling between Madden-Julian Oscillation and SST anomaly (Liang et al. 2021). Diverse background ocean conditions would lead to diverse WWB induced changes in SST anomaly response in the equatorial Pacific (Fedorov et al. 2015;Puy et al. 2019). The WWB dependence on background state was demonstrated in forecasting the 2014/15 and 2015/16 El Niños with the seasonal Met Office Glo-Sea5 system and ECMWF System 4 (Ineson et al. 2018). Hence assessing the ocean conditions in seasonal prediction systems is crucial for diagnosing sources of forecast errors in El Niño. WWB-forced oceanic Kelvin waves suppressed the equatorial thermocline, and enhanced the anomalous warming during 2015/16 El Niño event (Abellán et al. 2018;Chen et al. 2017;Gasparin and Roemmich 2017;Hong et al. 2017). The propagating Kelvin waves can also modulate the Equatorial Undercurrent (EUC), as evidenced from Tropical Atmosphere Ocean (TAO) moorings (Roundy and Kiladis 2006), which in turn affects the upper ocean current shear in the eastern equatorial Pacific, where EUC shoals.
Oceanic advection also played a key role in the growth of the 2015/16 El Niño warming (Abellán et al. 2018;Rudnick et al. 2021). A strong freshwater anomaly in the equatorial Pacific was observed during this event driven by both surface forcings and zonal advection (Gasparin and Roemmich 2016). In general, the upper ocean freshening during El Niño events in western equatorial Pacific is driven through formation of a thicker barrier layer that is advected eastward via the zonal current anomaly initiated by WWBs (Gasparin and Roemmich 2017;Wei et al. 2021;Zhang et al. 2021). Biases in barrier layers can lead to mixed layer cooling biases that affect El Niño growth in coupled models (Wei et al. 2021). Salinity stratification was a key reason for the sudden termination of a developing El Niño in the summer of 2014 (Chi et al. 2019;Corbett et al. 2017). A strong westward surface current during the 2015/16 El Niño (Rudnick et al. 2021) drove strong off-shore currents in the eastern equatorial Pacific. The observed weak warming in this region was mostly due to the westward advection via these surface currents that are in turn driven by strong southeasterly wind anomalies along the Peru coast (Paek et al. 2017;Rudnick et al. 2021;Zhong et al. 2019) facilitating coastal cooling.
Climatological biases, especially the excess equatorial cold and dry bias in the central Pacific and warm coastal bias along Peru coast, affect seasonal predictions of El Niño (Guilyardi et al. 2009;Ding et al. 2020;Ray et al. 2018;Vannière et al. 2013;Wu et al. 2022;Bellenger et al. 2014;Graham et al. 2017). Studies have shown that such bias affects CFSv2 forecasts too (Xue et al. 2013;Pillai et al. 2017). In the 2015/16 El Niño, models failed to capture the stronger SST anomalies in the central Pacific and weaker anomalies in the east (Jong et al. 2018;Siler et al. 2017). Although the location of SST anomaly during El Niño varies continuously from event to event across equatorial Pacific (Capotondi et al. 2015;Giese and Ray 2011;Johnson 2013;Ray and Giese 2012), relatively weaker warmings in the east have been linked to an overall drying effect on the US winter, particularly in Southern California (Johnson and Kosaka 2016;Yu and Zou 2013). This indicates that forecasting the location of the anomaly is crucial not only for the El Niño evolution in the equatorial Pacific, but also for its diverse teleconnection patterns.
As a step towards improving NOAA's seasonal forecasting capabilities, a new coupled system within the Unified Forecast System (UFS) framework is being developed through a community-wide effort led by NOAA's Environmental Modeling Center, targeting the configuration of a future operational Seasonal Forecast System (SFS v1). The current study begins with a brief comparison to the Climate Forecast System, version 2 (CFSv2), to illustrate the improvements or degradations in forecasts of SSTs in the Niño index regions and forecasts of other surface fields. Next, the study assesses the experimental ensemble system's prediction of the 2015/16 El Niño event in the equatorial Pacific, with the goal of identifying forecast deficiencies, addressing which may help improve the forecast system. It is encouraging that UFS's ensemble spread encompasses the observations with a reduced ensemble mean forecast error in the equatorial Pacific. Ongoing efforts target further improving the reliability of the UFS ensemble forecasts by bringing the ensemble spread closer to the root-mean-square of the ensemble errors.
The remainder of the paper is organized as follows: Sect. 2 provides a description of the coupled model, ensemble generation and analysis method, Sect. 3 provides a brief description of the reference data sets (observation-based analyses and reanalyses) used to validate the model, Sect. 4 presents the results, and Sect. 5 summarizes the results and discusses the next steps.

Model
The Unified Forecast System (UFS) Coupled Model Prototype 5, a subseasonal forecast system developed by National Centers for Environmental Prediction (NCEP), has been used and documented in earlier studies (Krishnamurthy and Stan 2022;Stefanova et al. 2022). The same system is now extended to 9 months for the current seasonal forecast experiment with added modifications, as follows. The coupled model application consists of an atmospheric component employing the FV3 dynamical core with the CCPP at C192 ( ∼ 50 km) in the horizontal and L64 in vertical resolutions with GFSv15.2 physics package; an ocean component (GFDL MOM6, (Adcroft et al. 2019)) at 0.5 • grid; and a sea ice component (Los Alamos CICE6) at 0.5 • resolution on the same grid as MOM6. Atmospheric initial conditions are from the Global Forecast System version 15 (GFSv15) retrospective analysis. Ocean initial conditions are from the Climate Prediction Center (CPC) ocean data assimilation (Sluka 2018) based on MOM6 at 0.25 • resolution and then interpolated to 0.5 • for the current experiment. This product was generated using a 3D-Var assimilation scheme ingesting insitu temperature and salinity from the World Ocean Database and satellite-derived SST. The 9-month long forecast run was initialized on 1st June 2015, and consists of one control-run and 40 ensemble members.
Among the 80 ensemble members of the uncoupled GFSv15 generating 6-hour forecasts by ensemble Kalman filtering (EnKF), 40 members are selected to generate the initial perturbations for the atmosphere. The GFSv15 analysis field is then added to each of the perturbations to create the atmospheric initial condition for each ensemble member. Two stochastic schemes-Stochastic Perturbed Physics Tendency (SPPT) and Stochastic Kinetic Energy Backscatter (SKEB) are applied in the atmosphere model at each time step as in the uncoupled operational Global Ensemble Forecast System version 12 (GEFSv12) (Guan et al. 2022). However, the magnitude in SPPT is reduced by 30% in amplitude and that in SKEB is increased to 0.7 (instead of 0.6) compared with the GEFSv12 configuration. SPPT scheme perturbs the total parametrized tendencies in temperature, horizontal winds and specific humidity representing the uncertainty associated with the model physical parameterization. It adjusts El Niño forecast errors by reducing excess convective activity in the western Pacific, and increases the ensemble spread to a better calibrated system in ECMWF seasonal forecasting system (Weisheimer et al. 2014). In NCAR-CAM4, SPPT improved the zonal wind variability and distribution of WWBs, which further improved the variability in SST and tropical convection, and broadened the Niño 3.4 power spectrum, reducing the overly periodic nature of ENSO (Christensen et al. 2017). Both the ocean and ice initial conditions remain unperturbed in this UFS experimental seasonal forecast. However separate stochastic perturbations are implemented in the ocean-one on the energetics based Planetary Boundary Layer (ePBL) that perturbs the Total Kinetic Energy (TKE) generation, and the other involves perturbing the potential temperature, salinity and layer thickness.
To assess improvements in the UFS ensemble forecast of 2015/16 El Niño, specifically in the surface fields, the real-time forecast from NCEP operational CFSv2 is used for comparison. Table 1 provides a comparison of specific features in the two systems. The initial conditions for CFSv2 for the atmosphere, land, ocean and sea ice are from the Coupled Data Assimilation System (Saha et al. 2010). In the forecasts the observed greenhouse CO 2 concentration is specified, hence the increasing trend from 1982 is included. The real-time forecasts are initialized at 0000, 0600, 1200, 1800 UTC every day, i.e., 4 members each day. For this study, an ensemble of 40 members was created by taking all members initialized 4 times daily over a span of 10 days between 23rd May-1st June, 2015. This ensemble's forecasts of sea surface temperature (SST), zonal wind stress and sea surface height (SSH) in the equatorial Pacific are assessed and compared to reanalyses and UFS ensemble forecasts.
While a long reforecast period allowing for the calculation of model climatology is available for the CFSv2, no such reforecast has yet been produced for the seasonal UFS, and thus its model climatology is not available. In order to keep the evaluation of the two systems at an equal footing, the observed climatology for the period 2011-2018 is subtracted from each system's forecast to define its anomalies. It is acknowledged that this definition of anomalies (i.e., relative to the observed climatology instead of to that of the corresponding model) means that model biases are playing a role in the resulting assessment.

Verification datasets
A diverse set of observation-based products for the period June 2015-February 2016 are used to validate the seasonal forecasts. Daily SST field is obtained from Operational Sea Surface Temperature and Sea Ice Analysis (OSTIA, (Stark et al. 2007)) that is produced on an operational basis at the UK Met Office using satellite data from sensors that include the Advanced Very High Resolution Radiometer (AVHRR), the Advance Along Track Scanning Radiometer (AATSR), the Spinning Enhanced Visible and Infrared Imager (SEVIRI), the Advanced Microwave Scanning Radiometer-EOS (AMSRE), the Tropical Rainfall Measuring Mission Microwave Imager (TMI), and in situ data from drifting and moored buoys. The daily sea level anomaly obtained from Archiving, Validation and Interpretation of Satellite Oceanographic (AVISO, (Ducet et al. 2000)) altimetry data is used to validate the SSH. The latest climate reanalysis ERA5 provided by the European Centre for Medium-Range Weather Forecasts (ECMWF; (Hersbach et al. 2020)) is used to validate the daily windstress. Daily sea surface salinity (SSS) fields from Seamless System for Prediction and EArth System Research developed Ensemble Coupled Data Assimilation system (SPEAR-ECDA reanalysis; (Lu et al. 2020)) are used. Pentad current fields from OSCAR (Bonjean and Lagerloef 2002) are used to validate sea surface zonal current (SSU). Monthly temperature fields from ECMWF-ORAS5 (Zuo et al. 2018) are used to assess the initial conditions, as well as lead month averaged subsurface temperature forecasts. Thermocline depth is calculated as the depth of the 20 • C isotherm (D20) in forecasts and compared to that in ORAS5. Shortwave and longwave radiative fluxes are obtained from CERES SYN1 • (Doelling et al. 2013), and latent heat flux is obtained from NASA Langley product (Rutan et al. 2015).

Initial conditions
The ocean initial conditions (ICs) in the equatorial Pacific used in the experimental seasonal ensemble forecast from June 2015 are assessed here. Figure 1 shows the equatorial Pacific surface zonal windstress ( x ) anomaly, and thermocline depth (D20), along with the upper 200 m temperature (T) anomaly and zonal current (U) in the ICs and in the corresponding observation-based products for the month of June 2015. The westerly x anomaly in the ICs in western to central equatorial Pacific resembles that in ERA5 although the peak strength is shifted slightly eastward. In contrast to westerly x anomaly decreasing monotonously eastward in ERA5, a secondary peak appears in the ICs around central equatorial Pacific. Further eastward, overly strong easterly x anomalies dominate the equatorial surface in the ICs. The difference in equatorial zonal x affects the slope of the equatorial D20, which is too shallow in the ICs compared to that in ORAS5, particularly between central to eastern regions. The contrast in zonal subsurface temperature anomalies across the equatorial Pacific-cooler in the west, and warmer in the east-is clearly captured in the ICs as in ORAS5, although the warm anomaly near D20 in the central equatorial Pacific is too strong in the IC's. The Equatorial Undercurrent (EUC) and the westward near surface U in the eastern equatorial Pacific are too strong in the ICs, which strengthens the vertical current shear, thus enhancing shear-driven cooling at depths where the EUC shoals along D20 in the eastern equatorial Pacific. The strong easterly x anomalies driven off-shore currents likely enhance upwelling of colder waters along the eastern equatorial coastal region. This is consistent with shoaling of the equatorial D20, and, with overly strong vertical current shear drives a mixing driven cooling of the upper ocean in the region. In the far western equatorial Pacific subsurface, westward U mostly dominates the subsurface above the D20 in the ICs. This is in contrast to that in ORAS5, where eastward U pushes anomalously cooler and fresher waters eastward, thereby freshening the central equatorial subsurface. Given the long memory of the ocean, the effect from the ocean ICs is expected to be retained in the forecasts, at least for the 1st month, and possibly extending to longer lead times.

Ensemble spread and RMSE
The ensemble spread and root mean square error (RMSE) of the ensemble members are calculated as follows: where N = number of spatial grid points, M = number of ensemble members (40), x m,n = forecast member, x n = ensemble mean of forecasts, obs n = observation. An ideal forecast ensemble system is characterized by comparable values among the ensemble spread and RMSE. Using the above formulas the ensemble spread of Niño 3.4 SST grows rapidly with time and beyond the 2 nd forecast month becomes larger than the RMSE (Fig. 2a), which indicates that the uncertainty from the model internal variability is too large. The ensemble spread in the coastal Niño region (Niño 1+2) SST also increases with time but has a larger RMSE than Niño 3.4 SST from the starting month, indicating larger biases coming from ICs. The difference between the ensemble spread and RMSE for the Niño 1+2 SST is smaller than for the Niño 3.4 SST (Fig. 2b). Compared to the SST in Niño 3.4 and Niño 1+2 domain, the Niño 4 averaged x shows a relatively better representation of the uncertainty spread (Fig. 2c). (1)

Forecast of Niño SST anomaly
The UFS forecast of Niño 3.4 SST anomaly shows remarkable improvement over CFSv2 forecasts, with much reduced forecast error in the ensemble mean (Fig. 3a). Ensemble means of both systems start off warmer than observed in the first lead month, but CFSv2 shows too cold forecasts evolving rapidly from July through December, thereby underestimating the strength of the developing El Niño. The UFS ensemble mean also has a cold forecast error during this period, but is of much smaller magnitude than CFSv2 (Fig. 3a). It is worth noting that Niño 3.4 covers the region notorious for excess equatorial cold tongue bias, a common bias among climate models affecting development of El Niño in the models (Guilyardi et al. 2009;Jiang et al. 2020). This bias affects the seasonal forecast of Niño 3.4 SST anomaly in seasonal prediction models (Ding et al. 2020;Vannière et al. 2013;Wu et al. 2022). The cold forecast error in UFS persists from August until December, and decreases from December onward, with a distinct improvement over its single member control-run. The CFSv2 performs better than UFS during the last two forecast months -January and February, 2016. Until December 2015, UFS outperforms the CFSv2 in forecasting the Niño 3.4 SST anomaly, and with a larger ensemble spread encompasses the observed SST, unlike CFSv2. However, further work is needed to optimize the uncertainty spread in the ensemble forecasts (e.g. by tuning the stochastic physics perturbations) to bring the ensemble spread closer to the RMSE of the ensemble mean. Unlike the Niño 3.4 region, the UFS ensemble mean SST anomaly in the coastal Niño (Niño 1+2) region does not show any improvement over its CFSv2 counterpart (Fig. 3b).  (blue), and UFS (red) ensemble forecast systems. The ± ensemble spread is shaded and the ensemble mean is a line for both the CFSv2 (blue) and UFS (red) systems. Also plotted are the UFS control forecast (single member) in dash (red), and the observed SST anomaly from OSTIA (black) in a, b, and the observed x anomaly from ERA5 (black) in c. A boxcar smoothing of 5 days is applied to SST, with no smoothing in x . Observed climatology of SST and x for the period 2011-2018 is used to calculate both the forecast and observed anomalies. This provides a uniform baseline for comparing the three products Although the UFS ensemble mean anomaly starts off overly warm from the 1st month, the warm forecast error grows slowly. By contrast, the CFSv2 ensemble mean starts off relatively closer to observations, but the warm forecast error rapidly increases (and exceeds that of UFS) from July onwards. This indicates that initial bias in the UFS from ICs is potentially persisting during the forecast period in the first few months, gradually growing with longer lead time. However, in spite of the relatively better initial anomaly, CFSv2 does not retain its advantage, as the forecast error rapidly increases until September. The magnitude of the warm forecast error changes in CFSv2, with larger errors in July-August, January-February, and minimal error in October, as opposed to minimal changes in the UFS. However, neither system's ensemble spread clearly encompasses the observed SST anomaly throughout the 9 months forecast period. As in Niño 3.4, the Niño 1+2 observed SST anomaly in CFSv2 is clearly outside its 1 ensemble spread. However it straddles at the lower edge of the uncertainty spread in UFS. Unlike in Niño 3.4, the UFS shows no improvement in Niño 1+2 SST forecast compared to its single member control forecast. The physical processes driving this consistently warmer forecast error in UFS is investigated further at the surface and subsurface in later sections.
The Niño 4 averaged zonal x anomaly in UFS shows remarkable improvement from CFSv2. Located in the western equatorial Pacific, Niño 4 averaged x provides a measure of the wind anomaly that eventually initiates, feeds into the growth, and subsequently initiates decay in the El Niño anomaly across the equatorial Pacific. In June 2015, the UFS forecast x starts off much closer to the reanalysis than does the CFSv2. The UFS ensemble spread throughout most of the forecast period clearly encompasses the reanalysis, ERA5 (Fig. 3c). By contrast, CFSv2 produces an overly strong westerly x anomaly which persists throughout the forecast period except for the period from October through December. UFS ensemble mean captures the low-frequency variation in x closer to the reanalysis, although the highfrequency variations seems underestimated, which likely is an artifact of smoothening in ensemble average. The highfrequency component of x consists of the short-lived strong westerly wind events that are crucial for El Niño growth. Deficiencies in simulating the spatial location and/or timings of these events could affect the El Niño forecasts (Fedorov et al. 2015). Although the reanalysis is generally well within the 1 ensemble spread of the UFS system, occasionally the strong westerly wind events are at the upper edge of the ensemble spread, indicating underestimation of these events in most of the ensemble members. The short-lived peak x anomaly in January 2016 in reanalysis is not captured by the UFS ensemble mean, and probably feeds into the prolonged warming of the equatorial Pacific beyond February 2016 in observations (Santoso et al. 2017).

SST and surface fluxes
An analysis of the temporal evolution of the forecast across the equatorial Pacific facilitates understanding the source of forecast error during El Niño. Warmer SST forecasts start off from the 1st lead month across the equatorial Pacific in UFS, indicating too warm SST from ICs. Warmer forecast errors persist throughout the 9-month forecast duration along the western and eastern coastal equatorial Pacific in the UFS control run, indicating the lasting effect from ICs. In particular, the warm forecast error between 140 • E-160 • E increases between July and October. However these forecast errors are substantially reduced in UFS ensemble mean, but are not completely eliminated (Fig. 4b, c) and are larger than in CFSv2. On the other hand, a large cold forecast error in the central-eastern equatorial Pacific in the CFSv2 is noticeably reduced in UFS (Fig. 4). Unlike UFS, the CFSv2 starts off with a cooler SST in the eastern equatorial Pacific in the 1st lead month. This is of importance as the region is known for excess equatorial cold tongue bias in coupled climate models. The UFS forecasts in the region start off from warmer ICs and do not develop a cold error of comparable magnitude, possibly suggesting a smaller cold tongue bias.
In the western Pacific, alternating instances of excess heating and cooling from net surface heat flux are seen (Fig. 5d). A relatively large positive forecast error in downward shortwave fluxes (Fig. 5a), along with a negative forecast error in downward longwave fluxes (Fig. 5b) and a relatively smaller amplitude positive forecast error in outgoing longwave fluxes (not shown), indicate that during most of the forecast period there is excess radiative heating, which acts in support of the warm SST forecast error. The opposite signs in the downward shortwave and downward longwave fluxes forecast errors indicate errors from cloud cover. For example, less cloud-cover could cause a positive forecast error in downward shortwave fluxes in June, 2015 along with a negative forecast error in downward longwave fluxes. Instances of increased cloud cover (suggested by negative downward shortwave flux error and positive downward longwave flux error) are also seen. The positive forecast error in latent heat flux (Fig. 5c) during most of the forecast duration suggest excess evaporation has the counteracting effect of dampening the SST warming. The sign of the net heat flux error in the western Pacific is largely determined by the relative magnitudes of the generally excessive radiative heating and the generally excessive evaporational cooling described above. Neither error is dominant for prolonged periods, and in the long run the two roughly compensate for each other, consistent with the overall persistence of the initial warm SST bias in the region. In the eastern equatorial coastal regions, the overly warm SSTs at the beginning of the forecast indicate an initial warm error brought in from the ICs too. The latent heat flux error is positive, counteracting the excess heating from the downward shortwave fluxes, thus making little change in the net surface heat flux. Looking closely in the first two weeks of June, a positive forecast error in downward shortwave fluxes briefly changes to negative along with a opposite change in sign of downward longwave fluxes. A warmer SST persists from the ICs as the forecasts start off warmer in June, and the net surface heat flux contributes little to its change throughout the forecast duration. From July onward, a steady positive downward shortwave radiation error, suggests underestimation of cloud cover. Forecast errors consistent with underestimated cloud cover are predominant through the forecast duration, both in the eastern and western coastal regions and indicate the need to improve the parameterization of clouds, radiation and coupled feedback processes.
The relatively weaker cold SST forecast error (compared to CFSv2) in the central-eastern equatorial Pacific increases in July, peaks in September, and then subsides by November (Fig. 4b). This forecast error is consistent with the net surface heat flux error dominated by the latent heat flux component. A positive latent heat flux forecast error (cooling the surface) leads to a cold SST forecast, followed by a negative forecast error (heating the surface) as the cold forecast error decreases (Fig. 5c) in UFS. This suggests excess evaporative cooling driven by either dryness from low level humidity and/or too strong low level winds, both of which are part of the ICs, leading the system to a cold SST forecast. Additionally, a shallower forecast error in equatorial D20 sustains the colder SST forecast error (Fig. 8), as discussed later.

Surface windstress and thermocline depth
Equatorial westerly x anomalies play a key role in exciting bursts of equatorial Kelvin waves that depress the equatorial D20 and amplify the growth of an El Niño event. UFS control forecast shows weaker equatorial westerly x anomaly, which further weakens in the ensemble mean, likely due to averaging the 40 ensemble members (Fig. 6b-d). However, westerly x anomalies in CFSv2 ensemble mean are too strong across the equatorial Pacific (Fig. 6a), with easterly anomalies in the west seen during the peak of the event. Thus UFS does a good job in bringing the anomalies closer to ERA5.
The episodic westerly wind events in UFS around the dateline (crucial for exciting Kelvin waves) are too weakly captured in both the control forecast and ensemble mean. In fact, the wind events are not spatially and temporally coherent to the reanalysis (Fig. 6) as also evidenced in Niño 4 averaged x anomaly (Fig. 3c). These episodic westerly wind events lead to eastward propagating equatorial Kelvin waves that increase the SSH and deepen the equatorial D20 (Figs. 7, 8). In the reanalysis, a strong westerly wind burst occurs around the dateline in October, 2015 and January, 2016 (Fig. 6) and consequently the observed SSH from AVISO increases (Fig. 7) eastward along the path of the propagating Kelvin waves.
Such a chain of events is visible in both the UFS control forecast and ensemble mean, although with a negative forecast error in SSH and D20, mostly due to weak westerly wind events displaced to the east and also inconsistent timings compared to that in reanalysis (Fig. 6). The SSH forecast across the equatorial Pacific is slightly higher (by 5m) in UFS in June, likely from the ICs (Fig. 7). From July onwards the SSH is consistently too high in the western equatorial Pacific but not as much in the east. This indicates a developing error in the east-west surface pressure gradient within the forecast period, which does not relax as much as that observed during El Niño. This would likely influence the forecasts of zonal winds, the equatorial D20 slope and the depth of the EUC core. The improvement in ensemble mean SSH evolution is however evident from the reduced forecast error compared to the control run (Fig. 7b,c). When compared to the evolution of SSH forecast error in CFSv2 (Fig. 7a), although the spatial and temporal pattern of the forecast error remains unaltered, the amplitude of the error is substantially reduced in UFS (Fig. 7b).
Comparing the equatorial D20 anomaly forecast across the equatorial Pacific clearly shows a delay in D20 deepening in the east in UFS compared to that in ORAS5 (Fig. 8). This results in a shallower forecast error in the central to eastern equatorial Pacific in July through October, which facilitates a cooler SST forecast error (Fig. 3). Although the anomalous D20 deepening in UFS starts later than in the reanalysis, it persists longer, resulting in a deeper forecast error in January and February, 2016-which contributes towards the warmer SST forecast error (Fig. 3) in the far eastern equatorial Pacific.

SSS and zonal current
During anomalous warming in the equatorial Pacific, surface freshening occurs due to both precipitation anomaly moving eastward and fresher waters being pushed eastward, thereby creating a negative salinity anomaly in the central Pacific. This is evident in observed SSS from Argo array (Gasparin and Roemmich 2016), as also in the observation based daily estimation of SPEAR-ECDA reanalysis (Fig. 9). Anomalous precipitation in the UFS, calculated using Tropical Rainfall Measuring Mission (TRMM) climatology, indeed shows an eastward shift across the equatorial Pacific (not shown). However, the freshening in the central equatorial Pacific is too weak throughout the forecast period in both the UFS control run and ensemble mean. At the same time weaker eastward surface currents in both the control run and ensemble mean relative to that in OSCAR datasets (Fig. 9) indicate weaker advection of fresher waters eastward. In the east the westward surface currents are too strong seasonally in the ensemble mean, which likely affects coastal upwelling too.

Equatorial subsurface anomaly
The magnitudes of subsurface temperature forecast errors are larger than at the surface, and, due to the long memory of the ocean, biases from ICs likely persist into the forecasts. Warm anomalies in the ICs across the equatorial Pacific (Fig. 2), with a maximum at D20 depth grow too slowly in the UFS forecasts (Fig. 10). Warmer peak anomalies at D20 depth in ORAS5 rapidly grow in the next few months. The forecasts of both anomalous temperature and EUC in UFS are much weaker than those in ORAS5.
Specifically, in ORAS5, anomalies greater than 2 • C outcrop to the surface in the eastern Pacific by July 2015, and spread westward, warming the central-eastern equatorial subsurface by more than 2 • C in October 2015. However, in the UFS ensemble mean, a delay in westward spreading of surface warm anomalies is seen, as anomalies greater than 2 • C outcrop to the surface only in November (6th forecast month). This likely explains the warmer surface in the eastern coastal regions in the UFS ensemble mean, compared to that in reanalyses in December 2015. And thereby serves a reason for differences in peak El Niño spatial warm anomaly in UFS, with too strong anomaly along the coastal waters. A weaker subsurface warming in the central-eastern equatorial Pacific until November 2015 is also consistent with a shallow forecast error in D20 (Fig. 8). The EUC core along D20 in UFS shoals eastward in the ensemble mean as in ORAS5, but is too strong in the eastern equatorial subsurface closer to coast at longer lead months (Fig. 10). EUC strengthening is delayed by a month in the UFS, however currents greater than 0.2 m/s span a wider range of depths than in ORAS5 (Fig. 10). By September, when EUC in ORAS5 shoals, almost reaching the surface east of 105 • W , the EUC forecast in UFS is still deeper, with westward currents dominating near the surface.
The too strong westward surface currents in the UFS ensemble mean in eastern equatorial coastal regions are Fig. 8 Hovmoller diagram of equatorial averaged ( 2 • S-2 • N ) monthly thermocline depth (D20) anomaly across equatorial Pacific in a UFS ensemble mean; b corresponding field from ECMWF-ORAS5 and c forecast error in UFS ensemble mean Fig. 9 Hovmoller diagram of daily sea surface salinity (SSS) anomaly (a-c) and surface zonal current (SSU) (d-f) across the equatorial Pacific in a, d UFS control forecast, b, e UFS ensemble forecast mean and c, f concurrent fields in verification dataset. SSS is averaged in 5 • S-5 • N across equatorial Pacific to reflect the wide meridional freshening and SSU is averaged in 2 • S-2 • N . Verification dataset for SSS is daily fields in SPEAR-ECDA reanalysis, and for SSU is pentads from OSCAR

Fig. 10
Equatorial Pacific subsurface temperature ( • C ) forecast anomaly in shading, overlaid by the zonal current (m/s) forecast in contour (black) and D20 (m, blue line) forecast from UFS ensemble mean (a, c, e) and the concurrent fields from ECMWF-ORAS5 reanalysis (b, d, f). The ensemble mean averaged for the a 1st, c 4th, and e 7th month leads, and b June, d September, and f December monthly fields from ECMWF-ORAS5 are shown 1 3 coming from the ICs, as evident from the vertical U-profiles in regions close to the coast (Fig. 11). This likely enhances off-shore movement of water and increases coastal upwelling of colder waters, which could contribute to the delayed anomalous warming as discussed above (Fig. 10). From the 1st month forecast, overly strong westward currents (U in Fig. 11) dominate the upper 40 m. The weak near-surface westward U in ORAS5 and the overly strong near surface westward U in UFS from June remain almost unaltered during the forecast period. A too warm upper 15 m temperature is counteracted by a too shallow D20 and too strong upper ocean stratification extending to 50 m (dT/dz in Fig. 11) relative to ORAS5. Thus, in layers shallower than 20 m, the effect from ICs likely sustains the warm SST forecast error (see Sect. 4.4.1). This warming at the surface could not penetrate through deeper layers due to too strong stratification, resulting in the too warm surface persisting along the eastern coastal region at all lead months (Fig. 11).
At depths below the surface layer (between 20 and 100 m), cooler forecasts are seen in the eastern equatorial region (T in Fig. 11), relative to ORAS5. The upper ocean stratification forecast is too strong in the upper 40 m but weakens slightly at depths below during the entire forecast period. This likely enhances the strong mixing-driven cooling in subsurface temperatures. The eastward EUC forecast at depth strengthens more than in the reanalysis during the Fig. 11 Vertical profiles of temperature (T, top), stratification (dT/dz, middle), and zonal current (U, bottom) forecasts averaged over 90 • W :80 • W for the 1st (left), 4th (center) and 7th (right) lead months from UFS ensemble mean (red) and the concurrent profiles from ECMWF-ORAS5 (black). Meridional averages of T and dT/dz span over 5 • S :5 • N , whereas U spans over 2 • S:2 • N to reflect the core of the EUC closer to the equator. Also shown in each plot the ± ensemble spread from the UFS ensemble forecasts is shaded red. The peak magnitude in dT/ dz represents the depth of D20, and the maximum U represents the core of EUC forecast period and likely comes from ICs. This strengthens the vertical current shear, as the stronger near surface westward U strengthens in UFS forecasts relative to ORAS5. This likely enhances the shear driven cooling of the subsurface too.
The ensemble spread of T, dT/dz and U forecasts in the eastern equatorial Pacific increases slightly with depth (Fig. 11). However, the spread of T and dT/dz forecasts barely encompasses ORAS5 fields until the 7th month, indicating that biases from ICs likely dominate the forecasts even at longer leads. Ensemble spreads of T and U forecasts increase with time and encompass the observed T and U profiles by the 7th lead month. However such improvement in forecasts of dT/dz is not apparent in the upper 80 m, although the reanalysis is well-centered in the ensembles below 80 m from 7th lead month onwards. This likely indicates that biases in dT/dz from ICs are more likely to persist and need to be addressed by improving the ocean initialization.
Moving further westward away from the coastal effect into the region where cold SST forecast error develops, between 140 • W-100 • W , the ICs do not show large biases. Instead, forecast errors begin appearing from the 3rd forecast month (August 2015) onwards. As shown in Fig. 12, T, dT/dz and U are well simulated relative to the reanalysis in the 1st forecast month. However by the 4th forecast Fig. 12 As in Fig. 11, but for regions averaged over 140 • W :100 • W . The deepening of D20 (dT/dz peak) and EUC core (maximum U) are well captured in the forecasts month the ensemble mean forecasts are too cold, driven by too strong dT/dz, shallower D20 and too strong vertical current shear driven upper ocean cooling (Fig. 12). The near surface westward currents are too strong in UFS thereby excessively strengthening the vertical current shear. Thus, although excess evaporative cooling at the surface triggers the cold SST forecast error in the UFS ensemble mean, a shallower D20 and upper ocean processes further contributes to the cold forecast error. The ensemble spread is too narrow for T, dT/dz, and U in the 1st forecast month, although it increases with forecast leads (Fig. 12). By the 4th forecast month the ensembles are too cold in the subsurface, with the reanalysis T straddling the upper edge of the ensemble spread. By the time the cold SST forecast error starts to subside (Fig. 4) in the 5th forecast month (October), D20 anomalously deepens and EUC strengthens excessively (not shown) reducing the vertical current shear in ensemble mean, thereby decreasing the upper ocean cooling. Although the forecast of upper ocean stratification is still reasonably captured, as evident from the forecast spread (dT/dz) encompassing the reanalysis dT/dz, the EUC is too poorly simulated (overestimated), as the reanalysis EUC straddles the lower edge of the forecast spread in U. However, the deeper D20 and reduced evaporating cooling works in tandem to further decrease the developing cold SST forecast error from October onwards in the ensemble mean.

Summary and discussion
An assessment of the 2015/16 El Niño forecasts in a seasonal ensemble forecasting system currently being developed within NOAA's UFS, targeting NOAA's future SFS v1 is presented here. In this experiment the ensemble is generated by perturbing the atmospheric initial conditions only, in an ocean-ice-atmosphere coupled modeling system initialized on the 1st June, 2015. The addition of run-time stochastic physics-based perturbations in the atmosphere and ocean further adds sampling variability in the generation of 40 ensemble members.
Compared to NOAA's current operational seasonal forecasting system, CFSv2, the experimental UFS system shows promising progress in forecasting the 2015 El Niño event, with prominent improvement in SST, zonal windstress, and SSH evolution in the equatorial Pacific. Although a definitive attribution of the improvement in this experimental UFS over the CFSv2 forecasts is beyond the scope of the study, an attempt has been made to identify possible contributing factors. In particular, the Niño 3.4 SST anomaly shows much reduced forecast error in UFS compared to that in the CFSv2 system. With a broader ensemble spread encompassing the observed SST, the UFS significantly improves the reliability of the Niño 3.4 SST seasonal forecast. The low-frequency variation of Niño 4 x forecasts is also closer to ERA5, although the high-frequency variation is underestimated. However, the east-west surface pressure gradient (as indicated by SSH) does not relax as much as that observed during El Niño, which likely affects the x , equatorial D20 slope and the depth of the EUC core.
A cold SST forecast error in the central equatorial Pacific develops and decays within the forecast period. It is likely triggered from excess evaporative cooling and sustained by a strengthened upper ocean stratification, shallower D20 due to weaker downwelling Kelvin waves and a strengthened vertical current shear-driven cooling due to strengthening of near surface westward currents. The gradual deepening of D20, weakening of upper ocean stratification, and strengthening of EUC around November 2015, likely reduce the cold SST forecast error, assisted by a reduced evaporative cooling. This is evident from the cooler forecast error of the Niño 3.4 SST anomaly as well. However, overly strong easterly x persisting from initialization and throughout the forecast in the eastern equatorial Pacific closer to the coast (Fig. 1) strengthen the westward near-surface currents and likely sustain the too shallow D20 forecast in the region in June 2015. Additionally, a shallower D20 forecast in the central equatorial Pacific, driven by weaker downwelling Kelvin waves, spreads eastward, which could likely delay the subsurface warming in the eastern coastal region, and the consequent outcropping and westward spreading of surface anomalies (warmer than 2 • C ). This delay leads to a deeper D20 forecast error in the eastern coastal region, late in January 2016, and likely feeds into the warm coastal forecast error in January 2016 (Fig. 8). Despite the strong westward near-surface currents (even stronger than the verifying reanalysis; see Fig. 11), and shallower D20 forecasts (between June 2015 and December 2015), the coastal SST forecast is excessively warm (Sect. 4.4). The overly warm initial SSTs persists to contribute in sustaining the overly warm SST for the remainder of the forecast. Reduced cloud cover and excess warming from surface fluxes appear to play minimal role in overall SST forecast evolution, due to the counteracting effect from excess evaporative cooling, yet needs to be corrected for improving the forecast system. Taken together, these findings indicate that improvements are needed in the SST, surface winds, ocean currents, temperature profiles in the ICs, along with improvements in the model's atmospheric physics, for the central-eastern equatorial Pacific region.
As illustrated in Fig. 3, the ensemble spread in the UFS ensemble system encompasses the observed SSTs and x in the Niño regions better than that in CFSv2. A possible contributing factor for this is the introduction of stochastic perturbations in the UFS ensemble system, which may lead to a better representation of WWBs in the western equatorial Pacific. Beside triggering downwelling eastward Kelvin waves that plays a key role in El Niño development (Chi 1 3 et al. 2019;Corbett et al. 2017;Wei et al. 2021), WWBs initiate a zonal current anomaly that pushes fresher water eastward (Gasparin and Roemmich 2017;Wei et al. 2021;Zhang et al. 2021). A delay in downwelling Kelvin waves likely played a role in a delayed D20 deepening in the eastern equatorial Pacific in UFS ensemble mean (Fig. 8). The freshening anomalies are widely dispersed eastward in UFS ensemble mean, although concentrated around the dateline in the UFS control run, as in SPEAR-ECDA reanalysis. The location of the maximum freshening anomalies is likely a concern in coupled models, as El Niño intensity has been suggested to be strongest when the location is around 170 • W (Guan et al. 2022). Thus, among other factors, temporal and spatial errors in WWBs in the UFS ensembles may be inhibiting the eastward advection of fresher anomalies in the western equatorial Pacific.
Including WWB parameterization scheme (Gebbie et al. 2007) in an ensemble system (CESM) improved the occurrence and amplitude of WWBs due to enhanced Bjerknes feedback, and along with an improved magnitude of the surface wind stress skill in El Niño prediction enhanced too (Tan et al. 2020). In this experimental study, the x in the UFS ensemble mean is nicely improved compared to CFSv2. One of the factors contributing to this improvement may be the inclusion of the SPPT scheme. SPPT is known to improve the temperature, horizontal winds and humidity fields that are stochastically perturbed at every timestep, leading to simulating large WWBs and triggering El Niño if the timing is right (Christensen et al. 2017). Further improvement in the occurrence, amplitude and location of WWBs in UFS would likely improve the spatial location of the El Niño anomaly in equatorial Pacific. Adopting a stochastic cellular automata (CA) parametrization in the UFS-based configuration of the Global Forecast System's (GFSv15.1) deep cumulus convection scheme (Bengtsson et al. 2021) organized the convective precipitation in space and time in the tropical Pacific, and provided a more systematic impact on simulated convectively coupled equatorial waves. It is worth testing whether implementing this parametrization also has the potential to capture a wider range of uncertainty associated with x variability. This could further improve the timing and spatial occurrence of WWBs in the ensemble, that is necessary for modulating the amplitude and location of a developing El Niño. However, CA alone cannot provide sufficient uncertainty spread from cumulus convection, as can be given by SPPT (Bengtsson et al. 2021). Thus implementing the CA in addition to SPPT in UFS is a possible next step to improve the triggering and location of WWBs in equatorial Pacific.
Although UFS ensemble forecasts show promising improvement, there are a few shortcomings in the forecasts. One of the prime concerns is the overly warm eastern coastal forecasts, including the overly warm Niño 1+2 SST in all the ensemble members likely from ICs as mentioned earlier. An updated data assimilation system under active development will likely address this concern. Another cause of concern is in regards to the ensemble spread in Niño 3.4 SST forecast, which is too wide, although the model internal variability encompasses the observed SST. In order to optimize the forecast uncertainty, the ensemble spread needs to be consistent with the RMSE, and this is an active area of research in UFS seasonal experiments. Experiments on improving physics parameterizations for the coupled system and for optimizing ensemble configurations are ongoing. Preliminary experiments in UFS with stochastic physics perturbation introduced only in the ocean component result in reduced ensemble spread, although still broader than optimal. Further experiments with tuning of the perturbations (including CA and SPPT), and inclusion of ocean perturbations are some of the directions being explored.
A caveat of this study is that observed climatology is used instead of the model climatology in defining the forecast anomalies. In other words, model biases are included in the forecast error. This likely amplifies the forecast error in CFSv2, which has a documented (Xue et al. 2013) cold summer SST bias in the equatorial Pacific. Model biases are similarly part of the analyzed forecast error in the UFS ensemble, although the sign and magnitude of these biases is not known since this early experimental system does not yet have a model climatology established. Additional runs, sufficient to calculate a model climatology for the seasonal UFS, are needed before bias-corrected errors can be compared between the two systems.

Conclusion
An experimental seasonal forecast system developed within the Unified Forecast System framework is demonstrated to show a reasonably improved forecast of the 2015/16 El Niño event in equatorial Pacific compared to the current operational forecast system. The failure of many operational systems in early predictions of the enormous 2015/16 El Niño event underscored the need to improve the then-current prediction systems and include the broad uncertainty associated with the internal variability of the system. With a reduced cold forecast error in Niño 3.4 SST anomaly, the ensemble spread from the experimental UFS set-up encompasses the observed SST throughout the 9 month forecast duration. Errors from initial conditions appear to persist in the eastern equatorial coastal regions, as illustrated by the persistent warm errors in the coastal Niño 1+2 SST forecast. The effect from excess evaporative cooling counteracts the reduced cloud cover and excess warming from surface fluxes, revealing canceling errors. The cold SST forecast error in the central-eastern equatorial Pacific on the other hand appears to be triggered by excess evaporative cooling and sustained by upper ocean cooling at thermocline depth. The temporal and spatial discrepancies in weaker westerly wind bursts in the western equatorial Pacific lead to too shallow equatorial thermocline, likely from weaker downwelling equatorial Kelvin waves. Thus better initialization of the ocean and atmospheric model components needs to go hand in hand with improvements in parameterization of clouds, radiation and coupled feedback processes for reducing errors in El Niño forecasts. Also, further work is needed to obtain an optimal ensemble spread by focusing on optimizing the application of stochastic physics in seasonal ensemble forecasts. This experiment with a single start date lays the foundation for further coupled-model development along with improving the seasonal forecast system towards achieving the target of operational Seasonal Forecasting System version 1 under the realm of the Unified Forecast System in the next few years.