Mechanisms driving ESM-based marine ecosystem predictive skill on the east African coast

The extension of seasonal to interannual prediction of the physical climate system to include the marine ecosystem has a great potential to inform marine resource management strategies. Along the east coast of Africa, recent findings suggest that skillful Earth system model (ESM)-based chlorophyll predictions may enable anticipation of fisheries fluctuations. The mechanisms underlying skillful chlorophyll predictions, however, were not identified, eroding confidence in potential adaptive management steps. This study demonstrates that skillful chlorophyll predictions up to two years in advance arise from the successful simulation of westward-propagating off-equatorial Rossby waves in the Indian ocean. Upwelling associated with these waves supplies nutrients to the surface layer for the large coastal areas by generating north- and southward propagating waves at the east African coast. Further analysis shows that the off-equatorial Rossby wave is initially excited by wind stress forcing caused by El Niño/Southern Oscillation-Indian Ocean teleconnections.


Introduction
Marine ecosystems and the living resources they sustain are subject to pronounced climate-driven fluctuations that have long challenged sustainable marine resource management efforts (Finney et al 2010, Costello et al 2016, Tommasi et al 2017. These fluctuations are compounded by ocean warming, changing ocean productivity baselines, deoxygenation, and acidification, all of which potentially influence the distribution of marine habitat types, phenology, and functioning of marine ecosystems (Doney et al 2012, Cheung et al 2013, Glibert et al 2014, Dussin et al 2019, Kwiatkowski et al 2020. Although century-scale climate change assessments provide a scientific foundation for developing policy options, such long-term projections are not sufficient for the tactical management decisions on seasonal to multiannual time-scales critical to long-term resilience. Thus, reliable tools to anticipate seasonal to multiannual and longer-term changes in the interacting physical, chemical, and biological processes are both required for holistic marine resource management. Earth system models (ESMs), the most comprehensive climate models incorporating Earth's biogeochemical cycles, have the potential to predict phenomena emerging from diverse physical and biogeochemical processes, including marine ecosystem responses to climate variability (Watanabe et al 2011, Dufresne et al 2013, Lindsay et al 2014, Bonan and Doney 2018. Several recent studies have suggested that many critical drivers of marine ecosystems, including temperature, nutrients, oxygen, acidity and ocean productivity, may be predictable on seasonal to multiannual time horizons (Seferian et al 2014, Taboada et al 2019, Frolicher et al 2020. This promise has been verified through the recent integration of biogeochemical dynamics with physical climate assimilation and prediction systems (Park et al 2019). Extensive retrospective forecast experiments using an ESM revealed verifiable seasonal to multi-annual chlorophyll predictions in many ocean areas, and promising relationships with fisheries fluctuations in several large coastal areas.
One of the strongest linkages between skillful chlorophyll predictions and reconstructed fisheries catch was found along the east coast of Africa (Park et al 2019). Skillful multi-year chlorophyll predictions within two large coastal areas spanning the east African coast were highly correlated with fluctuations of small pelagic fisheries up to three years in advance (after allowing for a one year lag for ocean productivity signals to manifest in the fishable stock). The mechanisms underlying skillful chlorophyll predictions, however, were not explored, eroding confidence in predictions in any adaptive management measures for the significant industrial and subsistence fisheries that they support. The present study addresses this by conducting a detailed analysis of prediction skill across 20 years of multi-annual reforecast experiments initialized every month. Given that the countries of these regions have started to develop a transboundary monitoring and assessment program to deal with marine environmental issues and regional management approaches (Ménard et al 2007, Hutchings et al 2009, Vousden 2016, the current study may provide a timely contribution to inform these efforts by characterizing the dynamical mechanisms underlying prediction skill of seasonal to multiyear chlorophyll variability.

Global marine biogeochemical prediction system
The global marine biogeochemical prediction system used in the present study is based on the ESM developed at the Geophysical Fluid Dynamics Laboratory (GFDL). The ESM is a fully coupled atmosphere-land-ocean-sea ice model integrated with the GFDL's marine ecosystem model, the Carbon, Ocean Biogeochemistry and Lower Trophics (COBALTs) (Stock et al 2014a(Stock et al , 2014b. COBALT, an intermediate complex biogeochemical model, simulates 33 tracers to resolve global-scale cycles of nitrogen, carbon, phosphorus, iron, oxygen, and silica with three explicit phytoplankton and three explicit zooplankton groups. The horizontal resolution of the ESM is 2.5 • longitude × 2 • latitude for the atmosphere and land, and 1 • × 1 • for the ocean, sea ice, and marine biogeochemistry. Each grid has 24 hybrid sigma/pressure vertical layers for the atmosphere and 50 vertical layers for the ocean. The initial condition for the marine biogeochemical prediction system is produced from GFDL's data assimilation system, the ensemble-coupled data assimilation (ECDA) system (Chang et al 2013). The ECDA system assimilates both observed atmosphere and ocean states, including winds and temperature from the National Centers for Environmental Prediction-Department of Energy Reanalysis 2 (Kanamitsu et al 2002), oceanic profiles from the World Ocean Database and Argo profiles, and the sea surface data from the National Oceanic and Atmospheric Administration's optimum interpolation SST v2 high resolution dataset (Reynolds et al 2007). The system employs an ensemble Kalman filter assimilation scheme and is integrated with the marine ecological model, COBALT, to reproduce historical ocean biogeochemical fields. The modeled oceanic and atmospheric fields are optimally constrained by observations to suppress the degradation of biogeochemical fields caused by momentum imbalances from equatorial data assimilation (Park et al 2018).

Retrospective forecast experiment and skill assessment
The retrospective forecasts are the same as the ones used in a recent study (Park et al 2019). The forecast experiments are initialized at the 1st day of every calendar month during 1991-2017 and run for two years with 12 ensemble members (supplementary figure 1). The initial conditions for this experiment are produced by data assimilative hindcasts from the ECDA-COBALT system as described above. A previous work has already shown that the global ocean biogeochemical variables produced by this assimilation system can capture well observed large-scale biogeochemical patterns and variability (Park et al 2018), providing the appropriate initial conditions for the biogeochemical prediction system. The anomalies of predicted biogeochemical variables are calculated relative to the lead-dependent climatology from the 27 year ensemble mean predictions for each initialization month, given that the full-field assimilation method used here generally leads to model drift toward its own preferred state once the prediction starts. All forecast data analyzed in this study are the ensemble mean of predicted variables.
Predicted chlorophyll anomalies are compared with satellite-retrieved chlorophyll concentrations. We used the satellite chlorophyll data from the two ocean color sensors, the Sea-viewing Wide Field-ofview Sensor and the Moderate Resolution Imaging Spectroradiometer (Esaias et al 1998, McClain et al 1998. The original chlorophyll data was binned to a grid of 9 km × 9 km, thus the data has been regridded onto a 1.0 • × 1.0 • grid using a bi-linear interpolation method for computational efficiency. The median value of chlorophyll in each 1.0 • grid is used in the interpolation process given the nearly lognormal distribution of ocean chlorophyll concentrations (Campbell 1995). Chlorophyll prediction skill is measured by the temporal anomaly correlation coefficient between the predicted and satellite-retrieved chlorophyll concentrations after taking logarithmic transformation. The significance test for the correlation skill uses the reduced number of effective degrees of freedom that is defined by the autocorrelations of a time-varying field, providing a conservative significance threshold (Bretherton et al 1999).
Chlorophyll prediction skill is evaluated at the large marine ecosystems (LMEs) scale, by focusing on the two east African coast LMEs, the Agulhas Current and the Somali Current. The LMEs are ocean areas defined by ecological criteria including bathymetry, hydrography, productivity, and trophodynamics (Sherman 1991). These large coastal areas encompass regions from estuaries to the seaward boundaries of the continental shelf or of coastal current systems, in which the adaptive management strategy to environmental changes is important to facilitate sustainable marine resource utilization given their substantial socio-economic benefits (Sherman 2014, Tommasi et al 2017.

Chlorophyll prediction skill in coastal LMEs
LME-scale marine biogeochemical prediction from the global marine biogeochemical prediction system is assessed by anomaly correlation coefficient between the predicted and satellite-retrieved chlorophyll during the period 1998-2018. Chlorophyll anomalies are averaged in each LME region after removing the climatological annual cycle. The correlation coefficient is calculated for each initialization month and lead time (see section 2).
The prediction system shows notable chlorophyll forecasting skill in the two African coastal LMEs, Agulhas current and Somali coastal current systems (figure 1; geographical locations of both LMEs are shown in supplementary figure 2). In the Agulhas Current system, significant chlorophyll prediction skill appears up to two years, with the alternating pattern of high-and low-prediction skills in diagonal bands corresponding to high predictability seasons (figure 1(a)). The skillful predictions are mostly possible in the austral winter and early spring regardless of initialization month, whereas in other seasons chlorophyll predictions are not possible, particularly for austral summer forecasts. The highest chlorophyll prediction skill in the austral winter and early spring correspond to the period of highest chlorophyll concentrations during which mixed layers are deepest in this predominantly nutrient-limited subtropical system (supplementary figure 3(a)). The high prediction skill in austral winter and its reappearance in the following winter are similar to a prediction pattern of extratropical sea surface temperature and chlorophyll, in which a subsurface signal remains and reemerges when the mixed layer deepens by strong seasonal wind (Alexander et al 1999, Stock et al 2015, Park et al 2019. Indeed, the temporal evolution of subsurface chlorophyll anomalies regressed onto surface chlorophyll shows that significant chlorophyll signals remain evident beneath the mixed layer during summer and reemerge when the mixed layer deepens during the subsequent fall and winter (supplementary figure 4).
In the Somali coast, significant chlorophyll prediction skill is generally limited up to lead time of 1.5 years ( figure 1(b)). Unlike Agulhas region, predictions of both austral summer and winter are possible within short lead time, and skill is relatively high for predictions of late austral summer and early fall (January-April). This coincides with two periods of peak chlorophyll concentration (i.e. austral summer and winter) in this region, in which primary productivity is controlled by wind-driven mixing and consequent nutrient entrainment from deeper water (supplementary figure 3(b)) (Veldhuis et al 1997). Note that detrending the data has little effect on the chlorophyll prediction skill in the two LMEs (supplementary figure 5).
The skill of chlorophyll predictions can also be seen in the interannual time series of satellite chlorophyll anomalies in Agulhas and Somali LMEs, where observed patterns are well predicted from Januaryinitialized prediction at longer leads, i.e. 21 months and 20 months leads, respectively (figures 1(c) and (d)). This temporal chlorophyll variability is found not to be dominated by a certain area of the LMEs, given strong covarying relationships of chlorophyll within the LME systems (supplementary figure 6). The long-lead prediction is largely dominated by the high chlorophyll in 2002-2003 that corresponds with increased small pelagic fish catches during this period (Hutchings et al 2009), implying the potential utility for marine resource application as seen in a recent study (Park et al 2019). The skillful LME-scale prediction is particularly encouraging, given the coarse resolution of ocean grids in this global ESM that limits the simulation of coastal circulation and ecosystem processes. This indicates that LME-scale chlorophyll variability in these regions is substantially controlled by large-scale physical and biogeochemical dynamics at least for the two year prediction horizon.

Mechanisms of chlorophyll prediction skill
To understand the dynamics underlying chlorophyll prediction skill, the temporal evolutions of upper ocean (0-200 m) nitrate anomaly patterns associated with chlorophyll variations in the two LMEs are examined (figure 2). Macronutrients such as nitrate, phosphate and silicate have been found to have an important influence on phytoplankton growths in these systems Codispoti 1980, Barlow et al 2020), thus we analyzed nitrate as a proxy for nutrient variability associated with chlorophyll prediction skill (supplementary figure 7). The basinscale patterns of predicted nitrate anomalies regressed onto Agulhas satellite chlorophyll show that significant positive nitrate anomalies propagate westward from the central Indian Ocean. The positive predicted nitrate anomalies centered around 80 • E, 20 months ahead of the September forecast window, propagate westward as the forecast window is approached (figures 2(a)-(c)), and then spread both north and south along the western boundary (figures 2(d) and (e)). This eventually increases upper ocean nitrate and chlorophyll in the Agulhas system, providing a key source of chlorophyll prediction skill.
Large-scale evolution of nitrate anomalies linked to Somali chlorophyll prediction exhibit similar propagation behavior (figures 2(f)-(j)). The zonally elongated positive nitrate anomaly appears in the central Indian Ocean at 19 months lag, and the center of nitrate anomaly propagates westward and then spread to Somali system during the 19 months. The concurrent connection between nitrate and chlorophyll anomalies in the Somali LME is lower than that in the Agulhas, which is potentially a reflection of the coastally-restricted Somali system and the prominence of monsoonal wind, and river discharge forcing which may partially confound predictable signals arriving from the ocean basin (Halpern andWoiceshyn 2001, Mutia et al 2021).
The westward moving nitrate signal across the Indian Ocean represents the off-equatorial oceanic Rossby wave slowly propagating from the eastern Indian Ocean to the west. This wave signal found in the biogeochemical variable is consistent with coupled Rossby waves observed in physical variables, such as sea level or upper ocean heat contents in the Indian Ocean (White 2000, Jury andHuang 2004). The slow phase speed of Rossby waves in this off-equatorial region is due to the inverse relationship between phase speed and latitude, thus it takes 2-3 years to cross the Indian Ocean at 15 • S (Perigaud and Delecluse 1993). Consistent with this, the propagation speed of nitrate anomalies found in figure 2 is approximately 0.08 m s −1 .
The Rossby wave signal shown in the nitrate anomalies in the Indian Ocean implies that this physical oceanic phenomenon provides a basis for predicting ocean biogeochemical variables in the two African coast LMEs. The upwelling Rossby wave supplies nutrient to the euphotic layer, and its slowly moving signal together with reemergence of subsurface nutrient anomalies during periods of peak mixing and chlorophyll in these tropical and subtropical LMEs allows the prediction of coastal chlorophyll at longer leads. Skill diminishes when surface chlorophyll is suppressed by strong stratification over the austral summer and adjacent months.

Rossby wave initiated by El Niño/Southern Oscillation (ENSO)-Indian Ocean teleconnection
What mechanisms are associated with the initiation of upwelling Rossby wave and associated nitrate anomalies in the Indian Ocean? To identify the physical process, the analysis for the temporal evolution of largescale dynamics is extended up to four years before the increase in Agulhas chlorophyll anomalies. For this, we used the reconstructed physical and biogeochemical fields from the ensemble-coupled data assimilation system integrating with the global marine biogeochemical model (i.e. ECDA-COBALT), which has been used for the initialization of our two yearlong prediction run (Park et al 2019).
The zonal propagation of nitrate anomalies in the central Indian Ocean are well depicted in Hovmöller diagram of nitrate anomalies averaged over the latitudinal band 20 • -10 • S ( figure 3(a)). The significant positive nitrate anomalies found in figure 2(a) started about 42 months before the increase in Agulhas chlorophyll. The westward moving nitrate anomalies take about 2-3 years to reach the western boundary and they experience disruptions and enhancements during the propagation period, presumably due to the interplay between the Rossby wave signal and other sources of nitrate variability.
The regressed fields of physical and biogeochemical variables onto the satellite chlorophyll in the Agulhas system show that the initial positive nitrate signals are triggered and amplified by two consecutive La Niña-like patterns at around lag 42 and 32 months (figures 3(b) and (c)), which is consistent with the result shown in figure 3(a). At the time lag of 42 months, a weak La Niña signal in the equatorial pacific is captured with cyclonic wind stress anomalies over the southeastern Indian Ocean. This low-level wind perturbation is consistent with a Gill-type response to diabatic warming caused by the eastern Indian Ocean precipitation anomalies (Gill 1980). This is further supported by an idealized model experiment using a linear baroclinic model (Watanabe and Kimoto 2000) forced by a prescribed diabatic heating over the eastern Indian Ocean at around 110 • E, 15 • S (supplementary figure 8). The cyclonic winds can, in turn, lead to the increase in upper ocean nitrate through the upwelling of nutrient-rich waters. A similar feature is also found at the time lag of 32 months ( figure 3(b)). That is, a strong La Nina signal increases the precipitation in the eastern Indian Ocean, and then it generates cyclonic flow in the southern subtropical Indian Ocean, leading to the increase in upwelling and nitrate. The nitrate anomalies propagate westward with time and finally reach the western boundary at 0 month lag (figure 3(f)). Such westward nitrate propagation pattern is commonly detected after major La Niña events despite minor discrepancies in its magnitude and structure, and also detected after El Niño events with a reversed sign of nitrate anomalies (supplementary figures 9 and 10). Overall inter-basin interacting processes associated with Agulhas chlorophyll prediction skill are also found in the observation-based datasets (supplementary figure 11).
Previous work has shown that Indian Ocean dipole (IOD) events are accompanied by a westward-propagating Rossby wave in the southern tropical Indian Ocean generated by wind stress forcing factors (Masumoto andMeyers 1998, Vinayachandran et al 2002), implying the potential role of IOD in chlorophyll prediction skill. The ESM-based prediction system used here shows significant IOD prediction skill within 1.5 year prediction horizon, however, chlorophyll variability in the Agulhas system generally shows a less significant correlation with IOD index compared to ENSO index (supplementary figure 12). The strong correlation between ENSO and Agulhas chlorophyll at long time lags indicates a primary role of ENSO-Indian Ocean teleconnection in initiating slowly moving subtropical Rossby wave and subsequent chlorophyll variability in the western boundary. In the Somali system, however, chlorophyll is significantly correlated with IOD particularly for short time lags, with a significant correspondence with ENSO for longer time lags. This implies that Somali chlorophyll is affected by both ENSO and IOD and/or their co-variability. Lastly, we examine whether the westward moving Rossby wave signal linked to Agulhas chlorophyll anomalies is identified in the observation. The observed ocean temperature analyzed here is from the EN4 objective analyses data (Good et al 2013). The ocean temperature averaged in the 100-200 m depth, i.e. oceanic heat content, is used to diagnose the upwelling Rossby wave that transports nitrate signal to the western Indian Ocean. In both model and observation, the heat content anomalies regressed onto satellite chlorophyll in the Agulhas region shows the significant negative anomalies initiated four years before the increase in Agulhas chlorophyll, and this negative signal moves westward afterward (figure 4). The westward-moving heat content anomalies linked to Agulhas chlorophyll are also present in the predicted ocean temperature (supplementary figure 13). The overall results here imply that despite some discrepancies between model and observations, the ESM used here captures the basin-wide physical processes and its accompanying biogeochemical responses, which is a key source of chlorophyll prediction skill in the western Indian Ocean LMEs.

Discussion and conclusion
The ESM-based marine biogeochemical prediction system used in the present study shows significant skill in predicting satellite-derived chlorophyll fluctuations in the east African coastal regions. The skillful chlorophyll prediction arises primarily from successfully simulating coupling processes between atmosphere, ocean, and marine biogeochemistry in the Indian Ocean. The inter-basin interacting process between atmosphere and ocean exemplifies the compelling predictive potential of global ESM even for regional biogeochemical prediction, while regional model-based prediction systems may be limited in predicting such inter-basin process. This result is similar to enhanced prediction of climate variability by incorporating inter-basin precursors (Ham et al 2013, Cai et al 2019.
The ESM-based biogeochemical prediction system can extend beyond chlorophyll such as oxygen, net primary production, and zooplanktonSuccessful prediction of these key biogeochemical variables may provide richer information on the linkage between bottom-up drivers and fisheries production and its mechanistic principles across globally distributed ecosystems,thus offering considerable potential for anticipatory dynamic management of marine resources (Tommasi et al 2017). Realizing the full potential of global earth system predictions for marine ecosystem and resource applications, however, will require a range of modeling and observational advances. The lack of global-scale biogeochemical observations often hinders the generation of optimal biogeochemical initialization data and the validation assessment of model products. In this regard, assimilating the remotely sensed ocean color and the recently available biogeochemical Argo floats data into the ESM can bring additional gains in biogeochemical prediction skill and may eventually provide a robust time series for reforecast verification. Recent studies support this point by showing the improvement of biogeochemical reanalysis from biogeochemical assimilation and the degradation of biogeochemical predictions from the perturbations in biogeochemical initial conditions (Ciavatta et al 2014, Ford and Barciela 2017, Salon et al 2019, Ford 2021. Future inclusion of dynamic river nutrient fluxes into the ESM may also improve the coastal biogeochemical prediction skill given that coastal biogeochemical cycles are substantially controlled by river nutrients (Walker andRabalais 2006, Sigleo andFrick 2007). For example, a recent observational study showed that chlorophyll concentration in the east African coast is largely controlled by river discharge during the rainy season (Mutia et al 2021). In addition to the dynamic coupling between river and coastal water biogeochemistry, much higher resolution than available in the present model (1 • , ∼100 km) would be preferable to adequately represent the complex coastal physics and biogeochemistry towards improvement of the coastal chlorophyll prediction. Such high resolution models are expected to better simulate coastal retention and residence time of oceanic tracers probably due to better representation of coastal bathymetry and complex fine-scale dynamics (Liu et al 2019) and should be a priority for future work.

Data availability statement
The data that support the findings of this study are available upon reasonable request from the authors.