The Early Life History Transitions of the Bivalve Aulacomya atra From the Humboldt Current System Off Peru Are Affected by Human Exploitation and Modulated by El Niño–La Niña Cycle

Introduction

The Humboldt Current upwelling system off Peru (hereafter HCS) is considered the most productive marine ecosystem of the world (Ryther, 1969). Among the four recognized upwelling centers in the HCS, the permanent upwelling area of Bahía Independencia is the most important, with extremely high biomass and productivity of shellfisheries (Mendo and Wolff, 2003; Riascos et al., 2016). A strong, coastal upwelling driven by trade winds brings deep, nutrient-rich water into the euphotic layer, fueling vast primary production that is then converted in high secondary productivity (Arntz and Fahrbach, 1996). However, environmental conditions in the HCS are strongly modified by El Niño-Southern Oscillation (ENSO). During El Niño, the warm phase of ENSO, the arrival and poleward propagation of coastal Kelvin waves along the Peruvian coast lead to strong abiotic changes, including the anomalous increase of sea surface temperature (SST) and a drastic deepening of the pycnocline and nutricline, which reduce the efficiency of the upwelling process (Barber and Chavez, 1983; Arntz et al., 2006). These changes are linked to dramatic and widespread biological effects, operating at different spatial scales and levels of organization (Arntz and Fahrbach, 1996; Tarazona and Arntz, 2001; Peña et al., 2005; Riascos et al., 2012).

Despite decades of effort, understanding the links between climate variability and population dynamics of exploited marine species remains elusive, largely because natural ecosystems operate as complex adaptive systems and demographic parameters are influenced by large-scale oceanographic processes (Levin, 1998; Bakun, 2010). Therefore, classical experimental approaches are generally not feasible to test environmental drivers of population dynamics at the appropriate scale. However, the ENSO cycle has been recognized as a natural experiment to study biotic responses at ecosystem and inter-annual scales (Trillmich, 1991; Arntz and Fahrbach, 1996). Moreover, the population dynamics of most coastal species in the HCS is controlled by mesoscale oceanographic processes, including upwelling events, frontal systems and eddies that influence variations in recruit supply and upwelling shadows in bay systems impacting larval retention (Thiel et al., 2007).

The reproduction, and the subsequent population processes leading to recruitment, are key in controlling long-term variability and persistence of marine bivalve populations (Bricelj et al., 1987; Lima et al., 2000; Beukema and Dekker, 2007). When oceanographic processes favor the local retention of early life stages, they represent the most important driving force in shaping year-to-year population variation (Sinclair, 2009; Bakun, 2010).

The oceanographic conditions at Bahía Independencia are broadly characterized by low surface temperatures and high nutrient levels due to permanent upwelling. The permanent circulation of water masses modulated by wind patterns determines the presence of upwelling shadow fronts, as observed in similar embayments along the HCS (Marín et al., 2003; Quispe et al., 2010). These mesoscale processes may strongly regulate larval movement in this nearshore retention area and hence the recruitment dynamics and other demographic features of benthic organisms (Escribano et al., 2003; Ashton et al., 2008; Riascos et al., 2011). Therefore, these are ideal systems to understand how climate variability affects the critical stage of early life stages that may control population dynamics. The ribbed mussel Aulacomya atra (Molina, 1782) is among the most important bivalve species for shellfisheries in Bahía Independencia. A. atra is distributed along the Pacific Ocean from El Callao (Peru) to Estrecho de Magallanes (Chile), along the Atlantic Ocean (from southern Argentina to southern Brazil) and on the South African coast (Castilla and Guiñez, 2000).

In this context, our work aims to evaluate the coupling between the production and release of gametes, the subsequent larval settling and survival until recruitment of the commercially- important bivalve A. atra in Bahía Independencia throughout a 9-year period that included major EN and LN episodes. Moreover, we wanted to evaluate the oceanographic drivers of changes of demographic processes of A. atra, from the production of gametes to the recruitment of individuals to the adult population.

Materials and Methods

Study Area and Sampling

Monthly samples were taken between January 1996 and September 2004 on the rocky subtidal area of Santa Rosa Island (−14.322167; −76.156833) in Bahía Independencia, part of the Paracas National Reserve (Peru; Figure 1). Six replicate samples were taken each month by divers using 25 × 25 cm quadrats randomly placed on the rocky substrate at about 15 m depth, all specimens of A. atra were collected and placed in plastic bags. In the laboratory, the resulting animals were counted, the shell length (SL, the anterior-posterior shell length) was measured with a Vernier caliper (precision: ±0.05 mm) and the individual wet mass was registered using a Mettler scale (precision: ±0.01 gr) after cleaning fouling organisms growing on the shell. Soft body parts were separated, weighed and oven-dried to constant mass and finally ash-free dry mass (AFDM) was obtained by ignition at 600°C for 3 h.

FIGURE 1

Figure 1. Map of the study area (Isla Santa Rosa-Bahía Independencia), with the sampling location for adult Aulacomya atra (black circle) and the location of the system to monitor larval settlement (red point) at La Vieja Island.

Spawning Stock Biomass

Changes in the spawning stock biomass (SSB) (g AFDM m^–2), defined as the combined weight of all individuals in a population that are capable of reproducing, was estimated as a proxy for changes in the capacity for propagule production using a two-step approach. First, the monthly length-mass relationship was estimated from monthly shell lengths and their corresponding body mass. For this, a linear regression analysis on log-transformed data was used to estimate the constants a and b of the relationship:

