Daily anchovy (Engraulis ringens) egg mortality rates in different spawning zones along the Chilean coast

Regional Studies in Marine Science 35 (2020) 101200

Contents lists available at ScienceDirect

Regional Studies in Marine Science journal homepage: www.elsevier.com/locate/rsma

Daily anchovy (Engraulis ringens) egg mortality rates in different spawning zones along the Chilean coast ∗

Blanca Bustos a,b , Luis A. Cubillos b , , Gabriel Claramunt c,d , Leonardo R. Castro e a

Programa Magister en Ciencias con mención en Pesquerías, Facultad de Ciencias Naturales y Oceanográficas, Universidad de Concepción, Casilla 160-C, Concepción, Chile b COPAS Sur-Austral, EPOMAR, Departamento de Oceanografía, Facultad de Ciencias Naturales y Oceanográficas, Universidad de Concepción, Casilla 160-C, Concepción, Chile c Facultad de Recursos Naturales Renovables, Universidad Arturo Prat, Casilla 121, Iquique, Chile d Instituto de Ciencias Exactas y Naturales, Universidad Arturo Prat, Iquique, Chile e COPAS Sur-Austral, Departamento de Oceanografía, Facultad de Ciencias Naturales y Oceanográficas, Universidad de Concepción, Chile

article

info

Article history: Received 7 September 2018 Received in revised form 26 February 2020 Accepted 26 February 2020 Available online 29 February 2020 Keywords: Egg survival Patchiness index Spawning area Dispersion Aggregation Small pelagic fish Northern Chile Central-southern Chile

a b s t r a c t The daily egg mortality (Z ) of small pelagic fish can be modulated by environmental variables as well as by population derived factors such as egg production and/or dispersion and aggregation processes. The anchovy daily egg mortality rate was compared across spawning areas distributed in the North (18◦ 25’S-26o 03’S), Central (34o 30’S-37o 10’S), and South (38◦ S-41◦ 20’S) zones along the Chilean coast. Daily egg production (P0 ) and Z were estimated through a generalized linear model (binomial negative) from egg count data obtained in yearly Daily Egg Production Method applications (2002– 2014). Co-variables included P0 , the spawning area (SA), sea surface temperature (SST ), phytoplankton biomass (Chl-a), alongshore wind (V ) and an index of the El Niño-Southern Oscillation (ENSO) during the spawning period (August–September).Egg mortality was significantly higher in the North zone, influencing spatial differences rather than yearly effects. Egg mortality was not associated with ENSO, but correlated positively with P0 , SA, SST, and V, and negatively with Chl-a. The best model explained 42.3% (P<0.01) of the egg mortality rates and included SST and Chl-a. In the South zone, the egg mortality was significantly correlated with P0 and V, which was interpreted as density-dependent effects and winds driving egg concentration processes rather than dispersion. The Lloyd’s patchiness index was computed as an index of egg aggregation by embryonic development stages (Group A: early development, Group B: middle development, and Group C: close to hatching). The spawning area was negatively correlated with the patchiness index of Groups A and B in the North zone. Overall results suggested that anchovy eggs survival would be dependent on dispersion in the North zone but concentration and retention in the South spawning area and thereafter coincident with reproductive tactics of adults in the North and South areas, respectively. © 2020 Published by Elsevier B.V.

1. Introduction The anchovy (Engraulis ringens Jenyns, 1842) has a wide distribution in the Southeastern Pacific Ocean, from northern Peru (4◦ 00′ S) to southern Chile (42◦ 00′ S) (Serra et al., 1979), even extending into fjords and inner waters of the Chilean Patagonia (Bustos et al., 2008). Four stocks of anchovy are identified along its distribution such as the north-central stock off Peru (4◦ S– 14◦ S), the shared stock between south Peru and northern Chile (16◦ S–24◦ S), the central-northern stock (25◦ S–32◦ 10′ S), and the central-southern stock (33◦ S–42◦ S) (Alheit and Ñiquen, 2004; Valdivia et al., 2007; Canales and Leal, 2009). ∗ Corresponding author. E-mail address: [email protected] (L.A. Cubillos). https://doi.org/10.1016/j.rsma.2020.101200 2352-4855/© 2020 Published by Elsevier B.V.

The anchovy stocks seem to fluctuate synchronously in the long term, which is associated with decadal to interdecadal variability occurring in the Humboldt Current System (HCS) (Chavez et al., 2003; Alheit and Ñiquen, 2004; Cubillos et al., 2007a). The HCS is one of the most productive upwelling ecosystem sustaining important fisheries (Chavez and Messié, 2009). The life cycle of anchovy occurs in habitats characterized by coastal upwelling (Silva et al., 2016), trapped waves, and active mesoscale eddies affecting the surface circulation field, sea surface temperature and distribution of the planktonic biomass (Correa-Ramirez et al., 2012). In addition, hypoxic waters associated with fluctuations of the oxygen minimum zone could be important for the anchovy population dynamics (Morales et al., 1996) and species alternations (Bertrand et al., 2011). In the North zone of Chile, the south winds are persistent year round and generate offshore

2

B. Bustos, L.A. Cubillos, G. Claramunt et al. / Regional Studies in Marine Science 35 (2020) 101200

