Key climate change stressors of marine ecosystems along the path of the East African coastal current

the 21st century, to identify key regionally important climate change stressors over the East African Coastal Current (EACC) that flows along the coasts of Kenya and Tanzania. We also discuss these stressors in the context of projections from lower resolution CMIP5 models. Our results indicate that the main drivers of dynamics and the associated ecosystem response in the TWIO are different between the two monsoon seasons . Our high resolution model projects weakening of the Northeast monsoon (December – February) winds and slight strength- ening of the Southeast monsoon (May – September) winds throughout the course of the 21st century, consistent with CMIP5 models. The projected shallower mixed layers and weaker upwelling during the Northeast Monsoon considerably reduce the availability of surface nutrients and primary production. Meanwhile, primary production during the Southeast monsoon is projected to be relatively stable until the end of the century. In parallel, a widespread warming of up to 5 ◦ C is projected year-round with extreme events such as marine heatwaves becoming more intense and prolonged, with the first year-long event projected to occur as early as the 2030s. This extreme warming will have significant consequences for both marine ecosystems and the coastal populations dependent on these marine resources. These region-specific stressors highlight the importance of dy- namic ocean features such as the upwelling systems associated with key ocean currents. This indicates the need to develop and implement a regional system that monitors the anomalous behaviour of such regionally important features. Additionally, this study draws attention to the importance of investment in decadal prediction methods, including high resolution modelling, that can provide information at time and space scales that are more directly relevant to regional management and policy making.


Introduction
Living marine resources are expected to respond to anthropogenic climate stressors in multiple and complex ways, often compounded by non-climatic pressures such as marine pollution and habitat degradation. Many aspects of this response are uncertain, with the majority of impacts neither spatially nor temporally homogeneous over the global ocean (e.g. Mizuta et al., 2014;Popova et al., 2016). The first responses of marine ecosystems to stressors have already been observed, with warming Sea Surface Temperatures (SSTs) leading to widespread coral bleaching (e.g. Ateweberhan et al., 2011), and changes in major currents and the intensification of extreme events such as marine heatwaves causing shifts in the distribution of marine species (e.g. Coleman et al., 2013;Cetina-Heredia et al., 2015).
The tropical Western Indian Ocean (TWIO), encompassing the region 15 o S-10 o N, 35-55 o E, has largely escaped the attention of studies investigating future climate change projections for the marine environment. The main focus of Indian Ocean regional studies is instead on projected changes around Australia (e.g. Jourdain et al., 2013;Bull et al., 2020), the Indian subcontinent (e.g. Parvathi et al., 2017;Mohan and Bhaskaran, 2019a) or further south, at the boundary with the Southern Ocean (e.g. Nevison et al., 2016;Schneider and Reusch, 2016). However, a third of the regional population lives within 100 km of the East African coastline  where marine living resources are vital for food security (Bell et al., 2016;Taylor et al., 2019). Food insecurity is likely to grow in this region (McClanahan et al., 2013) given the potential loss of suitable habitat conditions for species targeted by wild-capture fisheries, and their consequential poleward distributional range retreat under climate change (Pinksy et al., 2018). Unlike in other regions, these shifts are not accompanied by the arrival of a sufficient number of new species from lower latitudes, also experiencing poleward distributional shifts, leading to an overall loss of species from the regional pool (Oremus et al., 2020). These changes could exacerbate existing fisheries conflicts in East Africafrom 1990 to 2017 in Tanzania, 80% of disputes between fishers and the government were due to declining fish populations and ecosystem changes (Glaser et al., 2018). Consequently, understanding how the marine ecosystem in coastal East Africa will respond to climate change is vital if future management decisions are to ensure food security for the regional population.
The ocean dynamics that shape the marine ecosystems and fisheries of the TWIO are subject to considerable variability associated with the monsoonal winds (Jury et al., 2010;Jebri et al., 2020). During the Southeast monsoon (SEM) that prevails from May to September, strong southeasterly winds cause the major surface currents to reach their peak speeds (Schott and McCreary, 2001). These include the South Equatorial Current (SEC), the Northeast Madagascar Current (NEMC), the East African Coastal Current (EACC) and the Somali Current (SC) -see Jacobs et al. (2020a,b) for a schematic. During this season, the fast, northward flow of the EACC and the NEMC result in a dynamic uplift upwelling along the Tanzanian and Kenyan coastlines and a westward advection of nutrient-rich waters from the northern tip of Madagascar towards the East African coast, contributing to elevated primary production and increased fish catches (Jebri et al., 2020). In contrast, the northeasterly winds that prevail during the Northeast monsoon (NEM) induce Ekman-driven upwelling along the Tanzanian and Kenyan coasts (Varela et al., 2015;Jebri et al., 2020). Additionally, the EACC is weaker during this season while the SC reverses to flow southwards (Schott and McCreary, 2001). The convergence of these currents forms the Somali-Zanzibar Confluence Zone (SZCZ; see Fig. 1) near the North Kenya Banks (NKBs). The subsequent divergent flow away from the coast leads to the occurrence of a seasonal shelf-edge upwelling (Jacobs et al., 2020a,b), which likely sustains a productive, and currently expanding, offshore fishery (Maina and Osuka, 2014).
Ocean currents are critically important for fisheries in the TWIO (Jury et al., 2010), with distinct current-induced upwelling systems prominent during each monsoon season (Jacobs et al., 2020a,b;Jebri et al., 2020). Changes in ocean circulation under climate change are one of the key stressors for marine ecosystems with some of the fastest warming occurring in regions with shifting boundary currents (Wu et al., 2012;Popova et al., 2016). In turn, fast warming creates a two fold physiological challenge for regional species populations: whilst warming is a key life-history driver (Portner and Knust, 2007), and thus a major source of stress for species experiencing climate change, warming modifies the distribution and availability of oxygen and nutrients across the water column, further compressing habitats (Breitburg et al., 2018). In the Western Indian Ocean, the NEMC is projected to weaken by 14% by the end of the century across Coupled Model Intercomparison Project Phase 5 (CMIP5) models (Stellema et al., 2019), while Popova et al. (2016) find evidence for its potential offshore shift by mid-century. As the NEMC, which feeds the EACC, has been found to be important for advecting nutrients (Jebri et al., 2020), this could have implications for productivity, and hence fisheries, along the East African coast.
The projected changes in winds in the TWIO will likely also affect primary production via altered intensity of upper ocean mixing. For example, in the Arabian Sea the winds are projected to weaken by 6.5% by the end of the century during the NEM, which reduces deepening of the mixed layer and inhibits the seasonal chlorophyll bloom (Parvathi et al., 2017). Primary production in the tropics is expected to decline by up to 30% Kwiatkowski et al., 2017) with the northwest Indian Ocean experiencing the greatest projected declines, most notably in the Somali upwelling region (Steinacher et al., 2010;Cabre et al., 2015;Roxy et al., 2015;Nakamura and Oka, 2019). These declines are primarily due to increased stratification driven by rising SST which reduces vertical mixing and the supply of nutrients to the euphotic zone (Bopp et al., 2001;Steinacher et al., 2010;Fu et al., 2016), although modification of monsoonal winds and ocean currents may also change nutrient supply via upwelling. However, the typically low resolutions (~1 • ) of the ocean components of CMIP5 modelsand even CMIP6 modelsdo not yet adequately resolve changes in ocean currents and upwelling systems, making future projection of such changes problematic.
The continuous anthropogenic warming of the TWIO will further increase the risk of marine heatwaves (e.g. Hobday et al., 2016), with wide-ranging consequences for marine ecosystems, increasing the compounding environmental pressure on species targeted by fisheries and coral reef communities with high conservation and economic value (Jury et al., 2010;Bell et al., 2016). Furthermore, changes in currents and reduced productivity may lead to the alteration of marine food webs (Behrenfeld et al., 2006), causing changes in the geographical distributions of various species (Brander, 2010). As this region is particularly reliant on living marine resources for its food security (Taylor et al., 2019), it is essential to understand and predict the future changes to marine ecosystems so that appropriate adaptive measures can be put in place.
In this paper we focus on the TWIO between 15 o S and 10 o N to encompass the key dynamical features that shape the EACC regionthe key topic of this Special Issue. Using a high-resolution (1/4 • ) ocean model run under the high emissions Representative Concentration Pathway 8.5 (RCP8.5) scenario to the end of the 21st century, we synthesise the key climate change impacts on the marine ecosystems of the wider EACC region and put them into the context of CMIP5 projections. The model employed in this study is of higher resolution than currently available earth system models (Popova et al., 2016) and thus allows a more detailed investigation of regional features and their impact on marine ecosystems. The greater regional detail represented by this model thus allows us to propose key climate change indicators which can be monitored in order to identify and document early signs of ecosystem responses to the continuing anthropogenic climate change in this region.
The paper is organised as follows: section 2 presents the methods and model used in this study, including a limited validation with observational data of both this model and an ensemble of CMIP5 models; section 3 presents results from the future projection of the high-resolution global ocean model used here, including a Lagrangian analysis of circulation change; and section 4 presents discussion and conclusions, including the potential implications of the results on fisheries management.