$$\text{AFDM=}a{\text{SL}}^b$$

where AFDM is the ash-free dry mass (body mass) and SL is the shell length. Thereafter, the SSB was calculated as the product of the mean abundance of sexually mature animals and their corresponding body mass (Ec. 1). Mussels ≥ 30 mm in shell length were considered sexually mature (∼25% of the maximum length present in pooled length-frequency data) (see Eversole, 1989).

Larval Settlement Monitoring

To monitor monthly changes in larval settlement of A. atra, a collecting system was developed at La Vieja Island, near the sampling area (Figure 1) to avoid interferences with local small-scale fishing activities. The system used for larval settlement consisted of six replicated units of 65 g of filamentous cultch inside a 0.015 m³ metal frame box covered with fishing net, arranged in two rows oriented perpendicularly to the coast (see Figure 1). All cultch were extracted monthly, stored in plastic bags and replaced with a new cultch array by scuba divers. Samples were preserved using 10% formalin solution and transported to the laboratory. Thereafter, cultch were gently rinsed with filtered seawater over 250 and 710 μm-mesh sieves to retain a biological fraction that included 1-month old A. atra spat (Ramorino and Campos, 1983). This material was examined under a stereomicroscope and the diagnosis criteria reported by Ramorino and Campos (1983) and Chanley (1968) were used to identify A. atra. Finally, the monthly mean number of postlarvae per cultch was calculated and used as an indicator of changes in larval settlement.

Recruitment

We first defined “recruits” of A. atra as animals with a SL ≥ 10.9 mm, estimated from the mean of the smallest animals in each month during the study period. Thereafter, the density of recruits (ind. m^–2 month^–1), was obtained by the “recruitment” function of the TropFishR package (Mildenberger et al., 2017, 2018). This function calculates recruitment patterns of a stock by means of a hindcast projection of length frequency data using the seasonalized von Bertalanffy growth curve (VBGF; Pauly, 1982). For each length frequency bin, the time at which the growth curve equals a length of zero is recorded. The overall recruitment pattern simply reflects the frequency of crossings as weighted by each bin’s catch frequency. The method assumes that (i) all individuals belonging to a data set grow as described by a single set of growth parameters and (ii) the output of this function is simply count data and corresponding densities (ind m^–2 month^–1). Therefore, it is possible to have densities for each month (Gayanilo et al., 2005). Previously, the growth parameters (L∞, K) were calculated using the routine ELEFAN_SA in TropFishR (see Riascos et al., 2017) and the estimate of “t0” was obtained by the Pauly equation (Pauly, 1983):

$${\text{log}}_{10}\left(-\mathrm{t0}\right)=-0.3922-0.2752{\text{log}}_{10}L\mathrm{\infty }-1.038{\text{log}}_{10}K$$

where L∞ is the asymptotic length (cm), K is the growth coefficient (year^–1) and t0 is the theoretical time zero, when growth curve crosses length equaling zero.

Environmental Parameters

During the El Niño (EN)–La Niña (LN) cycle, sea life in Peruvian coastal waters may be affected by both local and remote factors, through teleconnections (Riascos et al., 2017). Therefore, we included variables describing environmental variability at macro, regional, and local spatial scales. The multivariate ENSO index (MEI), which combines multiple meteorological and oceanographic factors in a single index, was used as a descriptor of ocean-atmosphere weather variability in the tropical Pacific¹. The monthly-averaged SST at the region “Niño 1 + 2” (0°–10°S, 90°W–80°W) (SST Niño 1 + 2) was also used as a relevant variable describing large-scale thermal variability (data provided by the Climate Prediction Center²). The El Niño-Coastal Index (ENCI), an index based on SST anomalies in the El Niño 1 + 2 Region (see Takahashi et al., 2014), was taken as a descriptor of interannual anomalies of SST at a regional scale (i.e., Peruvian coastal waters). Additionally, daily measurements of SST taken at a fixed station (−13.7075; −76.2208) were provided by the Instituto del Mar del Peru (impSST); sea surface salinity (Sal) in the Pisco Region was taken from Morón (2011). Finally, in situ measurements of SST, sea bottom temperature (SBT), dissolved oxygen in the surface (SDO), dissolved oxygen in the bottom (BDO) and chlorophyll-a concentration were taken monthly using a SBE 19 CTD profiler in two areas: Santa Rosa Island and La Vieja Island (Figure 1).

Data Analysis

A first assessment of variability of the SSB, the number of settled postlarvae and recruits showed that two of these variables had a significant decreasing trend over time (Figure 2), which would obscure the stationary properties in our data series. Therefore, these series were first detrended, assuming a linear model as a reasonable model to estimate the trend. To obtain the detrended series we subtracted the linear trend component as x = x − (a + bt) and worked thereafter with the residuals of SSB and density of recruits (Figure 3). A generalized linear model was used to test the influence of SSB (independent variable) and period (LN, EN, Neutral; covariant categorical factor) on the subsequent larval settlement (Log-transformed number of settled postlarvae per cultch; dependent variable) as follows:

