The future of coastal upwelling in the Humboldt current from model projections

Abstract

The Humboldt coastal upwelling system in the eastern South Pacific ocean is one of the most productive marine ecosystems in the world. A weakening of the upwelling activity could lead to severe ecological impacts. As coastal upwelling in eastern boundary systems is mainly driven by wind stress, most studies so far have analysed wind patterns change through the 20th and 21st Centuries in order to understand and project the phenomenon under specific forcing scenarios. Mixed results have been reported, and analyses from General Circulation Models have suggested even contradictory trends of wind stress for the Humboldt system. In this study, we analyse the ocean upwelling directly in 13 models contributing to phase 5 of the Coupled Model Intercomparison Project (CMIP5) in both the historical simulations and an extreme climate change scenario (RCP8.5). The upwelling is represented by the upward ocean mass flux, a newly-included variable that represents the vertical water transport. Additionally, wind stress, ocean stratification, Ekman layer depth and thermocline depth were also analysed to explore their interactions with coastal upwelling throughout the period studied. The seasonal cycle of coastal upwelling differs between the Northern and Southern Humboldt areas. At lower latitudes, the upwelling season spans most of the autumn, winter and spring. However, in the Southern Humboldt area the upwelling season takes place in spring and the summertime with downwelling activity in winter. This persists throughout the Historical and RCP8.5 simulations. For both the Northern and Southern Humboldt areas an increasing wind stress is projected. However, different trends of upwelling intensity are observed away from the sea surface. Whereas wind stress will continue controlling the decadal variability of coastal upwelling on the whole ocean column analysed (surface to 300 m depth), an increasing disconnect with upwelling intensity is projected below 100 m depth throughout the 21st Century. This relates to an intensification of ocean stratification under global warming as shown by the sea water temperature profiles. Additionally, a divergence between the Ekman layer and thermocline depths is also evidenced. Given the interaction of upwelled nutrients and microscopic organisms essential for fish growth, a potential decline of coastal upwelling at depth could lead to unknown ecological and socio-economical effects.

1 Introduction

The eastern boundaries of the Pacific and Atlantic Oceans at subtropical latitudes encompass large areas of coastal upwelling (Sydeman et al. 2014). This upwelling activity is characterized by a relatively shallow vertical transport of mass of water sourced from about 300 m depth (Talley et al. 2011), which is dynamically independent from the ocean gyres systems governing subtropical ocean circulation. An offshore Ekman transport forced by equatorward upwelling-favorable winds has been long accepted as the main driver of the greatest eastern coastal upwelling systems on the Earth (Sverdrup 1938). The offshore surface transport is replaced by water upwelled from depth (Neelin 2009) leading to the emergence of cold, nutrient-rich waters at the surface (Bakun 1990; Warwick and Marjorie 2001; Diffenbaugh et al. 2003; Bakun et al. 2010, 2015; Sydeman et al. 2014). As these systems play a significant role for marine life and hence the fishing industry, variations of coastal upwelling intensity have direct ecological and economic impacts (Vecchi et al. 2006; Bakun et al. 2010; Sydeman et al. 2014; Rykaczewski et al. 2015).

Although the importance of coastal upwelling systems is widely recognised, their future behaviour is still uncertain (Bakun et al. 2010). The so-called Bakun (1990) hypothesis, based on observational wind data from ship records during the second half of the 20th Century, expects an intensification of coastal upwelling in the major coastal upwelling systems along the subtropical zones as global warming increases. The mechanism for this is a seasonal strengthening during the warm season of upwelling-favorable winds along-shore, due to an increasing atmospheric pressure gradient driven by a large-scale exacerbation of the land-sea thermal contrast. Sydeman et al. (2014) applied logistic regression methods on observational data from the 20th Century in order to analyse upwelling-favourable wind intensification, in terms of velocity, during the warm season. In accordance with the Bakun (1990) hypothesis, with the exception of the Iberian system, they observed an intensification of wind stress during the warm season in all the major coastal upwelling systems on the planet, collectively called the Eastern Boundary Current System (EBCS). This result is similar to other studies based on observed data and modelling outcomes throughout the 20th Century arguing a potential strengthening of upwelling-favourable winds in coastal upwelling areas under future climate change (Cane 1997; Diffenbaugh et al. 2003; Garreaud and Falvey 2009).

In this research, we focus on the Humboldt upwelling system, which is one of the most productive systems in the world (Daneri et al. 2000; Bakun et al. 2010; Taylor et al. 2008; Garreaud and Falvey 2009; Oerder et al. 2015) and provides between 15 and 20% of the world’s fish extraction (Sherman and Hempel 2008). Several studies have identified an intensification of the upwelling-favorable winds throughout the warm season in the Southern Hemisphere (Bakun et al. 2010; Belmadani et al. 2013; Sydeman et al. 2014; Rykaczewski et al. 2015). Based on 15 global climate models contributing to CMIP3, Garreaud and Falvey (2009) suggest a strengthening of the upwelling-favourable winds in the Humboldt system at the end of the 21st Century. A similar trend has also been reported by in climate model simulations in CMIP5 (Rykaczewski et al. 2015; Wang et al. 2015). However, different spatial patterns of this intensification were found during the warm season. Whereas wind speeds increase in the southern upwelling zone, a weakening of upwelling-favorable winds during summertime is projected nearer the Equator. The latitudinal response of upwelling-favourable winds has also been described at local scales in the Humboldt system (Aravena et al. 2014).

