A growing population on a planet with limited resources faces the challenge of food security through the sufficient supply of mankind with carbohydrates, fats and proteins (Prosekov and Ivanova, 2018). Especially with regard to proteins, fish is of key importance: for example, small epipelagic fish are used for the production of fishmeal and fish oil which are used to feed aquaculture animals, land stocks and to produce nutritional capsules for human consumption (Shepherd and Jackson, 2013). The averaged global fishmeal and fish oil production between 2001 and 2006 was 6.3 and 0.95 Mt per year, respectively (Péron et al., 2010). From these, 1.7 and 0.27 Mt came from small pelagic fish landed in Peru (Péron et al., 2010). Small epipelagic fish in the Northern Humboldt Current system (NHCS) –located in the Eastern Tropical South Pacific Ocean (ETSP)– represent 10% of the global fish landings (Chavez et al., 2008). However, the exploitation potential of these coastal stocks is limited (Tarazona and Arntz, 2001) and they have collapsed in the past due to overfishing and recruitment failure, impacting the ecosystem (Duffy, 1983; Tarazona and Arntz, 2001; Herling et al., 2005). Their susceptibility to high temporal variability in environmental conditions such as temperature (Chavez et al., 2003) and oxygen (Bertrand et al., 2011), in combination with the possible impacts of climate change on the NHCS, bring uncertainty for their exploitation in the upcoming decades (see Salvatteci et al., 2022). Alternative fish stocks may be necessary in the coming years to satisfy the demand for fishmeal and release the pressure on currently over-exploited epipelagic fish. Hence, fisheries are turning their perspectives towards mesopelagic fish, which have, so far, remained relatively unexploited (St. John et al., 2016). These may be used to support the supply of fishmeal and also as source for nutraceutical products (St. John et al., 2016). However, exploiting these resources without prior knowledge on their fundamental ecological and biogeochemical role poses threats for the mesopelagic community and, potentially, also for the ocean health and global climate (St. John et al., 2016; Martin et al., 2020).

Estimates regarding total standing stock of the mesopelagic community are scarce and uncertain (Davison et al., 2013; Belcher et al., 2019). While earlier estimates of the global mesopelagic biomass were around 1 Gt wet weight (Gjøsæter and Kawaguchi, 1980), a more recent estimate based on large scale echosounder surveys is an order of magnitude higher (11-25 Gt; Irigoien et al., 2014); yet, the latter high estimate is currently under debate (Davison et al., 2015b; Anderson et al., 2019). A different approach was followed by Anderson et al. (2019) in a flow based steady state analysis of the mesopelagic food web. They derived an estimate of 2.4 Gt which is significantly lower compared to recent estimates by Irigoien et al. (2014). Although the model by Anderson et al. (2019) derives the fluxes from observations, the food web they consider includes only the basic fluxes and neglects more complex interactions that could potentially also have an impact on the resulting biomass such as predation over different life stages. In the NHCS, Vinciguerria lucetia, also known as Panama lightfish, and myctophids, commonly known as lanternfish, have been reported as the main constituents of the mesopelagic fish community (Cornejo Urbina and Koppelmann, 2006; Marzloff et al., 2009). These are vertical migrants whose distribution has been reported on the upper 50 m depth during the night and between 200 and 400 m during the day (Cornejo Urbina and Koppelmann, 2006). Between 2.9 (Castillo Valderrama et al., 1999) and 11.1 Mt (Castillo Valderrama et al., 1998) of Vinciguerria sp. have been estimated for the upwelling region off Peru. Hence, even for this well-studied region, we find a large uncertainty for the biomass estimate of this group.