$$y_{ij}=\mathrm{\mu }+{\mathrm{\alpha }}_i+\mathrm{\beta }{x}_{ij}+e_{ij}$$

FIGURE 2

Figure 2. Temporal variability of (A) the spawning stock biomass and the larval settlement and (B) the recruitment of A. atra. Data on settled postlarvae are log-transformed (Log X + 1) to avoid the overweight of a few high values and allow an easier view of the variability in the dataset.

FIGURE 3

Figure 3. Residual variability of (A) the spawning stock biomass and the larval settlement and (B) the recruitment of A. atra. Red and blue background denote the timing of EN and LN episodes, respectively, according to the El Niño-Coastal Index.

where y_ij is the jth value in the group i for the number of postlarvae per cultch, μ the grand mean, α_i the effect of the ith group (SSB), x_ij is the jth value in the group i for the covariate (period), i = 1,2,3 are the groups, β the slope for the covariate and e_ij the associated error term for the jth observation in the ith group. The model used the Poisson distribution for the response variable and log as the link function.

A similar model was used to test the relationship between SSB as the independent variable, the period as a covariant and the density of recruits as the dependent variable. Based on these results, least-square fitting procedures followed by F tests were used to assess the significance of the relationships (i.e., describe the coupling) between SSB versus larval settlement and SSB versus density of recruits under EN, LN, and neutral conditions.

A two-step approach was used to assess the oceanographic factors affecting the demographic processes from the production of gametes to the recruitment to the adult population of A. atra. First, although the various indices and environmental parameters included in this study describe variability at different spatial scales (i.e., in situ, local, regional), it is expected that they concurrently- affect demographic processes of A. atra and could therefore be highly- correlated. This interdependence among environmental factors or multicollinearity, creates problems with the interpretation of a least-squares multiple regression analysis. To address this problem, the two primary data sets on environmental variability at Santa Rosa Island (local and remote environmental parameters) and La Vieja Island (local and remote environmental parameters) were first analyzed by Principal Components Analysis (PCA). Data were normalized (subtracting the mean and dividing by the standard deviation for each variable) to account for the different scales among variables before calculating Euclidean distances. Moreover, Draftsman plots (see Supplementary Figures S1, S2) were used to assess multicollinearity and check that data are roughly- symmetrically distributed across the range of each variable, thus achieving approximate (multivariate) normality underlying PCA analysis (Clarke and Gorley, 2006). This procedure reduces the volume of data to a handful of orthogonal (i.e., mutually independent) principal components that account for a proportion of the variation in the environmental factors. Vector plots were superimposed on the bi-dimensional PCA plots and the length and direction of vectors were evaluated to assess the contribution of each environmental variable to each PC. In a second step, we assessed the relationship between SSB, larval settlement and recruitment of A. atra and the four most important principal components (PCA scores), using a forward stepwise multiple regression analysis. Compared to traditional approaches looking for correlations with multiple environmental variables but not considering collinearity, this approach that includes multicollinearity among independent variables is more robust (Myers, 1986).

Results

Early Stage Transitions During the ENSO Cycle

The temporal variability depicted in Figure 3 shows that the SSB of A. atra follows a relatively consistent seasonal variability, with peaks during June to August and November to January and decreases between April to May and September to November. However, seasonal changes were not clear for the larval settlement and subsequent recruitment of A. atra, hinting at the influence of oceanographic conditions, interspecific interactions and stochastic processes on e.g., the settlement success or mortality. The results of the generalized linear model (Table 1A) showed that the SSB was a significant predictor of the larval settlement of A. atra, whereas the type of period (EN, LN, neutral) (Table 1A) was not significant as a covariant for this relationship. Figure 4A depicts how the number of settled postlarvae varied as a function of the SSB: the relationship was best described by a polynomial model, where the larval settlement decreased from low to average values of the SSB and steeply increased with increases in the SSB.

TABLE 1

Table 1. Results of the generalized linear model testing (A) the influence of the spawning stock biomass and the type of period (EN, LN, neutral) on the larval settlement of A. atra, and (B) the influence of the spawning stock biomass and the type of period (EN, LN, neutral) on the recruitment of A. atra.

FIGURE 4

Figure 4. (A) Assessment of the relationship between the spawning stock biomass (SSB) versus the log-transformed number of settled postlarvae of A. atra; broken lines represent the 95% confidence interval for the regression line. (B) Assessment of the relationship between the SSB and the density of recruits during EN, LN and neutral conditions. EN is denoted by red dots, LN by blue dots and neutral conditions by black dots. Red and gray shadows represent the 95% confidence intervals for the regression lines.

In contrast, the SSB and the type of period were significant predictors of the density of recruits (Table 1B). Figure 3B shows strong fluctuation of the recruitment during EN and LN, with higher recruitment during the transition to, and onset of, EN. In contrast, recruitment was consistently low and stable during neutral conditions. Figure 4B shows that the recruitment of A. atra was generally lower during LN than during EN and neutral periods. The relationship between the SSB and the recruitment density of A. atra was significant during EN and neutral conditions, but not significant during LN. Again, a polynomial model best described these relationships, indicating that lower levels of recruitment occurred in average values of the SSB.

Environmental Drivers of Population Processes