A debate about a possible weakening of the Humboldt system during the last century and future projections has also taken place. Whereas some early studies argued possible errors in wind estimation methods and the equipment used for direct measurements (Peterson and Hasse 1987; Ramage 1987; Wright 1988; Cardone et al. 1990), other research based on observed data and modelling outcomes (Hsieh and Boer 1992; Clarke and Lebedev 1996; Vecchi et al. 2006) directly suggest a potential decline of the Humboldt upwelling system. Clarke and Lebedev (1996), using observational data of surface atmospheric pressure as a proxy for wind stress, identified weakening episodes of the equatorial trade winds in the Pacific during the 20th Century, especially from 1960. These results were consistent with surface temperature data reported off the Peruvian coast during the same period. Similarly, Vecchi et al. (2006) carried out a study based on observational data and GCM outcomes from the GFDL CM2.1 model. The study suggested a weakening of trade winds in the equatorial Pacific as GHG emissions increased, with a consequent decline of the Pacific equatorial upwelling and a potential future decline of the Humboldt upwelling system. This last assertion has been under debate. A disconnect between the northern Humboldt upwelling system and the easterly trade wind dynamics in the tropical Pacific has been argued (Strub et al. 1998). Additionally, as Belmadani et al. (2013) pointed out, the apparent weak connection between these two systems could partially explain opposing trends in 20th Century observations of the intensity of tropical and subtropical upwelling-favorable winds reported by Bakun and Weeks (2008).

1.1 Coastal upwelling in coupled atmosphere-ocean GCM (AOGCM)

Some relevant processes related to coastal upwelling are not fully represented in AOGCMs, because of model resolution (Boville and Gent 1998; Mote and Mantua 2002; Burroughs 2007; Palmer 2014). For instance, the inaccurate representation of marine stratus clouds in models may lead to errors in wind stress simulations (Mote and Mantua 2002). In the ocean, mesoscale eddies have a direct effect on water properties and the local current dynamics on a spatial-scale smaller than regular AOGCM ocean grids (Griffies 2004; Flato et al. 2013). In addition, because of the interaction with complex mixing processes, the estimation of mesoscale eddies in boundary regions is especially hard (Griffies 2004).

Different spatial resolutions in AOGCMs have led to substantially different (and even contradictory) outcomes between models. The studies by Hsieh and Boer (1992) and Mitchell et al. (1990) based on wind velocity and pressure gradient outcomes suggested a weakening of the major coastal upwelling systems. However, differences have been obtained as the spatial resolution of GCM has improved. While the spatial resolution in the model used by Hsieh and Boer (1992) was $3.75^{\circ }$ (roughly 300–400 km), Mote and Mantua (2002) based on two GCMs with resolution of about 200 km reported no significant variations of wind stress for the EBCS. Nonetheless, recent studies based on modelling approaches at finer resolutions do not necessarily suggest an opposite trend to earlier models. Indeed Vecchi et al. (2006), based on a single GCM (GFDL CM2.1) with a resolution of $2.0^{\circ }$ latitude and $2.5^{\circ }$ longitude, suggest a potential decline of the Humboldt upwelling system. On the other hand, the 23 coupled GCMs from the CMIP3 dataset with varying resolutions between $0.3^{\circ }$ and $5.0^{\circ }$ reported mixed results for the Humboldt system without a clear trend (Wang et al. 2010). However, there is a scarcity of coastal upwelling studies based on GCMs (Miranda et al. 2012) and the apparently opposing results from models with different resolutions might be partly explained by the scale of analysis rather than being contradictory projections (Bakun et al. 2010).

1.2 Wind stress and upward ocean mass transport

There is consensus regarding the role of wind stress as a main driver of coastal upwelling. Therefore, modelling and observational data-based approaches developed to date have mainly been focused on upwelling-favorable wind pattern analyses. Different approaches to conceptualize wind stress in climate models have been adopted. Some studies have considered surface wind velocity (Garreaud and Falvey 2009; Echevin et al. 2011; Sydeman et al. 2014), upwelling index based on wind stress and the Ekman mass transport (Aravena et al. 2014), surface atmospheric pressure gradients as a proxy for winds stress (Clarke and Lebedev 1996; Vecchi et al. 2006) and also estimate of wind stress from the time-averaged eastward (u) and northward (v) velocity components of wind (Yelland and Taylor 1996; Mote and Mantua 2002; Bakun et al. 2010; Varela et al. 2015).

As upwelling dynamics is hard to measure practically, there are not enough data from observations and most of the upwelling phenomenon has been inferred from Global Circulation Models (Talley et al. 2011) or indirect estimations from climate-ocean variables, such as wind stress patterns. In addition, upwelling-related ocean variables such as the sea surface temperature (SST) do not necessarily represent a proper proxy for time variability of coastal upwelling (Lorenzo et al. 2005; Chen et al. 2012). Lorenzo et al. (2005) reported positive correlations of upwelling-favorable wind intensity and SST in the California upwelling system. According to the authors, this apparently unexpected relationship might be related to a large-scale heat flux from atmosphere-ocean interactions. Indeed SST is not only affected by ocean upwelling, but also depends on complex ocean-atmosphere interactions (McCabe et al. 2015; Tim et al. 2015).