Ekman transport and coastal upwelling events. In Central and South zones of Chile, winds favoring upwelling are seasonal from September to March; i.e., during the austral spring and summer (Arcos et al., 1987; Castro et al., 2002; Sobarzo et al., 2007). El Niño events affect the seasonal variation of the wind off central Chile, delaying the season of upwelling events and the distribution of the phytoplankton biomass (Gómez, 2007; Montecinos and Gómez, 2010). In the coastal zone, the kinetic energy associated with mesoscale eddies presents interannual variability linked to El Niño events (Hormazabal et al., 2004). The application of the Daily Egg Production Method (DEPM) for spawning biomass estimations has allowed verifying patterns of variability in the spawning distribution of the anchovy along the Chilean coast (Braun et al., 2005; Claramunt et al., 2012). The spawning distribution of the anchovy shows spatial and temporal variability in the North zone off Chile (Claramunt et al., 2012), which contrasts with more persistent locations of coastal egg aggregations in central-southern Chile (Cubillos et al., 2007b; Castillo-Jordán et al., 2007). In the central-southern stock, Castro et al. (1997) and Cubillos et al. (2007b, 2013) identified two main spawning zones of anchovy separated at Punta Lavapié (37◦ 10′ S): the central spawning zone (33◦ 30′ S–37◦ 10′ S) and the south spawning zone (38◦ 00–41◦ 30′ S). Comparative analyses of environmental conditions affecting early life-history traits of the anchovy in the North, Central and South spawning zones off Chile have revealed important latitudinal gradients in biological and oceanographic features (Parrish et al., 1983). Egg size and their lipidic content, and the yolk volume of recently hatched larvae tend to increase from north to south; in contrast, the yolk-sac larvae instantaneous growth rate varied among subpopulations as a function of temperature, which decreases towards higher latitudes (Llanos-Rivera and Castro, 2006; Castro et al., 2009). In addition, in the ovary, oocyte sizes are smaller and fecundity is higher in the north as compared with central-south female anchovies (Castro et al., 2009; Claramunt et al., 2012). These characteristics, fecundity, oocite size and early life stage’s traits are maternally influenced and ultimately driven by local environmental conditions at the spawning zones along the Chilean coast (Castro et al., 2009). Claramunt et al. (2012) concluded that the environmental conditions varied notably among spawning zones and seasons and that the reproductive traits (oocyte size and fecundity) responded to environmental conditions within a couple of weeks and were also manifested in the daily spawning fraction estimates. In their study, phytoplankton biomass (measured as concentration of chlorophyll-a) was the main environmental variable explaining changes in the location of the spawning areas rather than temperature (Claramunt et al., 2012). Different methods have been utilized to estimate anchovy daily early life stage mortality rates in the HCS. By following daily spawned cohorts at the beginning of the winter spawning season (Essig and Cole, 1986), Castro and Hernandez (2000) estimated mean seasonal daily egg and larval (combined) mortality rates of up to 96%–98% in a coastal location off central Chile. A second way to estimate the instantaneous daily egg mortality rate is by means of the mortality parameter estimated during the application of the Daily Egg Production Method (DEPM). This method has been utilized for spawning biomass anchovy estimations in the North (Claramunt et al., 2012), Central and South spawning zones off Chile for over a decade (Cubillos et al., 2015) and until now, no studies have been oriented to study the spatial and temporal variability in the daily egg mortality rate, neither a comparative study has addressed potential associations among egg mortality estimates (or survival) and environmental variables during the reproductive seasons either at each zone or among zones.

It has long been accepted that daily egg mortality rate of small pelagic fishes could be modulated either by environmentally driven or density-dependent factors (Rojas et al., 2011; Xu et al., 2015). In the HCS, El Niño-Southern Oscillation, Ekman transport associated to alongshore wind, as well as sea surface temperature are candidate environmental co-variables that could affect the fate of the eggs. Density-dependent variables such as the daily egg production, the spawning area covered by adults, and the food available for adults (phytoplankton biomass, Claramunt et al., 2012) might also affect egg chances of survival. In addition, variations in patchiness and dispersion of the egg distributions might have effects on predators and thereafter on egg mortality (McGurk, 1986) The objective of this study was to determine whether differences in the daily egg mortality rate of anchovy occur among three spawning zones along the Chilean coast, and their potential relationships with density-independent and density-dependent variables. The seasonal patterns of spawning in the studied zones tend to be similar, peaking from August to December (Claramunt et al., 2014). The selected variables were local temperature (because of its influence on egg development rates), phytoplankton biomass (as a proxy for adult food availability), alongshore wind (as a proxy of upwelling) during the reproductive peak season, as well as large-scale anomalies associated to El Niño-Southern Oscillation events. Associations with density-dependent biological processes in each zone included were the daily egg production, size of the spawning area, and egg patchiness. 2. Materials and methods 2.1. Study area and data sources We used the data of egg abundance per embryonic development stages, which were generated by application of the Daily Egg Production Method (DEPM) from 2002 to 2014 (Cubillos et al., 2007b; Claramunt et al., 2012; Cubillos et al., 2015). The data were provided by FIPA (Fondo de Investigación Pesquera) and IFOP (Instituto de Fomento Pesquero), and summarized in Table 1. The ichthyoplankton surveys were carried out in the North (18◦ 25′ S–26◦ 03′ S), Central (34◦ 30′ S–37◦ 10′ S) and South (38◦ S– 41◦ 20′ S) zones along the Chilean coast (Fig. 1). The North zone is large and characterized by a straight coastal-line north–south oriented and narrow continental shelf. The Central zone has a wide continental shelf interrupted by two submarine canyons that are associated with important rivers: Biobío and Itata, and an irregular coastal-line with three northward oriented bays (Sobarzo et al., 2007). The South zone is moderate in size (intermediate compared with the north and central zones), wide continental shelf and numerous rivers, providing runoff (Fig. 1). 2.2. Egg sampling in the daily egg production method surveys In the North zone, the survey design consisted of stations separated by 5 nautical miles (nmi) in the cross-shore direction over perpendicular transects 10 nmi apart. In the Central and South zones, the survey designs were based on transects separated by 5 nmi and stations separated by 4 nmi. The realized surveys, in terms of the number of stations, number of stations with presence of eggs (>1 egg 0.05 m−2 ), average egg density (egg 0.05 m−2 ), study area and spawning area are summarized in Table 1. Anchovy eggs were collected from vertical tows from a maximum depth of 70 m in deep stations on following conventional sampling design for anchovy’s DEPM (Smith and Hewitt, 1985; Santander et al., 1984; Cubillos et al., 2007b). In shallow waters, eggs were collected from 5 m above the sea bottom. A Pairovet net was used in the surveys, consisting of two CalVET

B. Bustos, L.A. Cubillos, G. Claramunt et al. / Regional Studies in Marine Science 35 (2020) 101200

3

Fig. 1. Study area showing the approximate boundaries of spawning areas in the North (18◦ 25′ S-26◦ 03′ S), Central (34◦ 30′ S–37◦ 10′ S) and South (38◦ 00′ S–41◦ 20′ S) zones off Chile. Table 1 Description of the egg surveys carried out in the three spawning zones of anchovy (Engraulis ringens) along the Chilean coast, from 2002 to 2014. Zona

Year

Stations (n)

Stations with eggs (n)

Abundance of eggs (n)

Average density (egg 0.05 m−2 )

Average density in stations with eggs (eggs 0.05 m−2 )

Survey area (km2 )

Spawning area (km2 )

North

2002 2004 2005 2007 2008 2009 2011 2012 2013 2014

589 649 658 799 717 773 730 564 770 632

310 211 215 113 136 41 292 211 201 128

18 273 12 367 17 059 4724 10 304 4702 13 914 10 779 6283 7458

31.0 19.1 25.9 5.9 14.4 6.1 19.1 19.1 8.2 11.8

58.9 58.6 79.3 41.8 75.8 114.7 47.7 51.1 31.3 58.3

102 419 109 464 111 427 130 668 119 221 130 764 119 813 119 098 122 693 101 167