The principal component analysis (Figure 5) showed that PC1 explained a high proportion of environmental variability in La Vieja Island (45.9%) and Santa Rosa Island (47.6%). In both cases, EN-months (red squares) were clearly separated from LN-months (blue squares) through PC1, hinting at the strong influence of ENSO on abiotic conditions in both places. Moreover, it was clear that while the 1997–1998 EN episode represented extremely different physical-chemical conditions as compared with neutral periods, oceanographic conditions during the 1996–1997 and 1999–2001 LN episodes were not so clearly different from neutral periods. The length and direction of vectors of environmental factors (Figure 5 and Table 2) show that thermal conditions (in situ sea surface and bottom temperature, and ENCI) had the strongest effect on the differences depicted in PC1. In turn, the PC2 explained 26.4 and 20.4% of the variability in La Vieja Island and Santa Rosa Island, respectively, with an influence of biochemical factors (i.e., chlorophyll-a concentration, salinity and dissolved oxygen) better explaining this variability. The stepwise regression analysis (Table 3) showed that a regression model including one or two PC scores explained almost half of the total variability of the SSB and the density of recruits, and one-third of the variability of the number of settled postlarvae of A. atra. The high power of PC1 scores to predict changes in the SSB and the density of recruits suggest that shifts in the thermal regimes associated with ENSO variability are critical factors for the changes in the body condition, the abundance of adults and the survival of recruits in the benthic zone. In contrast, the fact that PC2 scores had more explanatory power to predict the number of settled postlarvae suggests that biochemical and food conditions are key factors affecting the larval settlement of A. atra.

FIGURE 5

Figure 5. Biplots of Principal Components (PC1 and PC2) after the PCA analysis of regional and local environmental factors at La Vieja Island and Santa Rosa Island and the superimposed vectors (gray lines) of environmental factors (Chlor-a: chlorophyll-a concentration; impSST: sea surface temperature reported by IMARPE; SST Niño 1 + 2: sea surface temperature in the El Niño 1 + 2 Region; SST: in situ sea surface temperature; SBT: sea bottom temperature; ENCI: El Niño-Coastal Index; MEI: Multivariate ENSO Index; SDO: in situ dissolved oxygen; BDO: in situ dissolved oxygen bottom; Sal: salinity). Red and blue text/squares represent periods of EN and LN conditions, respectively.

TABLE 2

Table 2. Results of the principal components analysis of regional and local environmental variability at La Vieja Island and Santa Rosa Island.

TABLE 3

Table 3. Results of the stepwise regression analysis testing the influence of environmental variability (reduced to five PCA scores) on (A) the spawning stock biomass, (B) the number of settled postlarvae and (C) the density of recruits of Aulacomya atra.

Discussion

Human Exploitation and ENSO Influence the Early Stages of Aulacomya atra

The persistent decreasing trend observed for the SSB, the larval settlement and the recruitment of A. atra was an obvious feature of the variability of these population parameters. A significant decreasing trend of the mean abundance (see Supplementary Figure S3) of this heavily exploited bivalve seems to be related to this result. A strong drop in the mean abundance, and the fact that fishermen commonly select larger animals, implies a drop in the SSB, less larvae in the water column and possibly fewer individuals surviving into recruits. This persistent decreasing trend in the stock is in line with an independent study that monitored benthic communities from Bahía Independencia between 2004 and 2010 (Riascos et al., 2016). They found that the continuous exploitation of four species (including A. atra) was the most important driving force of changes in the benthic community structure, characterized by a ∼43% reduction in biomass of the whole benthic community. Moreover, as the most important species for shellfisheries in the study area [the Peruvian scallop Argopecten purpuratus (Lamarck, 1819)] has been heavily exploited, its fishery was regulated, resulting in an increased fishing pressure on A. atra after 1996 (Mendo, 1997). Therefore, the ecological effects of human exploitation in Bahía Independencia are seemingly overlaying the expected responses of A. atra populations to changes in climate patterns. This underscores the need to integrate the interaction between humans, biological communities and the particularities of the physical setting in ecological studies (Riascos et al., 2016).

Seasonal changes in the frequency and magnitude of upwelling pulses – with the concurrent changes in temperature and food availability – have been suggested to drive the chronology of the reproductive cycle of A. atra at several locations in Chile (Jaramillo and Navarro, 1995; Avendaño and Cantillánez, 2013, 2014). In our case, strong decreases in the SSB, which most likely reflect spawning episodes (Riascos et al., 2017), coincide with periods of drops in the SST in Bahía Independencia (e.g., Vélez et al., 2005). This suggests that the reductions in the SSB, reflecting periods of gamete release, are influenced by seasonal pulses of cold water. These are likely related to increased upwelling strength in a region characterized by year-round intense upwelling conditions (Calienes et al., 1985).

The described relationship between the larval settlement and the SSB most likely reflects the time lag between drops in the SSB, indicative of gamete release, and peaks in the larval settlement 2–3 months later, when the SSB has recovered. However, there is considerable variability not explained by the modeled relationship between the larval settlement and the SSB, hinting at the complex and interacting oceanographic, biological and stochastic processes affecting the temporal coupling between these parameters.