Ecosystem and fisheries models are valuable tools to understand the dynamics of the ecosystems and their potential response under certain scenarios. Numerous models exist to simulate either single fisheries or whole ecosystems (see Fulton, 2010; Tittensor et al., 2018; Tittensor et al., 2021). Ecosystem models tend to focus on commercially exploited fisheries, especially regarding the epipelagic community and demersal species (e.g., Rose et al., 2015; Watson et al., 2015; Carozza et al., 2016; Petrik et al., 2019). In recent years, the focus has started to shift from a mainly fisheries management perspective (see Colléter et al., 2015) to also consider the impact of fish and fisheries on biogeochemistry (e.g., Megrey et al., 2007; Travers-Trolet et al., 2014b; Getzlaff and Oschlies, 2017; Aumont et al., 2018; Bianchi et al., 2021). Only few models also consider mesopelagic organisms, aggregating them as a functional group (e.g., Ainsworth et al., 2015; Aumont et al., 2018; Anderson et al., 2019) or simulating individual species (e.g., Travers-Trolet et al., 2014b). A sustainable exploitation of mesopelagic fish requires not only to estimate its biomass but also to understand the population’s vital rates such as recruitment and growth (St. John et al., 2016) and their ecological and biogeochemical role (Martin et al., 2020). Hence, understanding the benefits and limitations of current approaches to model mesopelagic fish is essential to identify opportunities and requirements for the development of mesopelagic fish models.

In this study, we provide a new estimate of mesopelagic fish biomass in the ETSP, a region that is relatively well studied with regard to observations of fish biomass. To figure out potential sources of the large uncertainties in biomass estimates, we employ two different model types for higher trophic levels, which are both forced by the same physical–biogeochemical model. We explore the strengths and weaknesses of the two models and how their different complexities affect the estimates of mesopelagic fish biomass and their temporal variability. This is a first step towards understanding the trophic web of the mesopelagic zone and how it can be modelled. To do so, we coupled the regional physical–biogeochemical model CROCO–BioEBUS (Coastal and Regional Ocean COmmunity model, Shchepetkin and McWilliams, 2005; Biogeochemical model for Eastern Boundary Upwelling Systems, Gutknecht et al., 2013) with the two ecosystem models to simulate biomass of mesopleagic fish in the ETSP. The multispecies model OSMOSE (Object-oriented Simulator of Marine Ecosystems, Shin and Cury, 2004) is a dynamic individual based model that simulates the whole life cycle and trophic interactions of nine fish species, functional groups –such as the mesopelagic fish– and macroinvertebrates. OSMOSE receives integrated plankton estimates by CROCO–BioEBUS as food forcing for the fish. The second mesopelagic fish model that we used is the flow analysis model by (Anderson et al 2019, called A2019 in the following). This model takes primary production as input which, in the case of our study, is provided by CROCO–BioEBUS. The model calculates the pathways taken by the biomass, either sinking to the deep water or being consumed by zooplankton, until it reaches mesopelagic fish. The model provides mesopelagic fish biomass as output. A detailed description of the models employed in this study is provided in Section 2. In Section 3, we present the spatial and temporal trends in the mesopelagic fish estimated by OSMOSE and the model by A2019. We also explore how an explicit trophic web and life cycle in OSMOSE may impact its response to temporal variability in the plankton forcing in contrast to the model by A2019. Section 4 discusses differences and strengths of the two models based on our results and how they could benefit each other, followed by the conclusions of the study.

In this study, we calculated the abundance of mesopelagic fish in the Eastern Tropical South Pacific (ETSP) by using the physical–biogeochemical model CROCO–BioEBUS (Shchepetkin and McWilliams, 2005; Gutknecht et al., 2013) coupled to two different models of higher trophic levels (from now on called fish models): the simple food-web model by A2019 and the multispecies individual-based model OSMOSE (Shin and Cury, 2001; Shin and Cury, 2004). Primary production, phytoplankton and zooplankton were estimated using the physical–biogeochemical model and these in turn were used to force the two fish models.