The upward ocean mass transport variable ($\hbox {kg}\,\hbox {s}^{-1}$) included in some of the AOGCMs contributing to CMIP5 represents a direct output from models, which can be used as a measure of upwelling activity. As Griffies et al. (2010) stated, this variable is an adequate indicator for analysis in part because other physics variables such as vertical velocity can be also derived from it. However, in CMIP3 models only the vertical velocity was included and there is a scarcity in studies of upwelling based on the upward ocean mass transport variable and its relationship with other climate-ocean variables related to coastal upwelling. Although a recent study by Xu et al. (2014) in the Southeast Atlantic Ocean incorporated the analysis of the upward ocean mass transport, this was focused in sea-surface temperature biases rather than exploring upwelling.

This paper aims to explore the coastal upwelling activity in the Humboldt system through the analysis of the upward ocean mass transport variable and its relationship with wind stress and ocean stratification throughout the Historical experiment (1850–2005) and the Representative Concentration Pathways 8.5 (RCP8.5) climate change scenario (2006–2100). In addition, the relationship between Ekman layer and thermocline depths is also explored throughout the period of study.

2 Methods

The CMIP5 ensemble offers a substantial improvement on CMIP3 models regarding not only improved resolutions but also a richer set of output data, especially for ocean variables (Griffies et al. 2010). An example is the upward ocean mass transport variable ($\hbox {kg}\,\hbox {s}^{-1}$) explored in this study. Surface meridional wind stress ($\hbox {N}\,\hbox {m}^{-2}$) on the upper level of the ocean grid was used as a simplification of the resultant wind stress given the coastal geometry of South America’s Pacific Coast. Additionally, ocean stratification was analysed through the sea water temperature profiles. For the purposes of this study, output from 13 AOGCMs contributing to CMIP5 (Table 1; determined by availability of data for the variables of interest) were examined throughout the Historical (Taylor et al. 2012) and RCP8.5 (Riahi et al. 2011) long-term experiments. The period considered in the analysis spanned from January 1851 to December 2005 for the Historical experiment, and January 2006 to December 2100 for RCP8.5.

Given the differing expression of the seasonal cycle, the analysis was carried out separately for the higher and lower latitudes in the Humboldt System (Aravena et al. 2014; Rykaczewski et al. 2015; Varela et al. 2015). Two areas of interest were defined over the coastal upwelling region, the Northern Humboldt, spanning between $16^{\circ }\hbox {S}$ and $28^{\circ }\hbox {S}$ latitude, and the Southern Humboldt spanning between $29^{\circ }\hbox {S}$ and $40^{\circ }\hbox {S}$ latitude. In order to define the western limit, a latitudinal profile of upwelling was built for the mid-latitude of each area (not shown). According to the latitudinal profiles, for each model similar longitudes around $75^{\circ }\hbox {W}$ were identified as turning point for coastal upwelling in both areas. Thus, $75^{\circ }\hbox {W}$ longitude was selected as a western meridional limit for the Northern and Southern Humboldt areas. For these areas, the normalised area-average of the upward ocean mass transport (upwelling) and the spatial mean of wind stress and sea water potential temperature were computed. Upwelling and ocean temperature were interpolated to five common different depths: 30, 50, 100, 150 and 300 m. Ocean mixed layer thickness was also obtained from CMIP5 data set.

The estimation of upward transport is based on the formulation by Adcroft and Hallberg (2006) and summarized by Griffies et al. (2010) as follows:

$$\begin{aligned} V^{(s)}= \rho w^{(s)} dx dy. \end{aligned}$$

(1)

$V^{(s)}$ is the vertical mass transport at s-surface, $\rho$ is the in situ or reference density, $w^{(s)}$ is the surface flux in the vertical coordinate and dxdy represents the area of the grid cell. In terrain-following models, the upward mass of water transport is computed by mass-volume conservation, whereas in isopycnal models this is computed by the physics of the mass of water across an isopycnal coordinate. In the CMIP5 dataset, all the hybrid and isopycnic models report upwelling outputs on levels of constant depth remapped from their native coordinates (Griffies et al. 2010).

The formulation of the Ekman layer depth ($D_{ek}$), as described by Cushman-Roisin and Beckers (2011), is based on the turbulence theory where the Ekman number ($E_{k}$) on the order of one is assumed as function of $D_{ek}$, the Coriolis parameter (f) and the eddy viscosity $\nu _{E} =\mu _{*}D_{ek}$.

$$\begin{aligned} E_{k}=\frac{\nu _{E}}{fD_{ek}^2}\approx 1 \end{aligned}$$

(2)

As friction velocity $\mu _{*}$ is defined as $\sqrt{\frac{\tau }{\rho _{w}}}$, $D_{ek}$ can be formulated as follows:

$$\begin{aligned} D_{ek} = \gamma \frac{\sqrt{\frac{\tau }{\rho _{w}}}}{f} \end{aligned}$$

(3)

where $\tau$ is the wind stress and $\rho _{w}$ is water density. Parameter $\gamma$ is a constant of proportionality which according to observations normally adopts a value of 0.4 (Madsen 1977; Garratt 1992; Perlin et al. 2007; Cushman-Roisin and Beckers 2011) which is used in this study. Even though 0.4 is the most accepted value (Perlin et al. 2007; Cushman-Roisin and Beckers 2011), lower values have also been applied for specific ocean settings (Mercado and Leer 1976; Mofjeld and Lavelle 1984; Stigebrandt 1985).