El Niño-Southern Oscillation significantly changes abiotic conditions in the northern HCS, with clear consequences at the population organization level for coastal bivalves that range from transient changes in body conditions to local extinctions (e.g., Carstensen et al., 2010; Riascos et al., 2011; Uribe et al., 2012). Our results indicate that EN and LN episodes did not affect the relationship between SSB and the larval settlement, but they did affect the relationship between SSB and the recruitment. One could reasonably argue that the transition of gamete release to larval settlement is mostly determined by endogenous, physiological factors driving the spawning pulses, whereas the transition of gamete release to recruitment implies longer exposure to oceanographic and biotic variability associated with EN–LN episodes. However, we must also be aware that environmental anomalies associated with EN–LN and biotic interactions operating in the water column may be completely different from those operating on soft bottoms. This difference seems to explain why pulses in the density of recruitment were observed despite the lack of settled postlarvae on artificial collectors suspended in the water column (Figure 2).

Diverse studies have shown that coastal species in the HCS react in different – often opposing – ways to the oceanographic changes associated with EN and LN episodes, depending on their biogeographic origin (Riascos, 2009; Carstensen et al., 2010; Ortega et al., 2012; Uribe et al., 2012). In this context, it seems unsurprising that while the lowest values of SSB of A. atra were observed during the strong EN episode in 1997–1998, the highest values were observed during the previous and subsequent LN episodes. Most likely, this result has to do with its “biogeographic affinity,” i.e., the similarities in the flora and fauna of one location compared to those of another location, resulting from their common evolutionary history and/or their contemporary connectivity (sensu Loveland and Merchant, 2004). A. atra currently display a disjointed geographic distribution in the Southern Hemisphere (southern coast of South America and South Africa), most likely resulting from either vicariance or larval connectivity, and the fossil record strongly suggests that the genus Aulacomya originated in the Antarctic region (Castilla and Guiñez, 2000). Therefore, it is not surprising that the body condition and abundance of A. atra are negatively affected during warming EN episodes and positively affected or unaffected during cold LN episodes, as observed for other species with Antarctic or sub-Antarctic affinity (Riascos et al., 2009; Carstensen et al., 2010; Uribe et al., 2012). Moreover, this is in line with the negative population effects reported by Soenens (1985) and Díaz and Ortlieb (1993) for this species along the central Peruvian coast during previous EN episodes.

The recruitment process of A. atra showed a distinct response to EN–LN episodes – the density of recruits was higher during EN and lower during LN, as compared with neutral periods. It would not be reasonable to claim that the oceanographic anomalies associated with these episodes affect adults, larvae or recruits differently. A more plausible explanation for the increased recruitment during EN, despite a simultaneous drop in the SSB, seems to be a reduction in interspecific competitive interactions during the recruitment phase of A. atra on shallow benthic areas from Bahía Independencia that may overcompensate for the mortality associated with abiotic stress. This mechanism was described by Riascos et al. (2017) to explain the increased recruitment of the scallop A. purpuratus during EN in 1997–1998 in this bay; there was a decreased competition for food with filter-feeder polychaetes and a decreased predation by crabs, which were negatively affected during EN (Taylor, 2008).

Environmental Factors Affecting Early Stage Transitions of Aulacomya atra in Bahía Independencia

Bay systems in the HCS represent areas of propagule retention, prompted by upwelling flow structures such as filaments, squirts and shadow fronts (Marín et al., 2003; Thiel et al., 2007). Thus, circulation patterns and their effects on larval retention in Bahía Independencia are highly dependent on the temporal variability of the upwelling strength. The results of the PCA on the environmental variability at this bay indicates that the EN–LN cycle is the main source of variability in the oceanographic conditions, particularly thermal conditions. This result underscores one of the distinctive features of the HCS: unlike other coastal ecosystems at similar latitudes, the cold HCS exhibits stable oceanographic conditions, with little seasonal changes, that are heavily – often dramatically – modified during the EN–LN cycle (Arntz et al., 2006; Thiel et al., 2007; Carstensen et al., 2010).

Wind forcing in Bahía Independencia prompts the upwelling of equatorial subsurface waters characterized by high salinity, high nutrient levels and low dissolved oxygen. In the southern half of the bay (south of La Vieja Island) wind forcing generates a clockwise circulation pattern that favors shoreward water transport (Quispe et al., 2010). This circulation pattern seems to favor larval retention, as observed in other bay systems from the HCS (Thiel et al., 2007). With weakened upwelling-favorable winds during EN, there is an intrusion of subtropical, nutrient-depleted, warmer subtropical waters and a counter-clockwise circulation favoring offshore transport (Quispe et al., 2010). In this context we hypothesize that, when compared to normal or LN conditions, the effects of EN episodes for A. atra would be characterized by (i) nutrient-depleted waters that would reduce food availability, thus affecting body condition and the production of propagules, (ii) changed coastal circulation patterns that may result in offshore transport, thus reducing the larval pool available for settling and (iii) anomalous warmer waters that may change the balance of interspecific interactions that would ultimately increase the relative survival rate of mussels until successful recruitment.

Data Availability

All datasets generated for this study are included in the manuscript and/or the Supplementary Files.

Author Contributions

MS coordinated the monitory program. JR and MS created the databases, conceptualized the study, and analyzed the data. All authors drafted the first version of the manuscript, and revised and approved the final version.

Funding

MS was supported by several funding schemes including the United States Agency for International Development (USAID), the Alfred Wegener Institute for Polar and Marine Research (AWI), Germany, and the National University of San Marcos (Peru).

Conflict of Interest Statement