The spatial distribution of the mesopelagic fish simulated with OSMOSE and with the model by A2019 using BioEBUS primary production is similar ( Figure 1 ). The mesopelagic fish are absent in the coastal shallow water, which, for the model by A2019, was masked from the primary production input (Section 4.3), and their largest concentration is present off Peru. Independent of the choice of the mortality value, the model by A2019 produces a larger biomass than OSMOSE ( Figure 1 top). Applying the larger mortality of 3.02 y^–1 to the model by A2019 ( Figure 1 top middle), the resulting total biomass is about 26% larger than the biomass in OSMOSE, whereas applying the lower mortality of 1.65 y^–1 ( Figure 1 top right) results in a total biomass that is about 130% larger than in OSMOSE (0.0208, 0.0262 and 0.0476 Gt WW of biomass estimated by OSMOSE, A2019 with high mortality and A2019 with low mortality respectively, Figure 1 ). The results indicate that the representation of different life stages of mesopelagic fish is likely to impact the resulting biomass. In our case, the comparatively large mortality of the juvenile fish in OSMOSE seems to play an important role and needs to be taken into account by the flux representations of the model by A2019. In the following we analyze the simulation with the higher mortality only (3.02 y^–1), which is closer to the output from OSMOSE.

The seasonal cycle in both higher trophic level models has the same frequency, but the amplitude of the mesopelagic fish biomass is much larger in the model by A2019 compared to the one in OSMOSE ( Figure 1 bottom). In the model by A2019, mesopelagic fish biomass is a linear function of primary production, hence the anomaly of the interannual variability of the mesopelagic fish biomass in the model by A2019 is identical to the anomaly of the interannual variability of BioEBUS’ primary production ( Figure 2 ). In contrast to that, we find that the amplitude of the seasonal cycle in OSMOSE is muted whereas the amplitude of its interannual variability is larger compared to the one in the model by A2019 and does not follow the same trend as the plankton input ( Figure 2 ). In addition, the anomaly in mesopelagic fish biomass simulated by OSMOSE does not seem to follow the trend observed in BioEBUS’ plankton biomass. The decoupling is indicated by a weak cross correlation of maximal 0.197 (see Appendix).

While the temporal variability in the model by A2019 is solely driven by the primary production, in OSMOSE it originates from the interplay of many other factors. These include the explicit life cycle of the fish, the trophic interactions with predators and different sources of prey, the spatial distribution of the fish and the variability in different sources of mortality. The main predator of mesopelagic fish in OSMOSE is the Humboldt squid. The main prey of mesopelagic fish are euphausiids and then large zooplankton. Removing all fish groups from the OSMOSE configuration except for the mesopelagics ( Figure 2 , pink line) resulted in a decrease in the amplitude of the interannual anomaly of mesopelagic fish biomass. In this simulation, there is a maximum cross correlation of 0.65 between the 12-month running means of plankton forcing and mesopelagic biomass with a lag of 8 months (see Appendix). The lag indicates that the fish biomass responds to changes in plankton forcing after about 8 months. This points to an impact of simulating a complex life cycle in OSMOSE, including egg production, by delaying the biomass response to plankton abundance. In summary, both the model by A2019 and OSMOSE show a similar seasonal cycle following the trend of the plankton and primary production. However, the interannual variability in OSMOSE is stronger and follows a different pattern even if mesopelagic fish are the only group present in this model.

In our study, we observed a higher biomass in the model by A2019 than in OSMOSE of about 26 to 130%, depending on the mortality used. In general, mortality in numerical models is a highly uncertain parameter. As seen in this study, depending on the underlying assumptions, anticipate the mortality of juveniles or that of adult, the resulting biomass of mesopelagic fish can differ substantially.

OSMOSE model was calibrated to match observed biomasses of mesopelagic fish off Peru between 2000 and 2008 provided by the Instituto del Mar del Peru (IMARPE) as described in Hill Cruz et al. (2022). Nonetheless, acoustic estimates of mesopelagic fish are prone to high uncertainty (Marzloff et al., 2009; Davison et al., 2015a; Davison et al., 2015b). In addition, observations by IMARPE are done within the Exclusive Economic Zone of Peru (maximum distance to the coast of 100 to 200 nautical miles, see Castillo Valderrama and Gutiérrez Torero, 2001; Castillo Valderrama et al., 2004; Castillo Valderrama et al., 2009a; Castillo Valderrama et al., 2009b; Castillo Valderrama et al., 2009c). Castillo Valderrama et al. (1999) suggest that the abundance of Vinciguerria lucetia might be underestimated in a case where the transects were done closer to the shore since this is an oceanic species. Hence, the somewhat high mortality of mesopelagic fish could stem from the low bias of observations against which OSMOSE was calibrated.