Table 1 Coupled AOGCMs contributed to CMIP5 considered in this study

3 Coastal upwelling in the historical simulations

An important pre-requisite for any insight to be gleaned from future projections is that the models used must capture the necessary physical processes with some fidelity. A lack of confidence in climate models’ representation of upwelling is what has led to the longevity of studies inspired by the Bakun (1990) hypothesis. In this section, we present the upwelling behaviour simulated by models contributing to CMIP5 in the historical experiment.

Despite the importance of the Humboldt coastal upwelling, there appears to be an unfortunate dearth of observations that could be used to evaluate the simulations. Whilst data from many individual cruises exist, it is not clear how well they quantify the large spatial scale resolved by the climate models. One could use ocean reanalyses instead. But as they have an ocean dynamical model at their core, it is unclear whether they would just incorporate the same systematic biases. The highly variable nature of the upwelling, especially on interannual and longer timescales, does not aid a quantitative comparison with observations. We therefore constrict our discussion to the large visual features of the simulations.

All of the models simulate a zone of coastal upwelling, with a downward mass flux further out into the Pacific (Fig. 1). The width of the upwelling zone is strictly confined to the coast of Chile, but extends laterally moving Equatorward along the Peruvian coast. The differences in horizontal resolution of the models is clearly visible (Fig. 1). There is a variation in magnitude of the vertical mass transport which is logically connected to the horizontal resolution (Eq. 1). The resolution has little substantive impact of the spatial patterns however. Furthermore, differences in models formulations may also provide a partial explanation for dissimilar outcomes (Taylor et al. 2012). Despite the above, the spatial distribution of upwelling in all the models depicts a realistic behaviour for coastal upwelling, and similar spatial resolutions tend to represent similar spatial upwelling distribution. More extensive coastal upwelling areas are simulated in coarser models than in the higher resolution ones (Fig. 1).

Although some differences are visible in terms of the magnitude of coastal upwelling between the Northern and Southern Humboldt areas for some finer models, the vertical profile of the upward ocean mass transport (${\hbox {kg}}\, {\hbox {s}}^{-1}$) is consistent in both areas (Fig. 2). Finer models tend to represent magnitudes closer to the ensemble mean or with higher intensity in both regions. Similar to the value proposed by Talley et al. (2011), upwelling corresponds to a relatively shallower vertical flux of water sourced from about 300–9500 m, being the ensemble mean roughly 400 m in both areas. Downwelling is also observable from this depth to the deeper ocean with reduced vertical fluxes as depth increases. At about 50 m depth the magnitude of upwelling is maximum in both areas. However, whereas in the Northern Humboldt the ensemble mean reaches $1.1\, \times \, 10^{-3}\, \hbox {kg s}^{-1}\,\hbox {m}^{-2}$ at this depth in the Southern region the integrated upwelling flux reaches about $0.9\, \times \, 10^{-3}\,\hbox {kg s}^{-1}\, \hbox {m}^{-2}$. Similar magnitude of difference can be appreciable for most of the single models. This stronger coastal upwelling in the Northern Humboldt can be particularly explained by an enhanced Ekman transport in that region influenced by a lower Coriolis force ($f=0.545 \times 10^{-4}$) in comparison to the Southern Humboldt area ($f = 0.834 \times 10^{-4}$). Furthermore, the intensity at large scale of equtorward upwelling-favorable winds in the Humboldt system is affected by its latitudinal location in relation to the South Pacific Anticyclone (SPA), which induces more intense equatorward winds along the western coast of South America from about $35^{\circ }\,\hbox {S}$ latitude.

A different seasonal cycle (Aravena et al. 2014) is simulated for the Northern and Southern Humboldt areas throughout the Historical experiment for the 30-years period between 1971 and 2000 (Fig. 3). The maximum values of upwelling in the Northern Humboldt area (Fig. 3a1) take place between July and November, with peak values in September. On the other hand, minimum upwelling activity occurs in the warm season (DJF). This pattern of a extended coastal upwelling season most of the autumn, winter and spring in low subtropical latitudes in the Humboldt is consistent with wind stress patterns reported by Rykaczewski et al. (2015) and Varela et al. (2015), and the seasonal cycle of the upwelling index computed from wind patterns by Wang et al. (2010). The upwelling season at lower latitudes presents a more stable coastal upwelling activity throughout the year (Bakun and Nelson 1991; Schwing and Mendelssohn 1997; Rykaczewski et al. 2015). This pattern is especially observable at shallower water layers above 150 m depth in the Northern Humboldt area (Fig. 3a1).

In the Southern Humboldt area, coastal upwelling presents a distinctive seasonal behaviour (Fig. 3a2). The upwelling season spans most of the spring and summer (September to March), similar to the wind stress patterns described by Varela et al. (2015) based on reanalysis wind database with $0.3^{\circ }$ of spatial resolution. Maximum coastal upwelling takes place between December and January (except at 300 m depth where the maximum is in October). Minimum upwelling activity takes place throughout the cold season with negatives upward flux (downwelling) at the five depths analysed between May and July. This pattern of coastal downwelling in high latitudes during the cold season is consistent with previous studies based on observed and modelled winds patterns analysis (Bakun and Nelson 1991; Bakun et al. 2010; Rykaczewski et al. 2015; Wang et al. 2015).