The authors declare that the research was conducted in the absence of any commercial or financial relationships that could be construed as a potential conflict of interest.

Acknowledgments

We thank all the members of the research group DePSEA (Desarrollo Pesquero en áreas someras del Ecosistema de Afloramiento) at the Laboratory of Marine Ecology at the National University of San Marcos for their collaboration during sampling surveys. Especially to biologist Andrés Belapatiño for his valuable collaboration in obtaining larval data and who was part of the research project in Independence Bay. We also thank Dr. Octavio Morón from the Instituto del Mar del Peru for providing oceanographic data, and Dr. Martín Ríos Alcolea from the University of Barcelona for his valuable statistical advice. This work is part of a monitoring program, which was led by Dr. Juan Tarazona (deceased).

Supplementary Material

The Supplementary Material for this article can be found online at: https://www.frontiersin.org/articles/10.3389/fmars.2019.00496/full#supplementary-material

Footnotes

1. ^ https://www.esrl.noaa.gov/psd/enso/mei.old/table.html
2. ^ https://www.cpc.ncep.noaa.gov/data/indices/ersst3b.nino.mth.81-10.rt.ascii

References

Arntz, W. E., and Fahrbach, E. (1996). El Niño: Experimento climático de la naturaleza Causas físicas y efectos biológicos. Mexico: Fondo de Cultura Económica., 312.

Google Scholar

Arntz, W. E., Gallardo, V. A., Gutiérrez, D., Isla, E., Levin, L. A., Mendo, J., et al. (2006). El Niño and similar perturbation effects on the benthos of the Humboldt, California, and Benguela Current upwelling ecosystems. Adv. Geosci. 6, 243–265. doi: 10.5194/adgeo-6-243-2006

CrossRef Full Text | Google Scholar

Ashton, T., Riascos, J. M., and Pacheco, A. (2008). First record of Cymatium keenae Beu, 1970 (Mollusca: Ranellidae) from Antofagasta Bay, northern Chile, in connection with El Niño events. Helgol. Mar. Res. 62, 107–110. doi: 10.1007/s10152-007-0084-5

CrossRef Full Text | Google Scholar

Avendaño, M., and Cantillánez, M. (2013). Reproductive cycle, collection and early growth of Aulacomya ater, Molina 1782 (Bivalvia: Mytilidae) in northern Chile. Aquacult. Res. 44, 1327–1338. doi: 10.1111/j.1365-2109.2012.03149.x

CrossRef Full Text | Google Scholar

Avendaño, M., and Cantillánez, M. (2014). Reproductive cycle of Aulacomya ater [Bivalvia: Mytilidae (Molina 1782)] in Punta Arenas Cove (Antofagasta Region. Chile). Aquacult. Int. 22, 1229–1244. doi: 10.1007/s10499-013-9743-5

CrossRef Full Text | Google Scholar

Bakun, A. (2010). Linking climate to population variability in marine ecosystems characterized by non-simple dynamics: conceptual templates and schematic constructs. J. Mar. Syst. 79, 361–373. doi: 10.1016/j.jmarsys.2008.12.008

CrossRef Full Text | Google Scholar

Barber, R. T., and Chavez, F. P. (1983). Biological Consequences of El Niño. Science 222, 1203–1210. doi: 10.1126/science.222.4629.1203

PubMed Abstract | CrossRef Full Text | Google Scholar

Beukema, J. J., and Dekker, R. (2007). Variability in annual recruitment success as a determinant of long-term and large-scale variation in annual production of intertidal Wadden Sea mussels (Mytilus edulis). Helgol. Mar. Res. 61, 71–86. doi: 10.1007/s10152-006-0054-3

CrossRef Full Text | Google Scholar

Bricelj, V. M., Epp, J., and Malouf, R. E. (1987). Intraspecific variation in reproductive and somatic growth cycles of bay scallops Argopecten irradians. Mar. Ecol. Prog. Ser. 36, 123–137. doi: 10.3354/meps036123

CrossRef Full Text | Google Scholar

Calienes, R., Guillén, O., and Lostaunau, N. (1985). Variabilidad espacio-temporal de clorofila, producción primaria y nutrientes frente a la costa peruana. Bol. Inst. Mar Perú 10, 1–44.

Google Scholar

Carstensen, D., Riascos, J. M., Heilmayer, O., Arntz, W. E., and Laudien, J. (2010). Recurrent, thermally-induced shifts in species distribution range in the Humboldt current upwelling system. Mar. Environ. Res. 70, 293–299. doi: 10.1016/j.marenvres.2010.06.001

PubMed Abstract | CrossRef Full Text | Google Scholar

Castilla, J. C., and Guiñez, R. (2000). Disjoint geographical distribution of intertidal and nearshore benthic invertebrates in the Southern Hemisphere. Rev. Chil. Hist. Nat. 73, 1–23. doi: 10.4067/S0716-078X2000000400004

CrossRef Full Text | Google Scholar

Chanley, P. (1968). “Larval development in the Class Bivalvia,” in Proceedings of the Symposium on Mollusca, part II, Cochin, (Mandapam Camp: Marine Biological Association of India), 475–481.

Google Scholar

Clarke, K. R., and Gorley, R. N. (2006). Primer V6: User Manual - Tutorial. Plymouth: Plymouth Marine Laboratory.

Google Scholar