60 811 39 997 41 861 26 033 30 482 10 387 59 224 62 675 38 106 23 728

Central

2002 2004 2005 2007 2008 2009 2011 2012 2013 2014

308 137 130 180 179 164 158 174 178 149

82 49 39 62 9 35 50 6 40 43

1242 893 1815 1116 24 434 1187 58 593 1210

4.0 6.5 14.0 6.2 0.1 2.6 7.5 0.3 3.3 8.1

15.1 18.2 46.5 18.0 2.7 12.4 23.7 9.7 14.8 28.1

9175 7944 8678 12 000 12 334 11 246 10 431 11 642 12 023 9558

3491 2856 2585 4200 665 2466 3277 395 2681 2680

South

2002 2004 2005 2007 2008 2009 2011 2012 2013 2014

358 143 140 203 184 203 199 182 174 211

103 51 63 129 73 75 57 34 23 60

3552 2713 4856 19 499 6115 2560 1510 736 638 1496

9.9 19.0 34.7 96.1 33.2 12.6 7.6 4.0 3.7 7.1

34.5 53.2 77.1 151.2 83.8 34.1 26.5 21.6 27.7 24.9

13 775 9734 9057 13 457 12 152 13 550 13 619 12 406 11 721 14 232

5882 3194 4138 8582 4892 5169 3858 2312 1508 4157

net (California Vertical Egg Tow, Smith et al., 1985) with a mouth

Once the plankton sample was on board, the plankton collected

area of 0.05 m , 150 µm mesh and equipped with a depressor

were preserved in a 4% formaldehyde solution (10% formalin) in

of 36 kg. The velocity of each vertical tow was 70 m per minute.

seawater buffered with borax (tetraborate of sodium).

2

4

B. Bustos, L.A. Cubillos, G. Claramunt et al. / Regional Studies in Marine Science 35 (2020) 101200

2.3. Embryonic stages, aging, daily egg production and mortality The plankton samples were analyzed in the laboratory. Eggs and larvae of anchovy were sorted and counted under a stereoscopic microscope. The eggs were classified in 11 embryonic development stages according to Moser and Ahlstrom (1985) and Santander et al. (1984). The time of spawning was set as at midnight (00:00 h) (Claramunt et al., 2007). Aging of eggs was following Lo (1985) by using a stage-temperature dependent model obtained from incubation experiments (Claramunt et al., 2007). Assuming that the egg abundance declined exponentially, the daily egg production and mortality rate were estimated by the following expression: Pt = P0 e−Zt

(1)

where Pt is the egg abundance at age t (egg per 0.05 m2 per day), P0 is the daily egg production rate (egg 0.05 m−2 d−1 ) and Z is the daily instantaneous total mortality rate of eggs (d−1 ). Ages less than four hours after spawning were discarded, as well as older than 95% of hatching ages. P0 and Z were estimated through generalized linear models using negative binomial and link log, which was appropriated for variation of egg counts as a function of age (Bernal et al., 2001). The negative binomial distribution has a dispersion parameter, which was estimated iteratively according to fitting procedures described by Venables and Ripley (2002) for the package ‘MASS’ for R (R Core Team, 2019). The estimates of the daily egg production and mortality rate of anchovy were conditioned to the number of stations in the spawning area (SA). Stations with presence of at least one egg (egg 0.05 m−2 ) delimited the SA, which was estimated by searching a radius equivalent to the half of distance between stations. The radius allowed connecting neighboring stations by using the ‘geofun’ package of the ‘egg’ library (http://sourceforge.net/projects/ichthyoanalysis/), available for the language and software R (Stratoudakis et al., 2006). 2.4. Egg patchiness per stages In order to obtain a patchiness index for egg counting, the eggs number from eleven embryonic stages was grouped into three categories. Eggs in the stages I–IV (early development stages) were grouped as Group A, where most eggs were less than 24 h. Eggs in stages V–VIII (early development stages) was grouped as Group B and consist of egg less than 48 h and older than 24 h. The Group C consisted of eggs in stages IX to XI, which are close to hatching. The Lloyd index (Lloyd, 1967) of patchiness was computed according to McGurk (1986), i.e., p=1+(

σ2 − 1)/µ µ

(2)