Fig. 1

Spatial distribution along the Humboldt system of the upward ocean mass transport from CMIP5 data set during the spring season, 1971–2000. Black squares in each panel are the Northern Humboldt (top) and Southern Humboldt (bottom) areas. Data is vertically interpolated to 100 m depth in native units (${\hbox {kg}}\,{\hbox {s}}^{-1}$)

Fig. 2

Vertical profiles of normalized area-average upward ocean mass transport in the historical simulation. Black lines show the ensemble mean for the period 1971–2000. Individual models are coloured by their horizontal resolution, with red indicating finer and grey coarser models (Table 1)

Fig. 3

Model ensemble mean of the seasonal cycle of coastal upwelling at 30, 50, 100, 150, and 300 m depth for the Historical period (1971–2000) (a1, a2) and the climate change signal (b1, b2) in terms of absolute values computed for 2071–2100 with respect to 1971–2000. Red stars in b1 and b2 represent significant deviation from zero at 95% confidence of the climate change signal (two-tailed Student’s t test)

4 Coastal upwelling under climate change

Despite the dissimilar annual cycle of coastal upwelling in the Northern and Southern Humboldt areas (Fig. 3a1, a2), the climate change signal projected for 2071–2100 has similar patterns in both regions (Fig. 3b1, b2). Under the RCP8.5 (2071–2100) a gradual decline of coastal upwelling in the summertime and an increasing trend in winter is projected at all analysed depths, and the climate signal deviates significantly from zero (95% confidence) at most of the depths in both areas (Fig. 3b1, b2). This seasonal behaviour suggests that upwelling seasonality under climate change will be mainly driven by changes in the wind stress annual cycle. Previously, an increasing trend in the austral winter intensity and a future decline during the warm season has been described for the Humboldt system (Belmadani et al. 2013). However, even though similar patterns of the climate change signal are projected for the whole water column until 300 m depth, as depth increases a more pronounced annual amplitude of the climate change signal is projected throughout the year.

Similar annual amplitudes of the climate change signal are observed in both regions for shallower waters (e.g. 30 m depth). In the Northern Humboldt, whereas in the summertime (December-January) a decrease of around $-0.1 \times 10^{-3}\,{\mathrm {kgs}}^{-1}\,{\mathrm {m}}^{-2}$ is projected at the surface (30 m depth), during the cold season (June-July) an increase by $0.1 \times 10^{-3}\,{\mathrm {kgs}}^{-1}\,{\mathrm {m}}^{-2}$ has been projected. However, at 300 m depth, coastal upwelling would increase by $0.15 \times 10^{-3}\,{\mathrm {kgs}}^{-1}\,{\mathrm {m}}^{-2}$ in June, whereas in December would decrease by around $-0.4 \times 10^{-3}\,{\mathrm {kgs}}^{-1} \,{\mathrm {m}}^{-2}$. This higher sensitivity under an extreme climate change scenario is also observable in the Southern Humboldt area at these depths.

4.1 Wind stress and coastal upwelling connections

The seasonal cycles of coastal upwelling for both the Northern and Southern Humboldt areas are consistent with the seasonal patterns of wind stress reported in previous studies (Wang et al. 2010; Rykaczewski et al. 2015) in the 20th and 21st Centuries. Also, outcomes from modelling approaches and observational data in previous studies strongly suggest that wind stress has been the main driver of coastal upwelling during the 20th Century in the Humboldt system. To extend this analysis, we use the area average of simulated monthly means of the upward ocean mass transport variable contained in the CMIP5 data set (Table 1) to explore the connection of modelled wind stress and coastal upwelling activity at different depths during the 20th and 21st Centuries.

In the Northern Humboldt area, wind stress is projected to increase by 5.8% by 2071–2100 compared to 1971–2000 under the RCP8.5 scenario (Table 2). However, coastal upwelling is projected to rise at surface layers (i.e., 30 and 50 m depth), but to decline below about 100 m depth. This is particularly relevant looking at the consistent connection at these depths between simulated wind stress and coastal upwelling during the 20th Century. Despite the almost null (−0.5%) anomalous wind stress computed for 20th–19th C the computed anomalies of the upward water flux at 100, 150 and 300 m depth (−2.9, −5.4 and −30.2% respectively) are substantially less significant than the computed anomalies for 21st–20th C in the same area at the same depths (−5.8, −15.1 and −104.5% respectively). This, in spite of the projected increase of wind stress for 21st–20th C (5.8%) in the Northern Humboldt area. This pattern is similar for the Southern Humboldt area. The increase of wind stress in 20th–19th C was only 0.6%, yet was associated with negative anomalous upwelling. These changes are not statistically significant, so likely represent internal variability. Meanwhile the projected percentage change for 21st–20th C at 100, 150 and 300 m depth has some significance (Table 2), despite the noticeable increase of wind stress computed for 21st–20th C (19.6%).

This apparent disconnection between wind stress and coastal upwelling below 100 m depth is spatially presented for the spring season in Fig. 4. Whereas coarser models clearly represent a decline of upwelling covering most of the both the Northern and Southern Humboldt areas, finer models depicts a more complex pattern suggesting a less pronounced decline (Fig. 4). However, despite the differences given by models resolution, the spatial distribution of the climate change signal at 100 m depth suggest a clear disconnection from the increasing wind stress projected by 2071–2100.