Díaz, A., and Ortlieb, L. (1993). El fenómeno ‘El Niño’ y los moluscos de la costa peruana. Bull. Inst. Fr. Etudes Andines. 22, 159–177.

Google Scholar

Escribano, R., Fernández, M., and Aranís, A. (2003). Physical-Chemical processes and Patterns of Diversity of the chilean eastern boundary pelagic and benthic marine ecosystems: an overview. Gayana 67, 190–205.

Google Scholar

Eversole, A. G. (1989). “Chapter 3: Gametogenesis and Spawning in North American Clam Populations: Implications for Culture,” in Clam Mariculture in North America. Developments in Aquaculture and Fisheries Science, eds J. J. Mauzi and M. Costagna (Amsterdam: Elsevier), 75–110.

Google Scholar

Gayanilo, F. C., Sparre, P., and Pauly, D. (2005). ”FiSAT II - FAO-ICLARM stock assessment tools II: user’s guide”. Rev. version. Rome: Food and Agriculture Organization of the United Nations.

Google Scholar

Jaramillo, J. R., and Navarro, J. M. (1995). Reproductive cycle of the Chilean ribbed mussel Aulacomya ater (Molina, 1782). J. Shellfish Res. 14, 165–171.

Google Scholar

Levin, S. A. (1998). Ecosystems and the biosphere as complex adaptive systems. Ecosystems. 1, 431–436. doi: 10.1007/s100219900037

PubMed Abstract | CrossRef Full Text | Google Scholar

Lima, M., Brazeiro, A., and Defeo, O. (2000). Population dynamics of the yellow clam Mesodesma mactroides: recruitment variability, density-dependence and stochastic processes. Mar. Ecol. Prog. Ser. 207, 97–108. doi: 10.3354/meps207097

CrossRef Full Text | Google Scholar

Loveland, T. R., and Merchant, J. M. (2004). Ecoregions and ecoregionalization: geographical and ecological perspectives. J. Environ. Manage. 34, S1–S13. doi: 10.1007/s00267-003-5181-x

PubMed Abstract | CrossRef Full Text | Google Scholar

Marín, V. H., Delgado, L. E., and Escribano, R. (2003). Upwelling shadows at Mejillones Bay (northern Chilean coast): a remote sensing in situ analysis. Investig. mar. 31, 47–55. doi: 10.4067/S0717-71782003000200005

CrossRef Full Text | Google Scholar

Mendo, J. (1997). Estrategia para la Gestión Sustentable de Recursos pesqueros de la Bahía Independencia Pisco- Perú. IOC/INF 1046, 176–179.

Google Scholar

Mendo, J., and Wolff, M. (2003). El impacto de El Niño sobre la producción de concha de abanico (Argopecten purpuratus) en Bahía Independencia, Pisco, Perú. Ecol. Apl. 2, 51–57.

Google Scholar

Mildenberger, T. K., Taylor, M. H., and Wolff, M. (2017). TropFishR: an R package for fisheries analysis with length-frequency data. Methods Ecol. Evol. 8, 1520–1527. doi: 10.1111/2041-210X.12791

CrossRef Full Text | Google Scholar

Mildenberger, T. K., Taylor, M. H., and Wolff, M. (2018). TropFishR: Tropical Fisheries Analysis with R. R package version 1.2.1. Available at: https://CRAN.R-project.org/package=TropFishR (accessed April 2, 2019).

Google Scholar

Morón, O. (2011). Climatología de la Salinidad Superficial del Mar frente a la costa peruana. 1960 – 2008. Inf. Inst. Mar. Perú 38, 7–39.

Google Scholar

Myers, R. H. (1986). Classical and modern regression with applications. Boston, MA: Duxbury.

Google Scholar

Ortega, L., Castilla, J. C., Espino, M., Yamashiro, C., and Defeo, O. (2012). Effects of fishing, market price, and climate on two South American clam species. Mar. Ecol. Prog. Ser. 469, 71–85. doi: 10.3354/meps10016

CrossRef Full Text | Google Scholar

Pauly, D. (1982). “Studying single-species dynamics in a tropical multispecies context,” in proceedings of the ICLARM Conference on the Theory and management of tropical fisheries, (Manila: International Center for Living Aquatic Resources Management), 33–70.

Google Scholar

Pauly, D. (1983). Some simple methods for the assessment of tropical fish stocks. FAO Fisheries Technical Paper No. 234. Rome: FAO, 52.

Google Scholar

Peña, T. S., Johst, K., Grimm, V., Arntz, W., and Tarazona, J. (2005). Population dynamics of a polychaete during three El Niño events: disentangling biotic and abiotic factors. Oikos 111, 253–258. doi: 10.1111/j.0030-1299.2005.14067.x

CrossRef Full Text | Google Scholar

Quispe, D., Graco, M., Correa, D., Tam, J., Gutiérrez, D., Morón, O., et al. (2010). Variabilidad espacio-temporal de condiciones hidrofísicas en Bahía Independencia, Pisco-Perú, desde 1955 al 2004. Ecol. Apl. 9, 9–18.

Google Scholar

Ramorino, L., and Campos, B. (1983). Larvas y post larvas de Mytilidae de Chile (Mollusca: Bivalvia). Rev. Biol. Mar. Valparaíso. 19, 143–192.

Google Scholar