Determining mortality rate M , and its sources, of fish species in the ocean is not straightforward. A2019 derived the value of M = 0.67 y^–1 by calculating the inverse of a 1.5 year longevity ( t max ). While using longevity to estimate natural mortality is a convenient method when limited information is available, it might be an oversimplification and has been discouraged when more advanced methods are available (Kenchington, 2014). Often, longevity is determined from otoliths of caught fish (Takagi et al., 2006; Hosseini-Shekarabi et al., 2015). This reflects the age of the fish that survived until they were captured and provide information on how long a fish may live under certain circumstances, but does not distinguish between different processes such as predation, disease or starvation. There are a number of different methods to determine the natural mortality that involve longevity (see Kenchington, 2014; Then et al., 2014). The so called rule of thumb is also a function of longevity and estimates mortality as M ≈ 3 / t max (Hewitt and Hoenig, 2005; Maunder and Wong, 2011). Applying this rule of thumb to mesopelagic fish, assuming a longevity of 1.5 years as in A2019, would result in a mortality of 2 y^–1, which is larger compared to the value of 0.67 y^–1 used in A2019. Then et al. (2014) derived a different relation between mortality and longevity based on information from more than 200 fish species: M = 4.899 t m a x − 0.916 that yields a mortality of 3.38 y^–1 for our particular case. This is closer to our 3.02 y^–1 value estimated by OSMOSE for juvenile fish. Other approaches provide even higher values. For instance, the estimator by Hoenig (1983) for fish ( M = 6.99 t m a x − 1.22 , see Kenchington, 2014) applied to mesopelagic fish results in a mortality of 4.26 y^–1. In summary, our results suggest that the mortality estimated by OSMOSE is within the range of mortality estimates obtained with other methods based on longevity. Therefore, mortality rates diagnosed from OSMOSE can support the parameterization of less complex food web models since its higher complexity includes explicit and dynamic sources of mortality including starvation and predation for different life stages of the fish.

Both OSMOSE and the model by A2019 offer advantages and disadvantages. The model by A2019 represents the mass transfer from nutrients over plankton to mesopelagic fish in a strainghtforward way. Its simplicity allows us to carry out many sensitivity experiments that address the effects of a changing biogeochemistry on fish.

The comprehensive parameter uncertainty analysis carried out by A2019 provided low and high boundaries for the integrated biomass estimates between 40°N and 40°S. Such an approach becomes more complicated with OSMOSE since running this individual-based model requires higher computational resources; also the model species may collapse under certain conditions (see Supplement to Hill Cruz et al., 2022). Therefore, for unconstrained regions, especially in the high seas, the model by A2019 can be used for providing first estimates of potential mesopelagic fish biomass. Nevertheless, the model by A2019 does not explicitly simulate life cycles or a trophic web but only represents the (linear) biomass transfer from primary production to mesopelagic fish as a steady state. This might be an oversimplification for a highly dynamic region such as the ETSP where environmental variability has been seen to have strong impacts on small pelagic fish (Chavez et al., 2003; Alheit and Niquen, 2004). A dynamic model is necessary to study potential responses of this type in mesopelagic fish and to estimate a population growth rate and recruitment based on life traits.