Although wind stress is expected to intensify in both Humboldt areas, as Bakun (1990) suggests, a more pronounced increase is projected for the Southern Humboldt (Fig. 5a, b). This increasing wind stress trend during the 21st Century, with a more pronounced rise at higher latitudes than nearer the equator, was also reported by Rykaczewski et al. (2015). As Belmadani et al. (2013) argued, the increasing future trend of wind stress alongshore in the Southern Humboldt latitudes might be related to a strengthening and poleward displacement of the South Pacific anticyclone, which in turn would lead to an increase of the gradient of pressure alongshore (Bakun 1990). However, CMIP5 suggests that under substantial climate change, the increasing wind stress would lead to an increased upwelling at surface layers, but not at deeper depths. In the Northern Humboldt this disconnect becomes noticeable in the 1930s (Fig. 5c). From 100 to 300 m depth the vertical flux has gradually declined, and under the RCP8.5 this disconnection would increase significantly from 2020. Only at 30 m depth the vertical water mass transport would increase from about the year 2045 (Fig. 5c): coherent with the projected wind stress trend (Fig. 5a). In the Southern Humboldt, similar general trends are observed. However, the disconnection between wind stress and the intensity of coastal upwelling is appreciable from the 1970s, becoming more evident from the 1990s (Fig. 5d). In this region, the intensity of upwelling increases from 50 m depth above the surface from around 2025, which coincides nicely with the wind stress pattern computed for that region over the RCP8.5 scenario (Fig. 5b). Despite the disconnect between wind stress and coastal upwelling with depth, it is noticeable that wind stress drives the decadal variability even at 300 m depth and will continue setting it throughout the 21st Century.

Table 2 Coastal upwelling and wind stress percentage change in spring season for the Northern (NH) and Southern (SH) Humboldt areas

Fig. 4

Spatial distribution along the Humboldt system of the climate change signal computed for 2071–2100 with respect to 1971–2000 from CMIP5 data set in the spring season (September–October–November). Black squares in each panel are the Northern Humboldt (top) and Southern Humboldt (bottom) areas. Data is vertically interpolated to 100 m depth in native units (${\hbox {kg}}\,{\hbox {s}}^{-1}$)

Fig. 5

Ensemble mean time series of anomalies for wind stress (a, b) and coastal upwelling (c, d) at 30, 50, 100, 150 and 300 m depth for the Northern (a–c) and the Southern (b–d) Humboldt areas in spring season (September–October–November). Anomalies computed with respect to 1871–1900. Time series are presented for the spring season as this is an active coastal upwelling period in both regions along the whole Humboldt system

4.2 Stratification and thermocline depth implications

Our results indicate that wind stress is not the sole driver of change in coastal upwelling. We propose that ocean stratification due to global warming (Fig. 6) plays a significant role. The warming of the ocean has long been linked to greater stratification of the thermocline (Roemmich and McGowan 1995) and connections of stratification with upwelling have also been outlined (Chhak and Lorenzo 2007). In both the Northern and the Southern Humboldt areas, stratification would increase significantly under the RCP8.5 scenario as can be seen looking at the potential temperature trends (Fig. 6). Greater warming occurs nearer the surface of the water column during the 21st Century, which results in an increased ocean stratification. In the Northern region the projected increase of temperature by 2100 at surface sea water layers is slightly higher than in the Southern region, and ocean stratification started to progressively raise from the 1930 in both areas. As ocean stratification increases, upwelling of coastal sea waters sourced from 100 m to at least 300 m depth becomes less connected to the wind stress. This is being amplified as global warming increases. As the Bakun (1990) hypothesis proposes, a strengthening of upwelling-favorable winds is projected along the Humboldt system. However, even though the increasing wind stress under global warming would drive stronger upwelling at the surface (namely 30–50 m depth), its net effect on deeper ocean layers (100–300 m depth) would gradually decline because of the increasing stratification. It is noticeable that, unlike the wind stress, the change in the stratification is a gradual process which does not exhibit decadal variability in simulations. This may explain the predominant control of wind stress on the decadal variability of upwelling along the water column (Fig. 5). It might be possible that the increased stratification alters the expression of the wind stress-driven decadal variability. However, there are insufficient model realizations to test that hypothesis.

The effect of stratification is also appreciable from the analysis of the Ekman layer depths computed from a theoretical formulation (Eq. 3) and its correspondence with the ocean mixed layer thickness obtained from CMIP5 models for the 20th and 21st Centuries (Fig. 7). Ekman layer depth and the ocean mixed layer thickness start to diverge from the last decade of the 20th Century. The Ekman layer deepens slightly, as expected as wind stresses increase, and the ocean mixed layer becomes shallower. The surface mixed layer in the ocean is directly influenced by the wind and air-sea fluxes such as heat and freshwater exchange. Although the Ekman transport, mainly driven by wind stress, is not strong enough to influence the mixed layer by turbulence (Talley et al. 2011), these two apparently independent dynamical phenomena are simulated to have varied with some correspondence for most of the 20th Century under wind stress forcing. However from the end of the 20th Century, the gradually increasing stability of the upper ocean (because of the increased stratification) leads to a shallower surface mixed layer. This event is directly related with an increase of upwelling at shallow levels of the coastal sea water, where wind stress directly exert influence on the Ekman transport and layer depth, but a progressively decline of coastal upwelling intensity is projected at the thermocline where stratification significantly will increase (Fig. 10).