where p is the patchiness index, σ 2 is the variance of egg number and µ is the mean of egg number (egg 0.05 m−2 ). The index of patchiness was computed for each group and was used to establish a relationship with Z under the assumption that spatial aggregation explains the mortality rate of pelagic eggs following McGurk (1986, 1987). 2.5. Environmental variables In order to analyze remote effects of environmental variability on the daily egg mortality of anchovy, El Niño 3+4 SST anomaly index was utilized. This index is an indicator of central tropical Pacific El Niño conditions. It is computed from sea surface temperature data in the box between 170◦ W–120◦ W and 5◦ S–5◦ N, and was obtained from the NOAA Climate Prediction Center (http://www.cpc.ncep.noaa.gov/data/indices/). Local

environmental variables were sea surface temperature (SST ) and phytoplankton biomass (Chlorophyll-a, Chl-a), which were obtained from daily satellite images (resolution 4 × 4 km) of MODIS Aqua Level-3 (http://oceancolor.gsfc.nasa.gov) for the period that covered each survey. Satellite images were processed, scaled and atmospherically corrected to obtain geophysical units for temperature (◦ C) and phytoplankton biomass (mg m−3 ). Applying a spatiotemporal interpolation based on cokriging avoided the clouds present in satellite images (Marcotte, 1991; Navarro et al., 2004). Alongshore winds (V, m s−1 ) were based on re-analysis data, which are available in the NOAA website (https://www.esrl. noaa.gov/psd/). Unfortunately, the origin of wind data did not allow computing an upwelling index, but we utilized alongshore winds as a proxy for upwelling. Indeed, the alongshore wind is the main input for wind stress, and hence for an upwelling index (Bakun, 1973; Jacox et al., 2018). 7 2.6. Statistical analysis The available data for three spawning zones (north, central and south) were those for 2002, 2004, 2005, 2007–2009, and 2011–2014 (n = 10 years and 3 zones). The mean of each environmental variable was obtained during August–October for each zone, which represented the main spawning season in each zone (Claramunt et al., 2014). The daily egg mortality rate (Z) followed a lognormal distribution (Shapiro–Wilk’s W test: W = 0.936, P = 0.07), which is reasonable since Z values pertain to positive real numbers. We computed the daily survival rate (S) according to the expression S = e−Z for better interpretation of Z values. In order to evaluate spatial and temporal effects on logtransformed Z (logZ ), we used two-way ANOVA with the zone and years as main factors. After, the Tukey’s Honest Significant Difference method was applied for pair-wise comparison. Trend analysis was based on an ANCOVA, which allowed to comparing differences in the slope and intercept. In this case, the years were considered to be a continuous variable because the interaction between factors (Year-Zone) was not possible due to the reduced degree of freedom (10 years and 3 zones). As an alternative, we applied the Kruskal–Wallis test to evaluate the interaction Zone-Year, then the Dunn test for multiple comparisons (Zar, 2010). The association among variables was obtained through Pearson’s correlation matrix considering logZ as dependent variable and logPo, SA, ENSO, SST, Chl-a, and V as independent variables. Once identified significant correlations, we used linear regression to evaluate the effects of one or more environmental variables on logZ. Candidate models were sorted according to corrected Akaike’s Information Criterion (AICc) (Burnham and Anderson, 2002). In order to explain the effects of environmental variables on mortality, we studied the correlations between the Patchiness Index and spawning area (SA), as well as with the egg production (logPo) for each area. These analyses help us to interpret density-dependent effects on egg mortality. Statistical analyses were performed with R (R Core Team, 2019). In addition to specific packages already mentioned, the following packages were also utilized: ggplot2 (Wickham, 2016), ggpubr (Kassambara, 2018), FSA (Ogle et al., 2018) and rcompanion (Mangiafico, 2019). 3. Results The daily egg mortality rate (Z ) of anchovy fluctuated between 0.45 d−1 and 1.47 d−1 in the North zone, meaning daily egg survival rate (S) between 63.8 and 23.0%. Z was lower in the Central zone, ranging between 0.13 d−1 and 0.71 d−1 (survival

B. Bustos, L.A. Cubillos, G. Claramunt et al. / Regional Studies in Marine Science 35 (2020) 101200

5

Table 2 Estimates of the daily egg mortality rate of anchovy (Z ) and its standard error (Std. Err.), as well as the daily survival rate (S) across year (2002–2014) and spawning zone (North, Central, South). Year

2002 2004 2005 2007 2008 2009 2011 2012 2013 2014

North

Central

South

Z (d−1 )

Std. Err.

S (%)

Z (d−1 )

Std. Err

S (%)

Z (d−1 )

Std. Err.

S (%)

0.45 0.70 1.09 0.80 0.73 1.12 0.77 0.64 1.47 0.67

0.15 0.19 0.18 0.24 0.22 0.42 0.14 0.18 0.24 0.25

63.8 49.7 33.6 44.9 48.2 32.6 46.3 52.7 23.0 51.2

0.59 0.13 0.49 0.24 0.29 0.14 0.47 0.71 0.14 0.16

0.20 0.20 0.31 0.19 0.35 0.28 0.22 0.87 0.26 0.29

55.4 87.8 61.3 78.7 74.8 86.9 62.5 49.2 86.9 85.2

0.27 0.40 0.71 0.46 0.60 0.45 0.40 0.46 0.69 0.40

0.20 0.25 0.24 0.13 0.21 0.19 0.23 0.30 0.35 0.21

76.3 67.0 49.2 63.1 54.9 63.8 67.0 63.1 50.2 67.0

rate between 87.8 and 49.2%). In the South zone, the daily egg mortality rate ranged between 0.27 d−1 and 0.71 d−1 ; i.e., survival rate between 76.3 and 49.2% (Table 2). The daily egg mortality rate showed significant differences among zones (P < 0.05, df = 2, Fig. 2A), but not among years (P = 0.110, df = 9, Fig. 2B). The daily egg mortality rate was higher in the North zone (Fig. 2B, Table 2), and the Tukey’s post-hoc test showed that the main differences were due to the North zone (North-Central: P < 0.05, Central-South: P = 0.324, north–south: P < 0.05). Although the daily egg mortality rate showed variability in time (Fig. 2B), there were non-significant differences for the interaction Year-Zone (Kruskal–Wallis test: P = 0.465, df = 29). Indeed, the ANCOVA test for changes in the slope given the factor Zone confirmed that daily egg mortality rate did not show any trend in the time series here analyzed (cross factor Year * Zone: F = 1.92, P = 0.152), as well as significant differences in the intercepts associated to the factor Zone (F = 13.33, P < 0.001). In general, the daily egg production of anchovy was higher and less variable in the North zone (CV = 5.8%) as compared to the Central (CV = 20.9%) and South (CV = 12.7%) zones (Fig. 3A). The spawning area was bigger in the North zone also, with the smallest spawning areas in 2009 (Table 1 and Fig. 3A). In the Central and South zones, the spawning area was smaller and presented a slight trend to decrease from 2007 to 2012 (Fig. 3B). Phytoplankton biomass, during the reproductive period, was higher in the Central zone and variable in the North zone (Fig. 3C). The phytoplankton biomass showed a coincident decline in 2012 throughout the studied area (Fig. 3C). Inter-annual changes in SST were similar among zones, with warmer waters in the North zone and cooler in the South zone, and colder years were 2007 and 2013 (Fig. 3D). The alongshore wind (V ), favorable for upwelling, showed higher intensity in the North zone with a maximum in 2012. Instead, in the Central and South zones, the wind intensities varied similarly across the zones (Fig. 3E). The wind intensity was lower in the Central and South zones compared with the North one, because the alongshore wind was negative sometimes and reduced the mean during the reproductive season. Northern winds are non-favorable for upwelling, meaning convergence at the coast during the spawning period (Sobarzo et al., 2007). El Niño 3 + 4 showed warmer years in 2002, 2004, 2009, 2012 and 2014, while colder years were 2007 and 2011 (Fig. 3F). Without taking into account the zones, we found significant positive correlations between Z and the P0 (r = 0.751), between Z and the spawning area (r = 0.468), between Z and SST (r = 0.578), and between Z and alongshore wind (r = 0.379), but negative correlation between Z and the phytoplankton biomass (r = −0.532) (Table 3). Within each zone, Z was not associated with the environmental variables, excepting in the South zone in which Z showed a significant correlation with P0 and V (Table 3). The effects of the environmental variables (SST, Chl-a and V ) on logZ were evaluated by considering a set of linear models

Table 3 Matrix of Pearson’s correlation of bivariate relationship between the daily egg mortality rate of anchovy (Z, log-transformed) and covariables, where logPo: logarithm of daily egg production, SA: spawning area, ENSO: El Niño-Southern Oscillation variability index, SST : sea surface temperature, Chl-a: phytoplanktonic biomass, V : alongshore wind (Period 2002–2014, n = 10, df = degree of freedom). Zone

logPo

SA

ENSO

SST

Chl-a

V

Total (df = 29) North (df = 8) Central (df = 8) South (df = 8)

0.751***

0.468*

−0.137

0.578***

−0.532**

0.379**

0.467

−0.465

−0.353

0.007

0.039

−0.026

0.388

−0.225

−0.023

0.049

−0.104

−0.306

0.665*

−0.327

−0.455

−0.295

0.195

0.784**

*P < 0.05. **P < 0.01. ***P < 0.001.

(Table 4). SST and Chl-a explained 42.3% (P < 0.01) of the egg mortality variance (model M4), followed by model M6 in which Chl-a and V explained 39.3% (P = 0.001) of the egg mortality variance. SST was the single variable that explained 33.5% of the egg mortality variance (P < 0.001) (Table 4). Note that previous significant correlations between egg mortality rate and co-variables (excepting ENSO, Tables 3 and 4) were consequences of nested data that affected the correlations within each zone, particularly in the case of the spawning area (Fig. 4B), sea surface temperature, phytoplankton biomass and alongshore wind (Fig. 5). Nevertheless, daily egg production had the same positive relationship within each zone (Fig. 4A). The association between the spawning area (log SA) and the patchiness index, per embryonic development stages, showed a significant and negative correlation in the case of groups A and B in the North zone only. Contrarily, the correlation was significant but positive in the case of the egg production rate (Table 5). The relationship between the spawning area and the patchiness index of Group A showed that dispersion (low patchiness index) of early development egg stages was associated with a higher spawning area in the North zone only. On the contrary, higher egg patches tend to occur at lower spawning areas (Fig. 6). 4. Discussion Our objective was to evaluate differences in the daily anchovy egg mortality rates in three different spawning zones distributed along the Chilean coast and to assess the potential relationships with density-independent and density-dependent co-variables. The daily egg mortality rate can vary from year to year within a spawning area, but the effects associated with each spawning zone were the most significant. The daily egg mortality rate

6

B. Bustos, L.A. Cubillos, G. Claramunt et al. / Regional Studies in Marine Science 35 (2020) 101200

Fig. 2. Spatial and temporal comparisons of the daily egg mortality rate (Z) of anchovy among zones (A) and across years (B). In the boxplot, small vertical lines represent minimum and maximum values, the horizontal border of each box are the first and third quartile, and the horizontal line inside the box is the median. The times series shows Z values and the standard error of the estimates.

Fig. 3. Time series of the logarithm of egg production (logPo) (A), spawning area (SA) (B), phytoplanktonic biomass (Chl-a) (C) Sea surface temperature (SST ) (D); alongshore wind intensity (V ) (E) and ENSO index (F).

was computed from a standard methodology utilized in the application of the Daily Egg Production Method, and thereafter, variability associated with estimation error must be lower. Of

course, the temperature is an important factor affecting egg development (Lo, 1985), and higher egg mortality in the North zone

B. Bustos, L.A. Cubillos, G. Claramunt et al. / Regional Studies in Marine Science 35 (2020) 101200

7

Fig. 4. Bivariate relationships between logarithm of the daily egg mortality rate (logZ ) and density-dependent co-variables within each of the spawning zones (North, Central and South zones): logZ vs. logarithm of the daily egg production rate (A), Z vs. Spawning area (B). Table 4 Candidate linear models used to explain the relationship between daily egg mortality rate (logZ ) and the environmental variables sea surface temperature (SST ), phytoplanktonic biomass (Chl-a) and alongshore wind (V ). P: number of estimated parameters, df: degree of freedom, R2: coefficient of determination, F: Fisher test of analysis of variance, Pr(>F): probability of the F-test, AICc: corrected Akaike’s Information Criterion, ∆AICc: difference in AICc taken into account the best model in bold. Model

Expression

p

df

R2

F

Pr(>F)

AICc

∆AICc

M1 M2 M3 *M4 M5 M6 M7

logZ = −2.648 + 0.137 SST logZ = 0.211 − 0.343 Chl-a logZ = −1.133 + 0.097 V logZ = −1.526 + 0.100 SST − 0.216 Chl-a logZ = −3.142 + 0.194 SST − 0.075 V logZ = −0.175 − 0.323 Chl-a + 0.085 V logZ = f (SST,Chl-a,V )

2 2 2 3 3 3 4

28 28 28 27 27 27 26

0.335 0.284 0.144 0.423 0.363 0.393 0.423

14.08 11.08 4.696 9.891 7.680 8.733 6.35

<0.001

50.25 52.47 57.81 48.42 51.40 49.95 51.05

1.83 4,27 9,39 0 3,2 1,53 2,63

Table 5 Pearson’s correlation of bivariate relationships for the spawning area (log SA) and daily egg production (log P0 ) as a function of the Lloyd’s patchiness index of Group A (early development stages), Group B (middle development stages), and Group C (development stages close to hatching) (Period 2002–2014, n = 10). Variable

Zona

Group A (Stages I–IV)

Group B (Stages V–VIII)

Group C (Stages IX–XI)

Spawning area

North Central South

−0.907***

−0.808**

0.065 0.060

0.069 0.156

−0.416 0.385 0.507

0.672*

0.287 0.120 −0.199

0.497 −0.170 −0.292

Egg production

North Central South

−0.165 −0.519

*P < 0.05. **P < 0.01. ***P < 0.001.

would be a consequence of faster egg development rate than eggs in the Central and South spawning zones. In addition, the higher daily egg mortality rate in the North zone might be associated indirectly with egg traits such as volume (Llanos and Castro, 2004) or quality (e.g. biochemical composition, lipids) due to maternal influence (Castro et al., 2009). Indeed, Llanos and Castro (2004) found that anchovy eggs were smaller in the North zone than in the Central and South zones, probably due to warmer waters in the North. These intrinsic and environmental characteristics might explain the egg mortality differences that, along with the differences in development times,

0.002 0.039 <0.001 0.002 0.001 0.002

are affected also by temperature. In fact, sea surface temperature explained a significant fraction of the egg mortality variance, followed by phytoplankton biomass and alongshore wind as secondary co-variables explaining the egg mortality rate of anchovy across spawning areas. We expected a significant correlation between ENSO and Z since previous studies showed impacts of ENSO on the upwelling ecosystem and other traits of anchovy such as reproduction, growth and recruitment (Yáñez et al., 2008; Soto-Mendoza et al., 2010; Xu et al., 2015; Canales et al., 2019). The significant correlations between Z and environmental variables (sea surface temperature, phytoplankton biomass and alongshore wind) seem to be a pattern associated with the large latitudinal gradient rather than local effects. The pattern suggests that egg mortality is higher in the northern zone due to warmer waters, persistent upwelling, and higher phytoplanktonic biomass. The converse tends to occur in the southern zone. Also, higher egg mortality was associated with higher egg production and large spawning area. According to Claramunt et al. (2012), higher egg production and spawning area could be associated with smaller egg size and higher fecundity of anchovy in the North zone. Higher fecundity of anchovy in the North is a consequence of small egg size and quality (Castro et al., 2009), which would be linked with higher mortality rate. When nested-data effects associated with each spawning area were considered, correlations between Z and daily egg production and between Z and alongshore wind were significant only in the South zone. In this way, in the South zone, the anchovy spawning area seems to be configured spatially to other factors

8

B. Bustos, L.A. Cubillos, G. Claramunt et al. / Regional Studies in Marine Science 35 (2020) 101200

Fig. 5. Bivariate relationships between logarithm of the daily egg mortality rate (logZ ) and environmental variables within each of the study zones (North, Central and South zones): logZ vs. sea surface temperature (A), Z vs. Phytoplanktonic biomass (B), Z vs. alongshore wind (C).

such as local retention and/or concentration to favor local egg development and survival. Retention and concentration processes favor coastal aggregations of eggs (Canales and Leal, 2009, Rojas et al., 2011, Cubillos et al., 2013). Interestingly, the spawning area in the North zone showed a significant negative correlation with the patchiness index of early and middle development egg stages (Group A and B, respectively). Thus, while in the North zone factors favoring egg patchiness might be associated with high egg mortality, in the South zone spatial factors favoring egg concentration seem to be those associated with decreased egg mortality, through either less intense alongshore wind or northern alongshore wind. The overall results seem in agreement with observations on the contrasting distribution of eggs in different areas. The distribution of the spawning tends to occur offshore in the North zone (Claramunt et al., 2012), probably dispersed as a result of the constantly blowing upwelling-driver south winds. Instead, in the South zone, the distribution of the spawning tends to be coastal due to surface water convergence and almost inexistent upwelling during the spawning season (Castillo-Jordán et al., 2007; Cubillos et al., 2007b). Eggs are passive and distributed in the water column based on their buoyancy, and factors favoring either their concentration (fronts) or dispersion (transport) would determine their level of aggregation. These processes could be important for egg predation (McGurk, 1986; Krautz et al., 2003). According to the anchovy reproductive tactics of adults at each location, our results suggest that in the North zone, adults adopt a tactic that favors dispersion and select placing their eggs in oceanic waters to avoid predators (i.e., gelatinous macrozooplankton and euphausiids) (Pavez et al., 2006; Krautz et al., 2003). In the South zone, instead, female anchovy tends to place their eggs in coastal shallow areas where they concentrate and where euphausiids and the gelatinous macrozooplankton are absent in winter (Castro et al., 2000). Thus, the differences in the daily egg mortality of anchovy are associated with the differences in the reproductive tactics of adults, which in turn are conditioned also by maternal effects transferred to eggs (size and quality), and are tailored by local features among which aggregation or dispersion of eggs seem to play a major role. In the case of the anchovy Engraulis ringens along the Chilean coast, our results show that while the daily egg mortality might be dependent on dispersive processes in the North spawning area, in the South spawning area they are more associated with concentration and retention processes. Declaration of competing interest The authors declare that they have no known competing financial interests or personal relationships that could have appeared to influence the work reported in this paper. CRediT authorship contribution statement Blanca Bustos: Conceptualization, Methodology, Writing original draft, Data curation. Luis A. Cubillos: Visualization, Project administration, Writing - review & editing, Resources, Formal analysis. Gabriel Claramunt: Supervision, Resources. Leonardo R. Castro: Supervision, Resources. Acknowledgments

Fig. 6. Relationship between the spawning area and the patchiness index of early development embryonic stages of eggs (Group A) within each of the study zones: North, Central and South zones.

BB thanks COPAS South-Austral CONICYT PIA PFB31 for the scholarship to studying the Master in Fisheries, University of Concepcion. We thank FIPA (Fondo de Investigación Pesquera) and IFOP (Instituto de Fomento Pesquero) for providing the egg data here used. GC and LC were partially supported by FONDECYT, Chile 1161131. Thank you to referees, whose comments have helped us to improve this paper.

B. Bustos, L.A. Cubillos, G. Claramunt et al. / Regional Studies in Marine Science 35 (2020) 101200

References Alheit, J., Ñiquen, M., 2004. Regime shift in the Humboldt Current ecosystem. Prog. Oceanogr. 60, 201–222. Arcos, D., Núñez, S., Castro, L., Navarro, N., 1987. Variabilidad vertical de la clorofila en un área de surgencia frente a Chile central. Investig. Pesq. Valpso. 34, 47–55. Bakun, A., 1973. Coastal Upwelling Indices, West Coast of North America, 1946–71. US Department of Commerce, National Oceanic and Atmospheric Administration, National Marine Fisheries Service. Bernal, M., Borchers, D.L., Valdéz, L., Lanzós, A.L., Buckland, S.T., 2001. A new ageing method for eggs of fish species with daily spawning synchronicity. Can. J. Fish. Aquat. Sci. 58, 2330–2340. Bertrand, A., Chaigneau, A., Peraltilla, S., Ledesma, J., Graco, M., Monetti, F., Chávez, F., 2011. Oxygen: a fundamental property regulating pelagic ecosystem structure in the coastal southeastern tropical Pacific. PLoS One 6, e29558. Braun, M., Valenzuela, V., Claramunt, G., Reyes, H., Pizarro, M., Cataste, V., Herrera, G., Moreno, P., Gaspar, C., Díaz, E., 2005. Evaluación del stock desovante de anchoveta I y II regiones Año 2005. In: Informe FIP/2005-03. pp. 1–155. Burnham, K.P., Anderson, D.R., 2002. Model Selection and Multimodel Inference: A Practical Information-Theoretic Approach. Springer-Verlag. Bustos, C., Landaeta, M., Balbontín, F., 2008. Spawning and early nursery areas of anchoveta Engraulis ringens jenyns, 1842 in fjords of southern Chile. Rev. Biol. Mar. Oceanogr. 43 (2), 381–389. Canales, C.M., Adasme, N., Cubillos, L.A., Cuevas, M.J., Sánchez, N., 2019. Longtime spatio-temporal variations in anchovy (Engraulis ringens) biological traits off northern Chile: an adaptive response to long-term environmental change?. ICES J. Marine Sci. 75, 1908–1923. Canales, M., Leal, E., 2009. Parámetros de historia de vida de la anchoveta Engraulis ringens Jenyns, 1842, en la zona centro norte de Chile. Rev. Biol. Mar. Oceanogr. 44 (1), 173–179. Castillo-Jordán, C., Cubillos, L.A., Paramo, J., 2007. The spawning spatial structure of two co-occurring small pelagic fish off central southern Chile in 2005. Aquat. Living Resour. 10 (1), 77–84. http://dx.doi.org/10.1051/alr:2007018. Castro, L.R., Claramunt, G., Krautz, M.C., Llanos-Rivera, A., Moreno, P., 2009. Egg trait variations in anchoveta engraulis ringens: A maternal effect to changing environmental conditions in contrasting spawning habitats. Mar. Ecol. Prog. Ser. 381, 237–248. Castro, L.R., Hernandez, E., 2000. Early life survival of the anchoveta Engraulis ringens Off Central Chile during the 1995 and 1996 Winter Spawning Seasons. Trans. Amer. Fish. Soc. 129 (5), 1107–1117. Castro, L.R., Llanos, A., Blanco, J.L., Tarifeño, E., Escribano, R., Landaeta, M., 2002. Latitudinal variations in spawning habitat characteristics: influence on the early life history traits of the anchoveta, Engraulis ringens, off northern and central chile. In: GLOBEC Report 16. pp. 42–45. Castro, L.R., Quiñones, R., Arancibia, H., Figueroa, D., Roa, R., Sobarzo, M., Retamal, M., 1997. Áreas de desove de anchoveta y sardina común en la zona central. Fondo de Investigación Pesquera, Sub-Secretaria de Pesca, Informe final Proyecto FIP 96-11, Valparaíso, Chile. Castro, L., Salinas, G., Hernández, E., 2000. Environmental influences on winter spawning of the anchoveta Engraulis ringens off central Chile. Mar. Ecol. Prog. Ser. 197, 247–258. Chavez, F.P., Messié, M., 2009. A comparison of eastern boundary upwelling ecosystems. Prog. Oceanogr. 83, 80–96. Chavez, F.P., Ryan, J., Lluch-Cota, S.E., Ñiquen, M., 2003. From anchovies to sardines and back: Multidecadal change in the Pacific Ocean. Science 299, 217–221. Claramunt, G., Castro, L., Cubillos, L., Hirche, H., Perez, G., Braun, M., 2012. Inter-annual reproductive trait variation and spawning habitat preferences of Engraulis ringens off northern Chile. Rev. Biol. Mar. Oceanogr. 47 (2), 227–243. Claramunt, G., Cubillos, L., Braun, M., Serra, R., Canales, M., Sanchez, R., Flores, I., Moreno, G., Riquelme, K., Castillo, J., Valero, C., 2007. Mejoramiento del Método de Producción Diaria de Huevos aplicado en pelágicos pequeños. In: Informe FIP 2006-38. Claramunt, G., Cubillos, L.A., Castro, L., Hernández, C., Arteaga, M., 2014. Variation in the spawning periods of Engraulis ringens and Strangomera bentincki off the coasts of Chile: a quantitative analysis. Fish. Res. 160, 96–102. Correa-Ramirez, M., Hormazabal, S., Morales, C., 2012. Spatial patterns of annual and interannual Surface chlorophyll-a variability in the Peru-Chile system. Prog. Oceanogr. 92–95, 8–17. Cubillos, L.A., Castro, L., Claramunt, G., Navarro, E., 2013. Evaluación del stock desovante de anchoveta y sardina común entre la V y X Regiones. Año 2012. In: Informe FIP 2012-09. pp. 1–96. Cubillos, L.A., Castro, L., Claramunt, G., Navarro, E., 2015. Evaluación del stock desovante de anchoveta y sardina común entre la V y X Regiones. Año 2013. In: Informe FIP 2013-07.

9

Cubillos, L.A., Ruiz, P., Claramunt, G., Gacitúa, S., Núñez, S., Castro, L., Riquelme, K., Alarcón, C., Oyarzún, C., Sepúlveda, A., 2007b. Spawning, daily egg production, and spawning stock biomass estimation for common sardine (Strangomera bentincki) and anchovy (Engraulis ringens) off central southern Chile in 2002. Fish. Res. 86, 228–240. Cubillos, L.A., Serra, R., Fréon, P., 2007a. Synchronous pattern of fluctuation in three anchovy fisheries in the Humboldt Current System. Aquat. Living Resour. 20, 69–76. Essig, R.J., Cole, C.F., 1986. Methods of estimating larval fish mortality from daily increments in otoliths. Trans. Amer. Fish. Soc. 115, 34–40. Gómez, F., 2007. Variabilidad Ambiental Y Pequeños PeláGicos de la Zona Norte Y Centro-sur de Chile (Tesis Magíster en Ciencias mención Pesquerías). Universidad de Concepción, Concepción, p. 90. Hormazabal, S., Shaffer, G., Leth, O., 2004. Coastal transition zone off Chile. J. Geophys. Res.-Oceans 109, C01021. Jacox, M.G., Edwards, C.A., Hazen, E.L., Bograd, S.J., 2018. Coastal upwelling revisited: Ekman, Bakun, and improved upwelling indices for the U.S. West Coast. J. Geophys. Res.: Oceans 123. http://dx.doi.org/10.1029/2018JC014187. Kassambara, A., 2018. Ggpubr: ’ggplot2’ based publication ready plots. R package version 0.1.7. https://CRAN.R-project.org/package=ggpubr. Krautz, M.C., González, M., Castro, L.R., 2003. Detection of anchoveta (Engraulis ringens) eggs in euphausiid diets using immunoassays (ELYSA). J. Exp. Mar. Biol. Ecol. 294, 27–39. Llanos, A., Castro, L.R., 2004. Latitudinal and seasonal egg size variations of the anchoveta Engraulis ringens off the Chilean Coast. Fish. Bull. 102, 207–212. Llanos-Rivera, A., Castro, L.R., 2006. Inter-population differences in temperature effects on engraulis ringens yolk-sac larvae. Mar. Ecol. Prog. Ser. 312, 245–253. Lloyd, M., 1967. Mean crowding. J. Anim. Ecol. 36, 1–30. Lo, N.C.H., 1985. A model for temperature-dependent northern anchovy egg development. In: R., Lasker (Ed.), An Egg Production Method for Estimating Spawning Biomass of Pelagic Fish: Application To the Northern Anchovy, Engraulis Mordax. NOAA Tech Rep NMFS. 36, US Dep Commer, pp. 43–50. Mangiafico, S., 2019. Rcompanion: Functions to support extension education program evaluation. R package version 2.1.7. https://CRAN.R-project.org/ package=rcompanion. Marcotte, D., 1991. Cokrigeage with MATLAB. Comput. Geosci. 17, 1265–1280. McGurk, M., 1986. Natural mortality of marine pelagic fish eggs and larvae: the role of spatial patchiness. Mar. Ecol. Prog. Ser. 34, 227–242. McGurk, M., 1987. Natural mortality and spatial patchiness: reply to gulland. Mar. Ecol. Prog. Ser. 39, 201–206. Montecinos, A., Gómez, F., 2010. ENSO modulation of the upwelling season off southern-central Chile. Geophys. Res. Lett. 37, L02708 1-4. Morales, C.E., Braun, M., Reyes, H., Blanco, J.L., Davies, A.G., 1996. Anchovy larval distribution in the coastal zone off northern Chile: the effect of low dissolved oxygen concentrations and of a cold-warm sequence (1990-95). Investig. Mar. Valpso. 24, 77–96. Moser, H.G., Ahlstrom, E.H., 1985. Staging anchovy eggs. In: Lasker, R. (Ed.), An Egg Production Method for Estimating Spawning Biomass of Pelagic Fish: Application To the Northern Anchovy, Engraulis Mordax. NOAA Tech Rep NMFS. 36, US Dep Commer, pp. 37–41. Navarro, E., Schneider, W., Letelier, J., 2004. Estimation of onshore-o shore transport o central Chile by means of maximum cross-correlation using satellite derived SST. Gayana 68, 427–431. Ogle, D.H., Wheeler, P., Dinno, A., 2018. FSA: Fisheries stock analysis. r package version 0.8.22. https://github.com/droglenc/FSA. Parrish, R.H., Bakun, A., Husby, D.M., Nelson, C.S., 1983. Comparative climatology of select environmental processes in relation to eastern boundary currente pelagic fish reproduction. In: Sharp, G.D., Csirke, y.J. (Eds.), Proceedings of the Expert Consultation To Examine Changes in Abundance and Species Composition of Neritic Fish Resources. In: FAO Fishing Report 291, pp. 731–778. Pavez, M., Castro, L., Gonzalez, H.E., 2006. Across shelf predatory impact of pleurobrachia sp. (Ctenophora) on the small-copepods community in the coastal upwelling zone off mejillones, northern Chile: springs 2000–2002. J. Plankton Res. 28, 115–129. R Core Team, 2019. R: A Language and Environment for Statistical Computing. R Foundation for Statistical Computing, Vienna, Austria, URL https://www.Rproject.org/. Rojas, P., Landaeta, M., Ulloa, R., 2011. Eggs and larvae of anchoveta Engraulis ringens off northern Chile during the 1997-1998 El Niño event. Rev. Biol. Mar. Oceanogr. 46 (3), 405–419. Santander, H., Alheit, J., Smith, P.E., 1984. Estimación de la Biomasa de la Población Desovante de Anchoveta Peruana, Engraulis Ringens, en 1981 por Aplicación del ‘‘Método de Producción de Huevos’’, Vol. 8. Instituto del Mar del Perú, Callao, pp. 209–250. Serra, R., Aguayo, M., Rojas, O., Cañón, J., Inostroza, F., 1979. Anchoveta Engraulis ringens (Jenyns) teleostomi clupeiformes engraulidae. In: Estado Actual de Las Principales Pesquerías Nacionales. Bases Para Un Desarrollo Pesquero: I Peces. AP 79-18. En: CORFO-IFOP (eds.), pp. 1–52.

10

B. Bustos, L.A. Cubillos, G. Claramunt et al. / Regional Studies in Marine Science 35 (2020) 101200

Silva, C., Andrade, I., Yáñez, E., Hormazabal, H., Barbieri, M.A., Aranis, A., Böhm, G., 2016. Predicting habitat suitability and geographic distribution of anchovy (Engraulis ringens) due to climate change in the coastal areas off Chile. Prog. Oceanogr. 146, 159–174. Smith, P.E., Flerx, W., Hewitt, R., 1985. The CalCOFI vertical egg tow (CalVET) net. In: Lasker, R. (Ed.), An Egg Production Method for Estimating Spawning Biomass of Pelagic Fish: Application To the Northern Anchovy, Engraulis Mordax. NOAA Tech. Rep. NMFS 36, U.S. Dep. Commer, pp. 27–32. Smith, P.E., Hewitt, R., 1985. Sea survey design and analysis for an egg production method of anchovy biomass assessment. In: Lasker, R. (Ed.), An Egg Production Method for Estimating Spawning Biomass of Pelagic Fish: Application To the Northern Anchovy, Engraulis Mordax. NOAA Tech. Rep. NMFS 36, U.S. Dep. Commer, pp. 17–26. Sobarzo, M., Bravo, L., Donoso, D., Garcés-Vargas, J., Schneider, W., 2007. Coastal upwelling and seasonal cycles that influence the water column over the continental shelf off central Chile. Prog. Oceanogr. 75 (3), 363–382. Soto-Mendoza, S., Castro, L., Llanos-Rivera, A., 2010. Variabilidad espacial y temporal de huevos y larvas de Strangomera bentincki y Engraulis ringens, asociados a la desembocadura del río Itata, Chile. Rev. Biol. Mar. Oceanogr. 45 (3), 471–487.

Stratoudakis, Y., Bernal, M., Ganias, K., Uriarte, A., 2006. The daily egg production methods: recent advances, current applications and future challenges. Fish Fish. 7, 35–57. Valdivia, M., Chávez, R., Oliva, M., 2007. Metazoan parasites of engraulis ringens as tools for stock discrimination along the Chilean coast. J. Fish Biol. 70, 1504–1511. Venables, W.N., Ripley, B.D., 2002. Modern Applied Statistics with S, fourth ed. Springer-Verlag, NY. Wickham, H., 2016. Ggplot2: Elegant Graphics for Data Analysis. Springer-Verlag, New York. Xu, Y., Rose, K.A., Chai, F., Chavez, F., Ayón, P., 2015. Does spatial variation in environmental conditions affect recruitment?. A study using 3-D model of peruvian anchovy. Prog. Oceanogr. 138, 417–430. Yáñez, E., Hormazábal, S., Silva, C., Montecinos, A., Barbieri, M.A., Valdenegro, A., Órdenes, A., Gómez, F., 2008. Coupling between the environment and the pelagic resources exploited off northern Chile: ecosystem indicators and a conceptual model. Lat. Amer. J. Aquat. Res. 36 (2), 159–181. Zar, J.H., 2010. Biostatistical Analysis, fifth ed. Pearson Prentice Hall, Upper Saddle River, NJ.