Riascos, J. M. (2009). Thriving and declining: population dynamics of the macha (Mesodesma donacium, Bivalvia, Mesodesmatidae) along a latitudinal gradient of the Humboldt Current Upwelling System. Ph.D. thesis, University of Bremen: Bremen.

Google Scholar

Riascos, J. M., Avalos, C. M., Pacheco, A. S., and Heilmayer, O. (2012). Testing stress responses of the bivalve Protothaca thaca to El Niño-La Niña thermal conditions. Mar. Biol. Res. 8, 654–661. doi: 10.1080/17451000.2011.653367

CrossRef Full Text | Google Scholar

Riascos, J. M., Carstensen, D., Laudien, J., Arntz, W. E., Oliva, M. E., Güntner, A., et al. (2009). Thriving and declining: climate variability shaping life-history and population persistence of Mesodesma donacium in the Humboldt upwelling System. Mar. Ecol. Prog. Ser. 385, 151–163. doi: 10.3354/meps08042

CrossRef Full Text | Google Scholar

Riascos, J. M., Heilmayer, O., Oliva, M. E., and Laudien, J. (2011). Environmental stress and parasitism as drivers of population dynamics of Mesodesma donacium at its northern biogeographic range. ICES J. Mar. Sci. 68, 823–833. doi: 10.1093/icesjms/fsr026

CrossRef Full Text | Google Scholar

Riascos, J. M., Solís, M. A., Pacheco, A. S., and Ballesteros, M. (2017). Breaking out of the comfort zone: El Niño-Southern Oscillation as a driver of trophic flows in a benthic consumer of the Humboldt current ecosystem. Proc. R. Soc. B. 284:20170923. doi: 10.1098/rspb.2017.0923

PubMed Abstract | CrossRef Full Text | Google Scholar

Riascos, J. M., Uribe, R., Donayre, S., Flores, D., Galindo, O., Quispe, C., et al. (2016). Human footprints on benthic communities in marine reserves: a study case in the most productive upwelling system worldwide. Mar. Ecol. Prog. Ser. 557, 65–75. doi: 10.3354/meps11857

CrossRef Full Text | Google Scholar

Ryther, J. H. (1969). Photosynthesis and Fish Production in the Sea. Science 166, 72–76. doi: 10.1126/science.166.3901.72

PubMed Abstract | CrossRef Full Text | Google Scholar

Sinclair, M. (2009). Marine populations: an essay on population regulation and speciation. J. Mar. Biol. Assoc. UK. 69, 247–247. doi: 10.1017/S0025315400049353

CrossRef Full Text | Google Scholar

Soenens, P. (1985). “Estudios preliminares sobre el efecto del fenómeno “El Niño” 1982-1983 en comunidades de Aulacomya ater,” in El fenómeno El Niño Su impacto en la fauna marina, eds W. E. Arntz, A. Landa, and J. Tarazona (Callao: Insituto del mar del Peru), 51–53.

Google Scholar

Takahashi, K., Mosquera, K., and Reupo, J. (2014). “El Índice Costero El Niño (ICEN): historia y actualización,” in Boletín Técnico “Generación de modelos climáticos para el pronóstico de la ocurrencia del Fenómeno El Niño, ed. K. Takahashi (Lima: Instituto Geofísico del Perú), 8–9.

Google Scholar

Tarazona, J., and Arntz, W. (2001). “The Peruvian Coastal Upwelling System,” in Coastal Marine Ecosystems of Latin America, eds U. Seeliger and B. Kjerfve (Berlin: Springer Berlin Heidelberg), 229–244. doi: 10.1007/978-3-662-04482-7_17

CrossRef Full Text | Google Scholar

Taylor, M. H. (2008). The Northern Humboldt Current Ecosystem and its resource dynamics: insights from a trophic modeling and time series analysis. Ph.D. thesis, University of Bremen: Bremen.

Google Scholar

Thiel, M., Macaya, E. C., Acuna, E., Arntz, W. E., Bastias, H., Brokordt, K., et al. (2007). The Humboldt Current System of northern and central Chile: oceanographic processes, ecological interactions and socioeconomic feedback. Oceanogr. Mar. Biol. Ann. Rev. 45, 195–344. doi: 10.1201/9781420050943.ch6

CrossRef Full Text | Google Scholar

Trillmich, F. (1991). “El Niño in the Galapagos Islands: a natural experiment,” in Ecosystem Experiments, eds H. A. Mooney, E. Medina, D. W. Schindler, E.-D. Schulze, and B. H. Walker (Hoboken, NJ: John Wiley and Sons).

Google Scholar

Uribe, R. A., Oliva, M. E., Aguilar, S., Yamashiro, C., and Riascos, J. M. (2012). Latitudinal variation in the reproductive cycle of two bivalves with contrasting biogeographical origin along the Humboldt current upwelling ecosystem. Sci. Mar. 76, 713–720. doi: 10.3989/scimar.03480.26D

CrossRef Full Text | Google Scholar

Vélez, J. A., Watson, W., Arntz, W., Wolff, M., and Schnack-Schiel, S. B. (2005). Larval fish assemblages in Independencia Bay, Pisco, Peru: temporal and spatial relationships. Mar. Biol. 147, 77–91. doi: 10.1007/s00227-004-1535-z

CrossRef Full Text | Google Scholar