Ekman layer depth depends on wind stress, sea water density and the Coriolis force, with a constant of proportionality $\gamma$ (Eq. 3). However, accepted values in literature for parameter $\gamma$ come from observations, and are presumed constant in time. Our results suggest that stratification impacts the values of $\gamma$, and that it should vary with climate.

Fig. 6

Ocean Stratification. Ensemble mean time series of potential temperature anomalies at 30, 50, 100, 150 and 300 m depth for the Northern (a) and Southern (b) Humboldt areas in spring season (September–October–November). Anomalies computed with respect to 1871–1900

Fig. 7

Ensemble mean time series of the change of the ocean mixed layer thickness (blue), Ekman layer depth (black), and anomalous wind stress (green) with respect to 1871–1900 for the Northern (a) and Southern (b) Humboldt areas. Light blue and grey lines correspond to the ocean mixed layer and Ekman layer depth of each model, respectively

5 Model resolutions and limitations

It is clear from Fig. 1 that the representation of upwelling depends in part on model resolution. Finer models (Fig. 8; namely CNRN-CM5, GFDL-CM3, GFDL-ESM2M, MPI-ESM-MR, MRI-CGCM3 and MRI-ESM1) show a similar qualitative response at all depths in comparison to the full ensemble mean (Fig. 5). It is noticeable that the finer models ensemble mean simulate a higher magnitude of reduction in upwelling intensity, although this is probably a function of their larger baseline upwelling (Fig. 2). In accordance with the computed coastal upwelling trends (Figs. 5, 8) below 100 m depth most of the models project a pronounced upwelling decline in both the Northern (Fig. 9a) and Southern Humboldt areas (Fig. 9b). Although finer models tend to project slightly higher changes regarding absolute upwelling in the Southern region, no particular differences are visible between finer and coarser models regarding percentage change. The decrease in upwards mass transport below around 100 m is a robust result that is not dependent on the spatial resolution (Fig. 9). The precise depth at which the transition from more to less upwelling occurs is certainly model dependent, but resolution alone does not appear to be its driver.

Systematic SST biases have been seen in major coastal upwelling systems (Yang et al. 2015), with increasing resolution being one strategy to counter it (Gent et al. 2011). Whilst there clearly are SST biases in the simulation of the Humboldt Current, they are not predominantly associated with ocean model resolution (not shown) and are uncorrelated with the future response (also not shown). Nonetheless, we cannot preclude a systematic bias impacted the whole ensemble’s representation of upwelling. Regional ocean models can be used to overcome ocean-resolution dependencies (e.g. Oerder et al. 2015), unfortunately these cannot address issues associated with poor simulation of the overlying atmospheric conditions (Gent et al. 2011).

Fig. 8

Ensemble mean time series from finer models (Table 1) showing changes in wind stress (a, b) and coastal upwelling (c, d) at 30, 50, 100, 150 and 300 m depth for the Northern (a–c) and the Southern (b–d) Humboldt areas in spring season (September–October–November). Anomalies computed with respect to 1871–1900

Fig. 9

Vertical profiles of anomalous coastal upwelling in spring (September–October–November) computed for 2071–2100 with respect to 1971–2000. Black lines show the ensemble mean of the upward ocean mass transport change. Individual models are coloured by their resolution, with red indicating finer and grey coarser models (Table 1)

6 Considerations for marine life and fishery productivity

The Humboldt system sustains the largest yield of small pelagic fish on the planet. This marine ecosystem is highly susceptible to changes in upwelling patterns, such as alteration in frequency or magnitude (Aravena et al. 2014). Significant potential impacts due to climate change have been projected on the spawning rate of small pelagic fish, such as the Peruvian anchovy Engraulis ringens and the sardine Sardinops sagax (Brochier et al. 2013), despite the projected temperatures remaining within accepted ecological ranges (Bertrand et al. 2004). Small pelagic species strongly depend on zooplankton availability and the exceptional phytoplankton productivity in coastal upwelling zones (Ayón et al. 2011). Therefore the ecological impacts of a decline of upwelling below 100 m depth are not trivial to understand and project.

Our results, based on the upward ocean mass transport, suggest active upwelling until at least 300 m depth throughout the 20th and 21st Centuries. The decreasing intensity under 100 m depth is focused in coastal waters (Fig. 4), where daily zooplankton dynamics is active and fundamental for marine life. The macro-zooplankton is a fundamental dietary component for fish species, such as the Peruvian anchovy Engraulis ringens. The Humboldt System sustains a significant biomass of macro-zooplankton which usually migrates within a day from the surface to about 300 m depth playing an essential ecological role (Ballón et al. 2011). This is highly sensitive to climate oscillations and strongly depends on upwelled water (Mackas et al. 2006). Whereas in the Benguela system an increasing zooplankton abundance has been observed recently, for the California and Humboldt systems decreasing trends have been reported (Perry et al. 2004).