High resolution ocean projection
The physical framework used in this study is the Nucleus for European Modelling of the Ocean (NEMO) model (Madec, 2015). This is composed of an ocean general circulation model coupled to a sea-ice model, LIM2 (Timmermann et al., 2005). NEMO version 3.5 is used, with a horizontal spatial resolution of 1/4 • , which corresponds to grid cells of 27.8 km at the equator and 6 km at the poles. There are 75 vertical levels, which range in thickness from 1 m at the surface to 200 m in the deep ocean. Vertical mixing is parameterized using a turbulent kinetic energy scheme (Gaspar et al., 1990), which includes modifications from Madec (2008). Biogeochemistry is represented by the plankton ecosystem model MEDUSA-2 (Yool et al., 2013a,b). This is a size-based, intermediate complexity model that divides the plankton community into 'large' and 'small' components and resolves the elemental cycles of nitrogen, iron, carbon, silicon and oxygen.
The model is forced at the surface using output from a simulation of the Earth system model HadGEM2-ES (Jones et al., 2011), which includes representations of the terrestrial and oceanic carbon cycles, atmospheric chemistry and aerosols (Collins et al., 2011). This simulation was performed as part of the CMIP5 activity that contributed towards the Intergovernmental Panel on Climate Change (IPCC) Assessment Report 5 (AR5). To reduce the computational costs associated with this high-resolution model, its initial conditions are derived from a 'twin' 1 • model, which was spun up from 1860 to 1975 under the same forcing. The 1/4 • model is then run under historical atmospheric pCO 2 concentrations from 1975 to 2005 and then under the IPCC RCP 8.5 pathway from 2006 to 2099. This pathway is a scenario in which CO 2 emissions rise continuously over the 21st century, causing an additional 8.5 W m -2 of additional radiative forcing by 2100. It should be noted here that the overlap period  between the 1 • and 1/4 • simulations agree well across a broad range of physical and biogeochemical metrics (Yool et al., 2015). Further details of model implementation, forcing, equilibration and verification can be found in Yool et al. (2013a;2013b;2015) and Popova et al. (2010;. The projection is hereafter referred to as NEMO-MEDUSA. The spatial resolution of NEMO-MEDUSA is higher than that used in the CMIP5 models and shows good skill in reproducing the main dynamical and biogeochemical features; notably the ocean circulation, western boundary current systems and upper biogeochemistry and ecosystems, which makes it suitable for regional analysis (Popova et al., 2016).

Model validation
The skill of the 'present-day' state of the model is evaluated in the TWIO between 15 o S and 10 o N using in situ and satellite-derived observations for three properties directly relevant to this study: surface circulation, surface nutrient status and biological productivity.
Surface circulation is compared with absolute geostrophic velocities from the satellite altimetry product AVISO (distributed by the Copernicus Marine Environment Monitoring Service 1 ). The specific data used here is the delayed-time 'Update' DUACS-DT2018 version, with geostrophic zonal and meridional velocities computed weekly on a 1/4 • x 1/4 • grid for the period 2000 to 2009. Fig. 1 compares the modelled surface currents with observations as a decadal mean from 2000 to 09 for the two East African monsoon seasons. Decadal means capture the persistent features in the circulation whilst avoiding short-term variability dominated by mesoscale eddies. It shows good agreement between the observed and modelled surface circulation during both seasons with the position and monsoonal variability of the major currents captured well by the model. Specifically, the fast, northward flow of the EACC and SEC during the SEM and the widespread weakened flow and reversal of the SC to form the SZCZ 1 http://marine.copernicus.eu/services-portfolio/access-toproducts/?optio n=com_csw&view=details&product_id=SEALEVEL_GLO_PHY_L4_REP_OBSE RVATIONS_008_047.
during the NEM.
Surface dissolved inorganic nitrogen (DIN) is compared with the observational climatology  from the World Ocean Atlas 2013 (WOA13; Garcia et al., 2014). WOA13 is an objectively analysed product that includes all available observations from ship-based measurements and profiling floats. The model is able to simulate the spatial pattern of the climatological mean during each monsoon season (Supplementary Fig. 1) but generally underestimates the DIN concentrations in the TWIO. However, elevated concentrations near the East African coast in WOA13 are likely an artefact of extrapolation as the climatologies are based on a very small number of observations in this region.
Following Yool et al. (2013a, b) and Popova et al. (2016), observational primary production is estimated here using a simple average calculated from three satellite-derived estimates, the Vertically Generalized Production Model (VGPM; Behrenfeld and Falkowski, 1997), the Eppley-VGPM (Carr et al., 2006) and the Carbon-based Production Model (CbPM; Westberry et al., 2008) models. NEMO-MEDUSA accurately simulates the observed spatial pattern of depth-integrated primary production (IntPP) in the TWIO during both monsoon seasons ( Supplementary Fig. 2). Enhanced productivity occurs north of the equator during the SEM, while overall weaker IntPP exists during the NEM, except offshore of Kenya.
The next section assesses how well NEMO-MEDUSA represents the surface DIN concentration and primary production in the TWIO alongside the CMIP5 models. See Popova et al. (2016) for an extensive model validation.

IPCC AR5 model comparison
In order to assess where the NEMO-MEDUSA future projection belongs in the multi-model ensemble for the Indian Ocean, we assess the biogeochemistry of a suite of CMIP5 models (Taylor et al., 2012). We use CMIP5 models as opposed to those from CMIP6 to match the RCP8.5 scenario used with NEMO-MEDUSA. CMIP6 models use socioeconomic pathways as opposed to RCPs (Simpkins, 2017), and are therefore not directly comparable to NEMO-MEDUSA. First, we apply skill metrics to compare the models to observational data in order to identify an "inner ensemble" of models that simulate the observations best. This inner ensemble is then used to frame the assessment of future projections from 2006 to 2100 under the Representative Concentration Pathway 8.5 experiment (RCP8.5; i.e. a doubling of CO 2 emissions, and radiative forcing of 8.5 W m -2 by 2100). A list of the models used, in addition to NEMO-MEDUSA, is shown in Table 1. Our selection of CMIP5 models is based on those that include the variables DIN and IntPP.
We employ statistical metrics to score the models on how well they reproduce the spatial pattern of DIN over the tropical Indian Ocean. DIN is chosen because of good observational data coverage (here from the World Ocean Atlas, 2013; WOA13; Garcia et al., 2014) and its role as a key macronutrient driving marine primary productivity makes it a good indicator of the quality of a model's biogeochemistry. For consistency, all modelling projections analysed are linearly interpolated onto a standard 1 • x 1 • grid, and averaged over the first decade of the future projection run from 2006 to 2015. Over the wider Indian Ocean, most models represent the spatial pattern of observational data well, with greater concentrations in the Southern Ocean and northwest Indian Ocean ( Supplementary Fig. 3). However, DIN is generally underestimated in the tropics, which is a common feature across ESMs due to weaker upwelling in coarse resolution models (e.g. Seferian et al., 2013).
However, some models do not capture the spatial pattern of observational data well and/or overestimate DIN concentrations across the entire region ( Supplementary Fig. 3). To quantify these differences, skill metrics are calculated between the decadal mean of each model and WOA13 over the tropical Indian Ocean (20 o S-20 o N, 30-130 o E). It is important to assess a region larger than the TWIO, due to the coarse resolution of some models and importance of upstream processes. At the same time, we have chosen to restrict the area to 20 o S to avoid the large increase of nutrients towards the Southern Ocean, which may skew the metrics towards models correctly reproducing this gradient, rather than relatively low background values in the tropics. The metrics used to assess model skill are Pearson's correlation coefficient, r: √ where x and y are the observed and modelled data points respectively while x and y are the respective means; and the root mean square error (RMSE): where n is the number of grid points. These metrics estimate how much each model projection varies with the same tendency as the observational dataset, and the magnitude of the discrepancies between them, respectively (Stow et al., 2009). The analysis here complements prior studies to provide model assessment for the tropical Indian Ocean (as in Allen et al., 2007;Stow et al., 2009;Bopp et al., 2013;Ilyina et al., 2013;Seferian et al., 2013;Cabre et al., 2015;Rickard and Behrens, 2016;Bao and Li, 2016;Mohan and Bhaskaran, 2019b). The r and RMSE for each model are listed in Table 1. A large spread exists across both metrics with NEMO-MEDUSA achieving the highest r, 0.61, and CNRM-CM5 the lowest, 0.01. CNRM-CM5 also has the largest RMSE, 6.5 × 10 − 3 while MPI-ESM-MR and IPSL-CM5B-LR have the lowest RMSE, 0.45 × 10 − 3 . Models that fall outside the range of 2 standard deviations for each metric are considered the "outer ensemble" and include CanESM2, CMCC-CESM, CNRM-CM5, GFDL-ESM2G, GFDL-ESM2M and NorESM1-MR. The remaining models are characterised as the "inner ensemble" and are used to analyse the future projections of DIN and IntPP. The projected changes in DIN over the WIO (averaged from 15 o S-5 o N, 35-55 o E to avoid the Somali upwelling region) are shown in Fig. 2. Despite the models all having a negative bias, the inner ensemble best represent the spatial pattern of DIN across the tropical Indian Ocean ( Supplementary Fig. 3). The 3 IPSL models project increases of 84-158% by the end of the century, which is in stark contrast to the other 10 models, which project a decline ranging from 16 to 69%.
IntPP is shown as a decadal average over the first decade (2006-2015) of the future projection runs and from the observational Table 1 List of CMIP5 models (and NEMO-MEDUSA) and their corresponding skill metric scores based on the decadal mean period 2006-2015 (all months) when interpolated onto a 1 • grid; the correlation coefficient, r, and the root mean square error, RMSE, when compared to climatological DIN from WOA13. Models in italics are not in the inner ensemble and are excluded from further analysis.  (Seferian et al., 2013).
Observationally-estimated IntPP uses satellite-derived observations of chlorophyll and temperature together with productivity models rather than direct measurements, which adds an extra uncertainty when calculating skill metrics (Anav et al., 2013). Thus, we use the same "inner ensemble" based on the skill metrics calculated for DIN. The projected changes in IntPP over the WIO, are shown in Fig. 3. All inner ensemble models project a decline in productivity in this region ranging 19-40% (mean − 26%) by the end of the century. This decline parallels the projected global reduction of IntPP (Steinacher et al., 2010;Bopp et al., 2013;Cabre et al., 2015). Although there is evidence for greater projected declines of up to 67% for global equatorial regions, there is also a considerable model spread (Cabre et al., 2015).
This analysis reveals that future projections of annual mean IntPP and DIN in NEMO-MEDUSA are consistent with inner ensemble of the CMIP5 models in the WIO region. However, due to the considerable seasonal variability associated with the monsoonal winds, we further analyse NEMO-MEDUSA's future projection in the model on seasonal timescales to investigate the contrasting impacts in each monsoon season.

Lagrangian experiments
In order to assess projected changing upstream connectivity to the East African coastline, Lagrangian particle-tracking experiments are performed. These experiments use the Ariane (Blanke and Raynaud, 1997) particle-tracking software, forced using circulation output from NEMO-MEDUSA. Starting with a defined set of initial particle positions, Ariane uses an analytic method based on a bilinear interpolation of model velocity fields to calculate how each particle will move in time (Blanke and Raynaud, 1997). The velocity fields used as input by Ariane are pre-calculated 5-day mean fields from NEMO-MEDUSA. Particle transport can be calculated either forwards in time to investigate the downstream connectivity from defined source locations or backwards in time to investigate the upstream connectivity to defined sink locations  (as here).
Ariane has been used extensively (e.g. Van Gennip et al., 2017;Robinson et al., 2017;Jacobs et al., 2019;Kelly et al., 2020) in conjunction with NEMO, as it provides an efficient way of analysing advective pathways for novel experiments without requiring additional computationally-expensive runs of the full model itself. For a more in-depth discussion of Lagrangian analysis in general, the reader is referred to Van Sebille et al. (2018), and for Ariane specifically to Blanke and Raynaud (1997).
In these experiments, the backtracking mode of Ariane is used, analogously to the experiments in Popova et al. (2019). Here, particles are released annually at the beginning of each June over decadal periods from the start (2000-09) and end (2090-99) of the 21st century. June is selected to investigate projected changes in the surface circulation during the SEM (section 3.2). Each release is comprised of a cohort of 2183 particles initialised in a uniform distribution at the sea surface over the southern part of Tanzania's coastal zone (south of 8 o S, west of 40.5 o E, and north of the Mozambique border), which is where the NEMC reaches the coast and feeds the EACC. These particles are then backtracked upstream to their sources 30 days prior to identify interannual variability in sources, and how spatial patterns of these change under the emissions scenario examined.

Changes in East African monsoon winds
The mean wind field during the SEM and NEM, averaged from 2000 to 2009, is shown in Fig. 4a and b respectively. South of the equator, the prevailing south-easterly winds during the SEM are clearly visible, peaking at 9 m s − 1 north of the tip of Madagascar. North of the equator, the prevailing wind is south-westerly, blowing adjacent to the Somali coast, with speeds exceeding 10 m s − 1 . During the NEM, the wind reverses to become north-easterly, exceeding 10 m s − 1 off the coast of Somalia and 8 m s − 1 off the coast of Tanzania with comparatively lower speeds of 5 m s − 1 to the north of Madagascar. This is largely in agreement with CMIP5 models (Parvathi et al., 2017). However, despite CMIP5 models showing a bias of up to 1 m s − 1 during the NEM over the equatorial WIO, the biases are larger in the Arabian Sea and towards the Southern Ocean (Parvathi et al., 2017;Mohan and Bhaskaran, 2019a).
Projected changes in the wind speed by the end of the century are shown in Fig. 4c and d. First to note is the stark contrast in the sign of the anomaly between the monsoon seasons. During the SEM, there is a projected strengthening of the wind speed by up to 1 m s − 1 north of Madagascar, as well as in the equatorial region and adjacent to the Somali coast. These anomalies are generally found in the areas of the greatest wind speeds seen in Fig. 4a. Negative anomalies reaching − 0.6 m s − 1 are found in regions of slower mean winds, e.g. in the Mozambique Channel. Contrastingly, during the NEM, there is a stronger projected weakening across the entire region of up to 1.6 m s − 1 . The greatest projected reduction occurs adjacent to the Tanzanian coast and at the equator east of 50 o E.
Trends in the spatially averaged wind speed inside the EACC (see Fig. 4b), that encompasses the offshore waters of Tanzania and Kenya, during the SEM and NEM are shown in Fig. 5a and b. The contrasting trends between the monsoon seasons are further highlighted here, with an increase of 6% in mean wind speed over the whole region by the 2090s during the SEM, and a corresponding decline of 19% during the NEM, both consistent with CMIP5 models (McInnes et al., 2011;Parvathi et al., 2017). When analysing the future trends of climatic stressors, an important question to ask is when a stressor begins to fall outside of the range of its natural (before anthropogenic influence) or baseline (a fixed period) variability and if this constitutes a consistent trend. Following Popova et al. (2016), we consider 2000-19 as a period describing "baseline variability". Here, and in further analysis, 2 standard deviations are taken to be a representative range of this variability. The projected increase during the SEM occurs until the mid-century before it stabilises until the end of the century, with only a handful of years just exceeding its baseline range of variability (indicated by the thick black lines). In contrast, during the NEM, the wind speed undergoes a steady decline until 2100. Despite the larger range of baseline variability, during this season the wind is projected to regularly dip below this from the 2050s. The majority of the projected changes in overall wind speed are predominantly due to the meridional component with stronger northward winds projected during the SEM and weaker southward winds projected during the NEM ( Supplementary Fig. 5).
The projected trends in wind speed (from 2000 to 2030) are generally consistent with the past trends from 1980 to 2010 ( Supplementary  Fig. 6). This overlap indicates that the effects of climate change are already manifesting in this region. The strengthening of wind speed since 1980 during the SEM in the EACC region is projected to continue at around 0.01 m s − 1 y − 1 . During the NEM, the downward trend in wind speed since 1980 is projected to decline at a slightly faster rate of up to − 0.02 m s − 1 y − 1 in this region. Fig. 5c shows changes in the seasonal cycle between the first and last decade of the century. The reduced wind speeds are evident during the NEM with the greatest decline occurring in December (34%). There is also evidence for a slight strengthening (+19% in October) and lengthening of the strong SEM winds by the end of the century, which is in line with the latest IPCC findings (McInnes et al., 2011;IPCC, 2013). The increased duration of the SEM is due to an extension at the end of the season, with no changes observed at the beginning of the season, i.e. from May-July. Additionally, there are no projected changes in wind speed during the two inter-monsoon seasons, which suggests there will be minimal changes in the timing of the rainfall seasons on the East African coast. Hence, we will focus on the SEM and NEM seasons.

Changes to surface circulation
Changes in the wind field are likely to have consequences for the surface ocean circulation, one of the key climatic stressors of marine ecosystems (Van Gennip et al., 2017). The mean surface currents, averaged from 2000 to 09, during each monsoon season are shown in Fig. 6a and b. Fig. 6c shows the projected changes in the surface currents by the end of the century during the SEM. There is a projected strengthening of the EACC (north of 5 o S) and SC of up to 0.1 m s − 1 by the 2090s, consistent with a projected increase in wind speeds. South of 5 o S, the anomaly pattern of reduced flow adjacent to the coast and strengthened flow further offshore is indicative of an offshore shift of the EACC. This could cause a reduction in the dynamic uplift upwelling that occurs on the western side of the EACC along the Tanzanian coast (Jebri et al., 2020). This, combined with reduced flow through the Zanzibar and Pemba Channels, may lead to a reduction in primary production due to reduced availability of nutrients near the coast, potentially causing negative consequences for species targeted regionally by wild capture fisheries. There is also evidence for a weakened Southern Gyre (from 6 o S-1 o N,40-45 o E) and changes in the offshore veering of the SC (north of 5 o N) during the SEM. However, the most notable change in this region is the widespread negative anomaly (exceeding − 0.1 m s − 1 ) in the NEMC flow during the SEM, with evidence of a weakened SEC along 15 o S. The weakened NEMC could cause less nutrient-rich waters to be advected along Tanzanian and Kenyan coasts from the northern tip of Madagascar (Jebri et al., 2020), and hence a potential decrease in fish catches.
During the NEM (Fig. 6d), there is a slight strengthening of the northward-flowing EACC, of up to 0.05 m s − 1 , which meets a weaker, exceeding − 0.1 m s − 1 , southward-flowing SC. The projected strengthening of the northward flowing EACC and weakening of the southward flowing SC may cause the location of the SZCZ to reside in a more northward position towards the end of the century. As the position of the SZCZ is a major control on the upwelling and primary production over the NKBs (Jacobs et al., 2020a,b), the northward shift of the SZCZ could have a detrimental effect on the catch potential of the local offshore fishery during this time of year (Maina and Osuka, 2014). Further offshore in the tropical Indian Ocean (from 5 o S-5 o N), changes are found in the equatorial currents with evidence for a weaker eastward-flowing SECC and potential current shifts north of the equator.
The most notable change in this region is the weakened NEMC during the SEM (Fig. 6e). The weakening is projected to occur across all months but is greatest in May, with a change of − 26% by the 2090s (Fig. 6c). The reduced velocity is not as great during the NEM but there is also evidence that the NEMC has shifted offshore (Fig. 6d). Fig. 6f shows the time series of surface velocity across 47 o E during the SEM (green line in Fig. 6c), which reveals a stable westward flow until the 2040s whereafter it begins to frequently attain velocities outside of its baseline range of variability (indicated by the black lines; representing two standard deviations from 2000 to 2019). By the 2090s, the average NEMC velocity in this region is projected to weaken by − 16.5% during the SEM. This could have a negative impact on the productivity of the waters along the Tanzanian coast during the SEM as the NEMC has been found to be at the core of the advective supply of nutrients to the region (Jebri et al., 2020).

Changes in ocean connectivity
One of the important consequences of the changes in ocean circulation is the resulting change in connectivity between ocean regions. Ocean circulation is one of the underlying mechanisms of ecological connectivity between remote ecosystems (e.g. Popova et al., 2019) and in some regions it also acts as a mechanism of nutrient supply (e.g. Jebri et al., 2020). As these are likely to impact the productivity that underpins the occurrence and abundance of fisheries species (or species these depend upon) potential future changes are important to understand. To investigate potential connectivity change during the SEM, when considerable circulation changes are projected (Fig. 6), Lagrangian backtracking experiments are performed using the Ariane particle tracking tool (see Methods) for the beginning (2000s) and end (2090s) of the century. Virtual particles are initiated at the Tanzanian coast and backtracked for 30 days to find their source, i.e. where water arriving at the Tanzanian coast originates from. This experimental design is analogous to that of Jebri et al. (2020). The mean particle densities after 30 days of backtracking released each June from 2000 to 2009 and 2090-99 are shown in Fig. 7a and b respectively. June was selected on the basis of the substantial changes projected in the NEMC from May-June (Fig. 6e).
Numerous pathways are visible, including the NEMC that flows around the northern tip of Madagascar. The most obvious difference is the absence of particles off the NE coast of Madagascar in the 2090s. During the 2090s, the NEMC pathway originates at the tip of northern Madagascar with an eastern limit of 51 o E. In contrast, during the 2000s, this pathway originates further south on the east coast of Madagascar. This corroborates the Eulerian analysis (Fig. 6) that showed the NEMC is projected to be slower by the end of the century, reducing the connectivity time between this region and the Tanzanian coast. Advection of nutrients by the NEMC was found to be one of the main controls on interannual chlorophyll-a variability along the Tanzanian coast during the SEM (Jebri et al., 2020). The weakening of this current and its associated connectivity to the Tanzanian coast would imply a reduction in the delivery of nutrients to this region and hence, the primary production during this season. Further, the SEC and the NEMC are predicted to provide an important source of pelagic larvae, such as fish and corals, to the East African mainland from modelling pelagic larval duration (Crochelet et al., 2016), which is reflected in higher diversity and abundance of coral reef fauna in northern Mozambique (Obura 2012;Samoilys et al., 2019). Thus, weakening of these currents is likely to have a negative impact on the recruitment processes of coral reef fauna in mainland east Africa.

Thermodynamical and biogeochemical changes
In addition to ocean circulation changes, there are also projected changes in the thermodynamics and biogeochemistry of the TWIO. Fig. 8a and b shows that a widespread warming of up to 5 • C is projected by the end of the century under the RCP8.5 scenario during both monsoon seasons (decadal mean for 2000s is shown in Supplementary  Figs. 7a and 7b). This projected warming is consistent throughout the year with a mean SST increase of 4 • C in the EACC area during every month by the 2090s (Supplementary Fig. 8a). The time series showing the projected increase of SST during both monsoon seasons over the same area ( Fig. 9a and b) reveals a steady increase until the end of the century (>4 • C) with values exceeding the baseline range of variability regularly from the 2020s.
The projected changes in the mixed layer depth (MLD) are consistent with projected changes in the winds (Fig. 5). During the SEM, there is a mixed layer deepening of up to 10 m (on climatological values of up to 55 m; see Supplementary Fig. 7c) in a broad region that extends from the northern tip of Madagascar and up the coast of East Africa, with the greatest anomalies seen south of the equator (Fig. 8c). This pattern is consistent with the projected strengthening in the SEM wind speeds. By contrast, the projected widespread weakening of the winds during the NEM has led to the subsequent shoaling of the mixed layer by up to 20 m (Fig. 8d). Warmer surface waters may have also led to greater stratification, which may have also contributed to the shoaling of the mixed layer. The area of greatest shoaling occurs along the equator and to the south of this, west of 45 o E adjacent to the Tanzanian and Kenyan coastlines, where the MLD ranges from 40 to 70 m ( Supplementary  Fig. 7d). Supplementary Fig. 8b reveals that in the EACC area, the projected changes in the seasonal cycle of MLD are consistent with that of the wind speed in Fig. 5c. Specifically, the greatest change is during the NEM (− 39% in December), with minimal change in the intermonsoon periods and a slight strengthening and lengthening of the SEM season (+27% in September). The time series in Fig. 9c and d further highlight the differences between monsoon seasons with a weaker, but still positive trend during the SEM (+14% by the 2090s). This is accompanied by relatively low interannual variability, coherent with the regional wind field (Fig. 5b), which means the peaks of MLD do not exceed the baseline range of variability until the middle of the century. In contrast, there is a clear, negative trend (− 31% by the 2090s) during the NEM with much greater interannual variability evident, which leads to the 2050s being the first decade where some years exceed the baseline range of variability.
Projected changes in surface DIN and IntPP have very similar anomaly patterns. During the NEM, there are projected reductions in DIN (Fig. 8f) and associated reductions in primary production (Fig. 8h), of up to 0.2 mmol N m − 3 and 0.2 g C m − 2 day − 1 respectively, from 5 o S-5 o N along the coast of Tanzania and Kenya, and along the equator. This decline of DIN and IntPP can be attributed to two factors: MLD shoaling (Fig. 8d) and the weakening of upwelling-driven nutrient supply. The NEM is characterised by upwelling-favourable alongshore winds with coastal upwelling playing an important role in sustaining elevated primary production along the Kenyan and Tanzanian coast (Jebri et al., 2020). The Kenyan shelf is also the site of shelf-edge upwelling associated with the position of the SZCZ, located on the NKBs during this Fig. 9. Projected sea surface temperature ( o C), mixed layer depth (m), dissolved inorganic nitrogen (mmol N m − 3 ) and depth-integrated primary production (g C m − 2 day − 1 ) during the Southeast [a, c, e, g] and Northeast monsoon [b, d, f, h] in the EACC box (10-0 o S, 38-55 o E) over the 21st century in NEMO-MEDUSA. The cyan lines represent the decadal averages and the thick black lines represent the range of 2 standard deviations of the mean during the baseline period of 2000-19. (For interpretation of the references to colour in this figure legend, the reader is referred to the Web version of this article.) season (Jacobs et al., 2020a,b). The convergence of the EACC and SC then produces divergent flow away from the coast into the Indian Ocean, leading to the upwelling of nutrient-rich water and enhanced productivity (0.25 mmol N m − 3 and 0.55 g C m − 2 day − 1 respectively; see Supplementary Figs. 7f and 7h), which is thought to sustain offshore fisheries (Maina and Osuka, 2014). There are projected changes in the strength of these currents during the NEM (Fig. 6b), indicating a general weakening of the coastal upwelling along Kenyan and Tanzanian coastlines as well as a change in the strength of the divergent flow or a migration of the SZCZ. However, the exact mechanism remains uncertain due to still insufficient model resolution (of ¼ • ) (Fig. 4b) and models with resolution of ~1/12 • are required to adequately reproduce the dynamics of this feature (Jacobs et al., 2020a,b).
The projected reductions in DIN and IntPP in the EACC area is greatest during the NEM (− 49% and − 30% during February; Supplementary Fig. 8). There is also a 17% reduction in IntPP in July (− 26% DIN), which is lower during the other months of the SEM. Additionally, a decline in productivity is projected for the first inter-monsoon period (Mar-Apr), but the model projects a slight increase in productivity during the second inter-monsoon period (Oct). This is in accordance with the lengthening of the SEM period, as seen in Fig. 5c. A similar pattern is projected for DIN concentration in this region. Fig. 9e-h further illustrate the decline in productivity and DIN concentration during both monsoon seasons. However, during the SEM, the trend is weaker and by the end of the century DIN concentration and IntPP are still within their range in the baseline period (within 2 standard deviations of the mean; reductions of 12% for both by the 2090s). In contrast, the negative trends during the NEM are large enough so that the DIN concentration and IntPP exceed their baseline range of variability regularly from the 2050s, with a − 49% and − 27% change respectively by the end of the century. This is consistent with weaker winds, and consequently, less mixing and weaker upwelling, which reduces the transport of nutrients that into the upper mixed layer.
It should be noted that although the trends are weaker during the SEM, greater projected declines are found north of 5 o N and adjacent to the Somali coast with decreases in primary production of up to 0.2 g C m − 2 day − 1 (Figs. 8e) and 0.3 mmol N m − 3 in the DIN concentration (Fig. 8g). This projected decline in the major Somali upwelling system is outside the scope of this paper and a subject of future investigation.

Onset of rising temperature impacts on key commercial species
Both the estimated rise in SST and reduced productivity will affect marine ecosystems in the TWIO. Temperature will change suitable habitat distribution and productivity will alter trophic functioning and structure (Bryndum-Buchholz et al., 2019). The differing seasonal trends may also impact the life cycle of species, especially recruitment.
Of particular concern in the TWIO is the reaction of tuna, which use the this region for spawning and feeding (Dhurmeea et al., 2016). Reduced upwelling and associated productivity will likely result in a general decrease in availability, and therefore catch, of tuna in some areas. For example, Lan et al. (2013), in a study on the effects of climate variability on the distribution and fishing conditions of yellowfin tuna (Thunnus albacares) in the WIO, found that during positive Indian Ocean Dipole events, SSTs were relatively higher, productivity was lower, there was a drop in the catch per unit effort (CPUE) and the catch area was restricted to the northern and western margins of the WIO. Additionally, during a particularly strong event in 1997/98, purse-seine tuna fishing fleets underwent a major shift towards the eastern Indian Ocean (Le Blanc and Marsac, 1999), highlighting the vulnerability of this fishery to future changes in climate (Moustahfid et al., 2019).
Rising SST will likely cause major stress to many fish species, especially where it increases above their thermal tolerances (Dueri et al., 2016). Over time, rising SST is expected to cause the expansion of unsuitable habitats and shifts in the distributions of species both globally  and in the TWIO (Bell et al., 2016;Moustahfid et al., 2019). Increasing SST is expected to be one of the major drivers in the degradation of tuna habitat, particularly for skipjack tuna (Katsuwonus pelamis), which dominates tuna catches worldwide (Parker et al., 2015), with 20% of the total skipjack tuna catch coming from the Indian Ocean (Dueri et al., 2014). The equatorial Indian Ocean has been identified as one of the regions most affected by rising SST in CMIP5 models (Dueri et al., 2014) with considerable warming also found in this region in NEMO-MEDUSA (Figs. 8-10a, b). Although tuna are, to some extent, able to mitigate the higher SST by spending more time at greater depths, feeding will necessitate movement nearer to the surface. However, this may not apply to the abundant skipjack tuna, which have a much lower tolerance to the reduced oxygen levels associated with deeper water than other tuna species (Lehodey et al., 2011).
The upper estimated thermal tolerances of skipjack tuna and kawakawa (Euthynnus affinis) are 30 • C and 30.5 • C, respectively (Lehodey et al., 2011;Boyce et al., 2008). Tuna are large migratory species with extended geographical ranges, such that habitat compression, especially due to changes in SST, are likely to cause changes to navigation corridors and catch potential (Mugo et al., 2010). With a complex SST distribution in the TWIO, and heterogeneous warming, an exact estimation of how the distribution of suitable thermal habitat for tuna species is likely to change with continued greenhouse gas emissions is a challenge. To aid such estimates, Fig. 10 shows when SST is likely to exceed 30 • C for at least 6 months of the year, and for every month of the year. The majority of the surface waters of the TWIO become partially unsuitable thermal habitat for skipjack and kawakawa by the 2050s, with the Mozambique Channel and offshore equatorial regions being partially unsuitable even now. Given that the upper temperature limit for abundant occurrences of skipjack tuna in the Pacific Ocean is 29 • C (Lehodey et al., 2011), Fig. 10 also shows that reasonable catches of the commercially important skipjack tuna are likely to be disrupted by the warming ocean in the 2040s. In short, the expansion of unsuitable thermal habitat for skipjack and kawakawa in the coming years is expected to lead to significant changes in migration patterns, and a deepening of suitable habitats within the water column. These changes are expected to alter regional catch potential for these species and necessitate new strategies for their management and conservation (Bernal et al., 2017;Moustahfid et al., 2019).

Marine heat waves
In addition to the long-term annual average SST increase, marine ecosystems also suffer from short-term, localized extreme events (e.g. Hughes et al., 2017). Marine heat waves (MHWs) are periods where ocean temperatures are anomalously warm for an extended period of time, specifically when SST exceeds the 90th percentile of the climatological period for more than 5 days . Here, we examine the evolution of MHWs for the TWIO over the 21st century, as projected by NEMO-MEDUSA under RCP8.5. A fixed climatological baseline period is calculated from 1995 to 2014 using daily averages (interpolated from 5-day averaged model output). When using this approach, care must be taken when interpreting the results given the positive SST trend (Fig. 9), which makes the identification of MHWs sensitive to the choice of the baseline period. For reasons of consistency with other studies we have chosen a fixed baseline period, although an argument can be put forward that a moving baseline period may provide a more informative metric. What may count as a heatwave in a current decade is likely to become a new normal in the next one and re-definition of the baseline might need to be considered.
As shown by our model (Supplementary Fig. 9), the projected frequency and duration of MHWs in the TWIO increases in the first half of the 21st century until they converge to a continuous event lasting 1 year or more. This trend is consistent with the global projections from the CMIP5 archive (Oliver et al., 2019). One of the quantitative indicators of the evolution of the spatial patterns of MHWs is a date when MHWs are projected to last an entire year, i.e. becomes a continuous event. As shown in Fig. 11a, the onset of year-long heatwaves is heterogeneous. Along the Tanzanian and Kenyan coastlines, this occurs from 2035 to 2045 but for the offshore Northwest Indian Ocean it is projected to occur as early as the 2020s. For the East African coastlines, this will have severe consequences for the shallow coral reef ecosystems with acute thermal stress to corals likely to result in mass bleaching and mortality as documented in 2016 (Hughes et al., 2017). However, in the path of the SEC and NEMC the onset is delayed to the 2050s with a localized area at the northeast side of Madagascar showing onset as late as the 2060s. A delayed onset also evident in the upwelling systems such as Somali, NKBs and Tanzanian upwellings.

Ocean acidification and de-oxygenation
Two other major stressors of marine ecosystems are ocean acidification and deoxygenation which we only briefly mention here as the general patterns have been discussed elsewhere (Popova et al., 2016;Jiang et al., 2019). Ocean acidification arises from increased ocean uptake of carbon dioxide from the atmosphere, which has been exacerbated by man-made emissions (e.g. Gruber et al., 2012;Bopp et al., 2013). Much like climate change stressors, the effects of ocean acidification depend on region-specific life-history characteristics of species and populations, with species directly dependant on ocean carbonate chemistry for calcified skeleton maintenance (e.g. corals and shellfish), poor internal acid-base regulation ability, and long pelagic larval stages being the most vulnerable (Calosi et al., 2017;Kroeker et al., 2013;Widdicombe et al., 2008). In the WIO, decreased ocean pH is likely to affect important habitat forming, keystone species such as corals, especially under concurrent long-term warming and increased heat wave conditions (Kroeker et al., 2013). Resulting loss of suitable habitat conditions for carbonate reef maintenance, decreased growth and increased mortality could affect not only coral but associated reef communities, which support socially and culturally important, artisanal, reef dependent wild capture fisheries in the region (Samoilys et al., 2017). A general decline of surface pH by ~0.3 by the end of the century is projected by NEMO-MEDUSA under RCP8.5, in agreement with the CMIP5 models.
Ocean deoxygenation is the ongoing general decline in ocean oxygen concentrations Breitburg et al., (2018). It is driven primarily by lower oxygen solubility under higher surface temperatures manifesting itself as the expansion of oxygen minimum zones across the world's oceans as well as a general decline of the oxygen concentrations. The oxygen minimum zones are regions in which oxygen saturation in seawater is at its lowest. These zones occur at depths of between 200 and 1500 m, and are most pronounced along the western coasts of continents, mostly in the tropics (e.g. Karstensen et al., 2008). Unlike the Arabian Sea and the Bay of Bengal, the TWIO does not currently include any of the major oxygen minimum zones. However, it already has low oxygen concentrations in many coastal systems, in particular around the vicinity of the coast of Somalia. Expansion of these zones, as well as the general decline of oxygen concentrations, alongside ocean warming and acidification, will result in a major decline in the health of marine ecosystems, alter species distributions, modify key climate regulation pathways such as carbon sequestration, especially in the inner ocean, and thus potentially lead to serious impacts on the economies dependent on these ecosystems (Ravaglioli et al., 2019;Levin et al., 2020).

Discussion and conclusions
Although climate change is a global phenomenon, its manifestation in many regions is often unique and shaped by local factors. The TWIO region is especially prominent in this regard, with its distinctive dynamics being shaped by the unique geography with the presence of a landmass at the northern boundary driving seasonally reversible monsoonal winds and, consequently, seasonally reversible ocean currents (Schott and McCreary, 2001). In line with the topic of this special issue, we focus on the region influenced by the EACC flowing along the continental margins of Kenya and Tanzania. The dynamics of these tropical waters are influenced by three key currents: the NEMC, the EACC and the seasonally reversing Somali current with numerous upwelling sites along their advective pathways, which in turn influence catches in local fisheries (Jury et al., 2010;Jacobs et al., 2020a,b;Jebri et al., 2020;Painter et al., under review).
Monsoons are clearly the overall drivers of the dynamics of these currents and the associated ecosystem response, and are one of the key uncertainties in future climate projections (IPCC, 2013;Lee and Wang, 2014;Sharmila et al., 2015). AR5, and even AR6, models are still only at 1 • spatial resolution (Seferian et al., 2020) and cannot resolve intricacies of the dynamics driving key features of productivity in these waters. Similarly, the resolution of the ocean models in the climate downscaling effort for Africa is also about 1 • (e.g. de Castro et al., 2016). Here, we have used a high-resolution (1/4 • ) ocean model with surface forcing from a climate model to better resolve these key regional features under a relevant climate scenario. Such a model allows us to look beyond the large-scale patterns of temperature increase, stabilisation of stratification and decline of primary production identified by the CMIP5 analysis (Bopp et al., 2001Kwiatkowski et al., 2017).
Taking into account seasonal monsoonal variability, we have analysed the NEM and SEM periods separately. We find the following key projected changes, which are summarised in Fig. 12: • A weakening of NEM winds (− 19%) and a slight strengthening of SEM winds (+6%) by the end of the century, which is consistent with other CMIP5 model projections for the Indian Ocean. • Shallower mixed layers during the NEM (− 31%) and weakening of upwelling associated with the EACC, with a subsequent reduction in surface DIN concentrations (− 49%) and integrated primary production (− 27%). • Relatively stable dynamics of primary production during the SEM, showing only minor change throughout the course of the century. • Widespread warming of up to 5 • C year-round with an initial increased frequency and duration of marine heatwaves with the first year-long event projected to occur as early as the 2030s. • Weakening of the NEMC by up to − 26% in the early SEM (May-June) causing reduced connectivity and hence, reduced nutrient and pelagic larvae supply, between Madagascar and the East African coast.
These projected changes will have consequences for the marine ecosystems and socio-economics pertaining to the coastal population of Tanzania and Kenya who are dependent on marine resources for food security (Taylor et al., 2019). The considerable warming of this region projected under RCP8.5 is in line with global (e.g. Power et al., 2017) and tropical estimates (e.g. Mizuta et al., 2014). This will have significant negative impacts on the reef fisheries where widespread coral bleaching has already been reported for the region (e.g. Ateweberhan et al., 2011;Gudka et al., 2020) during recent decades. For example, SST exceeded the bleaching threshold for months at a time from 2014 to 17 at the Aldabra atoll in the TWIO, causing a >90% decline in soft coral cover (Cerutti et al., 2020). A steady rise in SST, combined with extreme events such as MHWs, will lead to negative impacts not only on the coral reef-based ecosystems and fisheries (Cowburn et al., 2018) but also on the wider marine ecosystems in the region (Frolicher and Laufkotter, 2018). There is already evidence that tuna species are shifting in response to climate stressors in the TWIO, which may have knock-on effects for the coral reef fisheries (Moustahfid et al., 2019).
Of great concern is also the projected increase in MHW events. These are generated by the combined effects of atmospheric and oceanic processes  and, if prevalent, can lead to potentially catastrophic consequences for marine ecosystems causing shifts in the geographic distribution of species, changes in phenology and widespread mortality (e.g. Wernberg et al., 2013;Moustahfid et al., 2019;Thomsen et al., 2019). MHWs have recently become a focus of research as the most visible signs of anthropogenic climate change impacts on ocean ecosystems. Here, we find that the majority of the TWIO will have experienced its first year-long heatwave by the 2030s. This is a significant result especially given the overall SST increase mentioned above. However, we also note that along the path of the SEC, the NEMC and key locations of regional upwelling along the shelf area, this could happen as late as the 2050s, suggesting that these regions could be a refuge for some species from MHWs. In the case of productivity, we find that the projected weakening of the NEMC during the early SEM, also projected by the CMIP5 models (Stellema et al., 2019), reduces the connectivity to the East African coast. This could reduce the amount of nutrients reaching the Tanzanian coast, and consequently, cause a reduction in primary production to the detriment of local pelagic fisheries (e.g. Anderson et al., 2016;Sekadende et al., 2020). Productivity in the EACC region is also dependent on the seasonal North Kenya Banks (NKBs) upwelling during the NEM (Jacobs et al., 2020a,b) and the dynamic uplift upwelling that occurs along the fast-flowing EACC during the SEM (Jebri et al., 2020). Here, we find the projected offshore shift of the EACC during the SEM would reduce the availability of nutrients, and hence limit primary production near the coast.
Shifting currents are known to cause shifts in the distribution of marine species (e.g. Coleman et al., 2013;Cetina-Heredia et al., 2015) and, in this case, could be especially detrimental for small pelagic fisheries in the TWIO (Sekadende et al., 2020). These fish are highly sensitive to environmental conditions and show specific tolerance to temperature, salinity, oxygen and pH (Alheit et al., 2009). However, due to the complex interaction of currents in this region, it is less clear how the upwelling on the NKBs, controlled by the position of the SZCZ, will be affected in future. Modelling this interaction will require another increase in resolution to 1/12 • . Any persistent shift of the SZCZ is likely to affect the offshore fishery (Maina and Osuka, 2014), which is believed to be sustained by the associated upwelling and enhanced productivity. Higher resolution model projections are required to fully unravel the changes in the local dynamics of this region, which is proposed as a new frontier for food security.
For the wider region, the projected decline in primary production using NEMO-MEDUSA is consistent with that found in the CMIP5 models, i.e. up to − 30% in the tropics (e.g. Bopp et al., 2013;Kwiatkowski et al., 2017). We note in the northwest Indian Ocean, a region with significant upwelling including that off Somalia, that primary production is projected in other studies to have greater declines than the TWIO (e.g. Roxy et al., 2015;Nakamura and Oka, 2019). This outcome is also indicated in our study (Fig. 8) but unraveling the mechanisms behind this decline is beyond the scope of this paper.
The contrasting NEMO-MEDUSA model projections for each monsoon season, associated with opposing trends in the regional wind field, highlight the need to investigate seasonal changes as opposed to solely studying annual and decadal averages for this region. In addition, the critical importance of regional upwelling system in maintaining primary production underscores the importance that future change may play in these highly dynamic systems. However, completely representing these features, as well as the wider complicated regional dynamics of the TWIO, requires the use of even higher resolution (1/12 • ) models than currently available, even in the case of ocean-only modelling systems such as that used here.

Implications for policy and management
Tropical marine fisheries are highly vulnerable to climate change impacts (Lam et al., 2020). This vulnerability drives an urgent need to understand the risk that climate change puts on the marine ecosystems, the associated fisheries and the socio-economic systems dependent on them. Understanding the risk will better inform climate change adaptation measures.
Although model projections of the key climate change stressors of marine ecosystems at regional scales are improving, understanding, documenting and predicting the response of the ecological and socioeconomic systems to these stressors is an equally challenging task. This requires specific observational and monitoring strategies and documentation of the system-level responses to climatically driven events that fall outside of natural range of variability. Projections of such system-level responses, which take into consideration specifics of the region, should be at the core of determing the climate change adaptation options.
An example of such a system already in place is the NOAA Coral Reef Watch (https://coralreefwatch.noaa.gov/). This is a coral bleaching alert system, based on the daily monitoring of satellite remote sensing SST data. The near-real time and long-term monitoring system assists the global community in the management, study, and assessment of warming events on coral reef ecosystems brought about by environmental change (Liu et al., 2013). Similarly, systems for monitoring, analysis and assisting responses to MHWs are also beginning to emerge (Holbrook et al., 2020). Nonetheless, similar systems for other climatic stressors are needed. For instance, a growing body of evidence suggests that upwelling systems and their variability play a key role in supporting regional ecosystems and fisheries (Jury et al., 2010;Jacobs et al., 2020a, b;Jebri et al., 2020;Kamau et al., 2021). For the EACC region in particular, near-real time monitoring of upwelling systems and the associated phytoplankton blooms (based on remotely sensed information in combination with improved fisheries catch statistics) is likely to provide a significant step towards understanding a system-level response to climate impacts on marine ecosystems.
We draw attention to our findings here and provide a summary of the primary EACC region-specific stressors pertaining to ocean currents and upwelling systems that could guide development of a regional monitoring system. Key components include: • Provision of near-real-time warning of anomalous behaviour of key regional upwelling systems (e.g. suppression or intensification of upwelling, high SST, low Chl-a via near real-time satellite remote sensing). • Monitoring of the dynamics of the NEMC and EACC using near realtime satellite altimetry data. • Documenting responses of ecosystems to events of anomalous strengthening or weakening of the coastal and shelf-break upwelling systems. • Putting in place a regional MHW alert system. With the exception of the coral-based ecosystems, the response of local fisheries to MHWs is poorly documented in the region, so observations of how fisheries respond to increasing MHW strength and duration will provide valuable information on the ecosystem response to the first manifestation of anthropogenic climate change. However, whether a stand-alone system or one embedded in existing platforms, like NOAA's Coral Reef Watch, would be beneficial to states in the TWIO, very much depends on the frequency of MHWs as well as the translation of scientific information into policy.

Modelling for climate risk management of fisheries: the next steps
Long-term climate change impacts are driving widespread changes in the distributions of marine species, with observed poleward migrations tracking change in the distribution of suitable habitats in many regions Polockzaska et al., 2013). Projections of climate change impacts on the marine environment similar to those described in this study, provide an understanding of the level of exposure of species to changing environmental conditions that, when combined with information on species sensitivity, aid th development of conservation and other ocean management policies that can help mitigate impacts on natural and associated human systems (IPCC, 2013). However, due to insufficient resolution of models and problems associated with the consequences of a long millennia-scale equilibration of the forcing fields (e.g. Yool et al., 2020), these projections are less useful at the spatial and temporal scales needed for management and policy decisions regarding climate change risk management for regional and local fisheries, and especially to artisanal fisheries at the core of communities with less adaptive capacity. Such decisions require analyses of climate change variations on timescales shorter than those often associated with climate projections (typically 100 years). Hobday et al. (2018) demonstrated a successful use of seasonal (up to three months) forecasting in decision-making for a range of fisheries and aquaculture-based industries in Australia. However, for developing climate change adaptation options this timescale is too short, and progress needs to be achieved through investment into decadal prediction methodsa time scale that is more relevant to regional policy making. Significant progress has been made recently in decadal prediction as seen in Smith et al. (2019). Particularly relevant to the TWIO is the finding that the Indian Ocean is the region of highest skill worldwide in decadal climate prediction in terms of SSTs (Guemas et al., 2013). However, this skill is related primarily to the long-term warming trend in the Indian Ocean (Guemas et al., 2013, Monerie et al., 2012. This means that once this trend is accounted for, the reliability of the decadal prediction decreases (Corti et al., 2012). Essentially this shows that it is the external forcing rather than internal variability that accounts for the skill in Indian Ocean decadal predictions. Nevertheless, this does suggest that decadal predictions of SST might help in determining the shorter-term behaviour of marine heat waves, and potentially their spatial variations. Therefore, decadal predictions could usefully contribute to policy and decision making for the TWIO with regard to protecting ecosystems that are at risk. In addition, recent progress has been made in the multi-year prediction of the Indian Ocean Dipole (IOD; Feba et al., 2020), which is a major influence on the SST and productivity of the TWIO.
Decadal prediction to-date has frequently focussed largely on the physical climate system (e.g. the Decadal Climate Prediction Project www.wcrp-climate.org/dcp-overview) and there remains the challenge of developing decadal prediction systems for biological and biogeochemical ocean processes more directly relevant to ecosystem management, and for the ecosystems themselves. Such systems are currently beginning to emerge both in the conservation and the fishing sector (Maxwell et al., 2015). For instance, Payne et al. (2017) reviewed progress has been made in producing ecological forecasts, but these forecasts are on near real time to seasonal time scales, rather than decadal ones. In terms of decadal predictions of biological and biogeochemical processes the first significant work was that of Séférian et al. (2014), who demonstrated multiyear predictability of net primary production (NPP) in the tropical Pacific, emphasising a strong link to large scale climate variations such as ENSO. The study further showed that the predictive skill for NPP was 3 years in contrast to that for SST of 1 year, the higher NPP skill being attributed to the advection of nutrient anomalies. Given that productivity in the WIO is closely linked to multi-decadal oscillations (Jury et al., 2010;Jebri et al., 2020) and that nutrient advection can play a part in productivity in the TWIO (Jebri et al., 2020) this suggests that NPP might be predictable here too at similar timescales. Tommasi et al. (2017) reviewed the role of seasonal to decadal climate forecasts in managing marine resources, but all the case studies they describe are at seasonal or shorter timescales. More recent work (e.g. Frölicher et al., 2020;Krumhardt et al., 2020) has shown that there is potential predictability in both NPP and other marine ecosystem drivers (such as temperature, pH, and oxygen) on time scales of the order of 3 years. Park et al. (2019), in a hindcast study, show that annual fish in the Agulhas and Somali large marine ecosystems can be predicted from chlorophyll and SST anomalies up to 2-3 years in advance (despite the 1 • model resolution being less than ideal for modelling coastal large marine ecosystems). These results show that decadal forecasts, together with the long-term projections of climate change on the marine environment, may in the near-future enable development of adaptation strategies for marine fisheries and communities dependent on the health and stability of local ecosystems in the TWIO.

Declaration of competing interest
The authors declare that they have no known competing financial interests or personal relationships that could have appeared to influence the work reported in this paper.