The individual-based model OSMOSE represents the whole life cycle of several fish and macroinvertebrate species. Since this is a dynamic model, it is a good alternative for representing non-linear changes in fish populations triggered by different pressures such as environmental variability. Our results show that the non-linearity affects how mesopelagic fish respond to variability in their food, which is produced by environmental changes. As a multispecies model, OSMOSE is also useful for exploring fishing strategies in an ecosystem-based fisheries management context (Briton et al., 2019; Fu et al., 2019; Guo et al., 2019; Fu et al., 2020). There is growing evidence of the importance of mesopelagic fish in deep water food chains (Mann, 1984; Davison et al., 2015a; Saunders et al., 2019). In our study, the main predator of mesopelagic fish is the Humboldt squid (see Appendix). This is a species of economic importance (Gilly et al., 2013) that feeds mainly on mesopelagic fish (Markaida and Sosa-Nishizaki, 2003). Therefore, any prospect on the exploitation of mesopelagic fish in the ETSP should consider the potential impacts on this species. Likewise, changes in the exploitation on the Humboldt squid should consider potential cascading effects on the mesopelagic fish and lower trophic levels. Furthermore, while the spatial distribution offshore in the model by A2019 was fully dependent on the primary production, in OSMOSE the distribution of fish schools is forced by habitat niche models. These maps vary in time and can be used to simulate horizontal migrations (see Grüss et al., 2014; Oliveros-Ramos, 2014) allowing for further refinement in the distribution of fish when this information is available. Nevertheless, our results show a similar spatial pattern in both OSMOSE and the model by A2019. This indicates that the food source has a strong impact on the distribution of mesopelagic fish. Hence, the approach used in the model by A2019 can provide a first idea of the distribution of mesopelagic fish where observations and habitat niche models are not available, for instance in a global setting.

Both OSMOSE and the model by A2019 lack a vertical dimension of fish distribution, but parameterize its effect on the trophic web. In the case of OSMOSE, this is done through a predatory accessibility matrix which constrains which fish groups can feed on others based on their hypothetical position in the water column. The model by A2019 goes a step further by also parameterizing vertical migrations of mesopelagic fish. It assumes in its trophic chain that migrant fish feed on epipelagic zooplankton for half of the day and on migrating zooplankton all day [A2019]. Vertical migrations are especially important in the context of mesopelagic fish due to their role in the active transport of organic matter to the deep ocean (Davison et al., 2013; Belcher et al., 2019; Hernández-León et al., 2019) and its implications for carbon capture, oxygen loss and nutrient cycles. However, vertical migrations are only implicitly represented in the model by A2019. A full coupling of the model by A2019 to the biogeochemical model would be necessary to capture the impact of vertical migrations on biogeochemical cycles.

Depending on the mortality rates applied, the biomass estimated by the model by A2019 was more or less (26% with high mortality and 130% with low mortality) higher than in OSMOSE. In both cases, the mortalities calculated by OSMOSE were larger than the original mortality estimated as the inverse of longevity by A2019. Although OSMOSE has been calibrated to match the limited amount of observations available for the coastal region of the ETSP, uncertainties are large with respect to the open ocean where mesopelagic fish are abundant. To reduce this uncertainty, more reliable estimates of mesopelagic fish biomass in the open ocean would be beneficial. OSMOSE exhibits a weaker seasonal cycle but a stronger interannual variability. On the contrary, the model by A2019 follows the same temporal pattern as the primary production forcing due to its linear nature. Looking ahead, both models have strengths and disadvantages. OSMOSE provides a richer representation of the life cycle and trophic interactions of the fish community. This is crucial for studies steered towards implementations of fisheries management. The strength of the model by A2019 lies in the direct, computationally efficient link between the biogeochemistry and the fish representation, which is useful for understanding potential boundaries for the available biomass of mesopelagic fish given certain biogeochemical conditions. The most comprehensive representation of the ecosystem could combine the strengths of both models. Being the most simple and computationally efficient approach, the model by A2019 could provide a first estimate for regions where data are limited, such as in the global ocean and the high seas. OSMOSE simulates a more comprehensive food web and is the preferred approach when observational data are available to calibrate the model. With its more detailed resolution of food web and ecological processes, it can also aid the validation of parameters in simpler models that do not have a direct observational counterpart such as mortality. In addition, model results can be assessed against many other quantities such as stomach content, trophic level, size structure and catches to provide a more accurate representation of the ecosystem.

The raw data supporting the conclusions of this article will be made available by the authors, without undue reservation.

MH and IK designed the experiments. IK set up the model by A2019. MH set up OSMOSE, carried out the simulations with both fish models coupled to the physical–biogeochemical model and analyzed the data. All authors discussed the results and wrote the manuscript. All authors contributed to the article and approved the submitted version.