A rise of nutrient abundance has commonly been associated with increased upwelling by increasing wind stress, yet increasing stratification might lead to a decline in nutrients supply (Bakun et al. 2015; Oerder et al. 2015). Observational data have shown a reduction in zooplankton abundance under stronger stratification (Roemmich and McGowan 1995). Whereas some modelling approaches argue a potential counteracting effect of zooplankton abundance by increasing wind stress and upwelling (Auad et al. 2006), our results suggest that stratification is already playing a significant role in a decreasing trend of subsurface upwelling intensity. As Bakun et al. (2015) state, the relationship between the upwelling behaviour throughout the water column, the role of stratification and wind stress intensification is still an open question. Although our results here do not refute the Bakun (1990) hypothesis, they expose the limit of its applicability. Furthermore, phytoplankton production might decrease in coastal upwelling ecosystems due to deeper and stronger surface mixing and light limitation associated with enhanced wind stress (Bakun et al. 2015). In effect, the optimal environmental range for primary productivity in coastal marine ecosystems has been associated with moderate upwelling intensity (Cury and Roy 1989; Bakun et al. 2015). Outcomes from a coupled physical-biochemical model by Brochier et al. (2013) project decreasing nursery areas for small pelagic fish under extreme climate change in the Humboldt ecosystem.

In the Humboldt system, nutrients are distributed throughout the water column, as shown by observational data obtained from oceanographic expeditions (Silva et al. 2009). A significant concentration of nitrate, associated with phytoplankton primary production (Fréon et al. 2009), and phosphate were measured under 200–1400 m depth. Whereas adverse impacts on nutrient abundance from a decline of coastal upwelling below 100 m depth seem plausible, additional research is required to confirm them. Reproduction success of small pelagic fish and changes under climate change also depend on additional factors not commonly considered in model formulations, such as ocean acidification (Bakun et al. 2015), ecosystem structure (Rykaczewski and Checkley 2008) or external dynamics associated to environmental factors and predators dynamics (Brochier et al. 2013). As Travers et al. (2007) suggest, these interactions should be addressed with marine ecosystem models comprising physical-biological integration and two-way interactions for the ecosystems components.

In addition to consequences on marine life and fishery at local or regional scales, changes in coastal upwelling patterns could lead to effects at global scale. Global food supply partially relies on fishery productivity in the major coastal upwelling systems (Fréon et al. 2008; Brochier et al. 2013). As Bakun et al. (2015) pointed out, coastal upwelling systems are resilient and robust regarding productivity under natural variability. However, the ecological impacts under enhanced anthropogenic climate change are still uncertain.

7 Conclusion

The Humboldt system is the most productive system in the world, and coastal upwelling plays a significant role in marine life and fishing along the western coast of South America. Therefore, climate impacts on this system could result in ecological and socioeconomic consequences not only locally but in the global food supply. As wind stress is a strong driver of coastal upwelling, most of the studies so far project its future behavior during the 21st Century from wind data, observed or modelled. We analysed coastal upwelling and related variables from 13 AOGCMs contributing to CMIP5. Despite the differences in spatial resolution, all the AOGCMs represent coastal upwelling, although finer models tend to have their coastal upwelling confined to the coast. The vertical profile of upwelling is consistent across the models, and no significant differences between finer and coarser AOGCMs are appreciable for either vertical distribution or the magnitude of the upward water mass flux.

Coastal upwelling in the Humboldt system exhibits distinctive seasonal cycles for the Northern and Southern regions. Whereas in the Northern Humboldt upwelling is present most of the year except summer, the upwelling season spans spring and summertime in the Southern Humboldt. This pattern is consistent with the seasonality of the upwelling-favorable winds previously described for both areas. According to this, we analysed future projections of coastal upwelling and its connection with wind stress during the spring season, which is an active upwelling period in the whole Humboldt system. Our results suggest that wind stress has not been the sole driver of change in coastal upwelling during the 20th Century, and ocean stratification has also exerted a significant effect on the upward water mass transport below 100 m depth throughout the second half of the past Century.

Wind stress is projected to increase during the 21st Century over the Humboldt Current, though more pronounced in the South. The wind stress increase leads to enhanced coastal upwelling at surface levels, namely 30-50 m depth. However, from the late 20th Century a disconnect between coastal upwelling from wind stress has taken place below 100 m depth, in spite of wind stress still controlling decadal variability at these depths during the 21st Century. This result is not dependent on model resolution. Ocean stratification is already partially obstructing the upward ocean mass transport at sea water layers under 100 m depth, and this trend will further accentuate as global warming increases. Additionally, increasing ocean stratification is leading to a shallower thermocline and hence a decreasing ocean mixed layer thickness (Fig. 10). These results strongly suggest that looking at wind stress and the Ekman depth dynamics are insufficient to fully describe changes in coastal upwelling for the 21st Century.

Fig. 10

Schematic representation of our suggested impact of climate change on coastal upwelling in the Humboldt upwelling system. During most of the 20th Century (left) wind stress has led to coastal upwelling activity between the surface and 300 m depth. However increasing ocean stratification has led to a shallower thermocline, and a disconnection of wind stress and coastal upwelling below 100 m depth has gradually taken place since the second half of the last Century. Although during the 21st Century (right) wind stress is projected to increase leading to an enhanced Ekman transport off-shore and an intensification of coastal upwelling at the surface (30–50 m depth), the effect of ocean stratification will progressively be strengthened as ocean temperatures increase. The potential ecological and socioeconomic effects of an eventual declining of coastal upwelling intensity below 100 m depth remain uncertain