# Reliable machine learning initialization methods for the calibration of Dynamic Energy Budget models

Diogo Oliveira, Gonçalo Marques, Filipe M.P. Santos, Laure Pecquerie, João M.C. Sousa, Tiago Domingos

► To cite this version:

Diogo Oliveira, Gonçalo Marques, Filipe M.P. Santos, Laure Pecquerie, João M.C. Sousa, et al.  
Reliable machine learning initialization methods for the calibration of Dynamic Energy Budget models.  
Ecological Informatics, 2026, 94, pp.103624. 10.1016/j.ecoinf.2026.103624. hal-05514637

HAL Id: hal-05514637

https://hal.science/hal-05514637v1

Submitted on 19 Feb 2026

HAL is a multi-disciplinary open access
archive for the deposit and dissemination of sci-
entific research documents, whether they are pub-
lished or not. The documents may come from
teaching and research institutions in France or
abroad, or from public or private research centers.

L'archive ouverte pluridisciplinaire HAL, est
destinée au dépôt et à la diffusion de documents
scientifiques de niveau recherche, publiés ou non,
émanant des établissements d'enseignement et de
recherche français ou étrangers, des laboratoires
publics ou privés.

<figure><img src="image_1.png" /></figure>

Distributed under a Creative Commons CC BY 4.0 - Attribution - International License
<figure><img src="image_2.png" /></figure>

# Reliable machine learning initialization methods for the calibration of Dynamic Energy Budget models

Diogo F. Oliveira a,b,c,*, Gonçalo M. Marques a,b, Filipe M.P. Santos c, Laure Pecquerie d,e,
João M.C. Sousa c, Tiago Domingos a,b

## ARTICLE INFO

**Keywords:**

Dynamic energy budget

Model calibration

Parametrization

Automated initialization

Neural networks

Nearest neighbors

## ABSTRACT

Dynamic Energy Budget (DEB) theory is a general theory that describes how organisms utilize the energy in food for maintenance, growth, development, and reproduction. DEB models have been widely applied in fields such as conservation biology, aquaculture and ecotoxicology, due to their ability to simulate how organisms respond to changing environmental conditions. To obtain a DEB model, the calibration problem must be solved: find the parameters that minimize the deviation between observed data and model predictions. While DEB model calibration is largely automated, the selection of initial parameters remains a key unresolved step, since the only automated method – the bijection method – often fails to produce a feasible initial parameter set. Consequently, modelers resort to trial-and-error to find parameters to seed the estimation. To bridge this gap, we propose using machine learning to initialize the calibration. We develop two models: a neural network and a 1-nearest-neighbor. Both models are built with a focus on feasibility, directly integrating parameter constraints into their structure. We train and evaluate our methods on the 5000+ DEB models in the Add-my-Pet database. Both methods generate feasible parameter sets in 99% of cases — compared to only 40% for the bijection method. The neural network initialization leads to improved DEB model calibration, achieving a calibration loss three times lower, on average, when compared to other methods. To support broader adoption, we have open-sourced our code and our models are available as initialization options within DEBtool, the primary software for parameter calibration.

## 1. Introduction

Mathematical modeling plays a fundamental role in biology, offering insights into complex systems that are otherwise difficult to study (Gunawardena, 2014; Kooijman et al., 2024). Biological models enable researchers to predict growth patterns, resource utilization, and the effects of environmental changes on organisms (Edelstein-Keshet, 2005; Grimm and Railsback, 2013). Such predictive capabilities are essential for fields like conservation, agriculture, and ecosystem management, where understanding biological responses can inform decision-making and policy (Elsawah et al., 2020; Jakeman et al., 2024).

Dynamic Energy Budget (DEB) theory provides a quantitative framework for understanding the metabolism of organisms across their life cycle (Kooijman, 2010; Sousa et al., 2008). Through a state-based model, DEB theory describes how organisms allocate the energy in food to growth, development, maintenance, and reproduction. DEB models respect the principles of conservation of mass and energy and can be used to predict the impacts of changes in temperature, food availability, and other environmental factors or stressors (Sousa et al., 2006). DEB theory has been applied in various fields such as conservation biology (Molnár et al., 2010; Haberle et al., 2020; Lavaud et al., 2021; Martins et al., 2024; Kim et al., 2019), ecotoxicology (Jager et al., 2006; Desforges et al., 2017; Gergs et al., 2021; Jager et al., 2023),

<footnote>ª Terraprima - Serviços Ambientais, Centro de Negócios do Porto Alto, Fracção S, Avenida das Nações Unidas 97, 2135-199 Samora Correia, Portugal</footnote>

<footnote>ᵇ MARETEC - Marine, Environment and Technology Centre, LARSyS, Instituto Superior Técnico, Universidade de Lisboa, Avenida Rovisco Pais 1, 1049-001 Lisboa, Portugal</footnote>

<footnote>ᶜ IDMEC, Instituto Superior Técnico, Universidade de Lisboa, Avenida Rovisco Pais 1, 1049-001 Lisboa, Portugal</footnote>

<footnote>ᵈ LEMAR, Institut de Recherche pour le Développement (IRD), UBO/IRD/CNRS/Ifremer, IUEM, rue Dumont d'Urville, 29280 Plouzané, France</footnote>

<footnote>ᵉ ISI, Computer Science Depart., FSSM, Cadi Ayyad University, Bd. Prince My Abdellah, 40000 Marrakesh, Morocco</footnote>

<footnote>* Corresponding author at: MARETEC - Marine, Environment and Technology Centre, LARSyS, Instituto Superior Técnico, Universidade de Lisboa, Avenida Rovisco Pais 1, 1049-001 Lisboa, Portugal.

*E-mail address:* diogo.miguel.oliveira@tecnico.ulisboa.pt (D.F. Oliveira).</footnote>
<figure><img src="image_3.png" /><figcaption>**Fig. 1.** DEB model calibration with machine learning initialization. A schematic of our proposed approach to solve the DEB model calibration problem. The machine learning models combine the species data with taxonomic classifications and habitat information to produce an initial guess of the parameters. This guess is first checked for feasibility and then is refined in the DEB model calibration procedure.</figcaption></figure>

aquaculture (Pouvreau et al., 2006; Giacoletti et al., 2024; Tan et al., 2022; Palmer et al., 2024), and livestock production (Marques et al., 2020; Oliveira et al., 2024). See the DEB Library<fnref n="1" /> for a list of more than 1200 publications that apply DEB theory.

The Add-my-Pet collection, an online repository of DEB models, currently houses models for more than 5500 species<fnref n="2" /> from all large animal phyla (Marques et al., 2018; Kooijman et al., 2024). For these models to be useful, they must be well calibrated with accurate parameter values. The calibration problem consists of finding the set of DEB parameters that minimize the error between the available data and model predictions, which is solved automatically by an optimization algorithm. In the end, the parameters are compared with those of similar species to verify consistency, a procedure called estimation-in-context, which follows from the principle that similar species should have similar parameter values (Marques et al., 2019; Kooijman et al., 2024). The available data is obtained from experiments, sampling and literature review, and therefore can differ significantly from species to species, both in characteristics and quality (Marques et al., 2019).

To start the calibration, an initial set of parameters must be chosen. Because DEB theory imposes constraints on parameter values to ensure that the model is realistic, the initial parameter set must respect several constraints, some of which require solving differential equations. If the parameter set is feasible but close to infeasibility, the estimation may halt due to the inability of the differential equation solver to converge. In addition, a poor initialization will lead to longer convergence times and may stop at a local minima with worse model fit (Li et al., 2020; Friedrich and Wagner, 2015).

Most modelers opt to choose parameters from a taxonomically similar species, following the estimation-in-context principle. The only automated method is the bijection method, an algorithm that predicts parameters from a few commonly available data types (Lika et al., 2014). Although it has not been previously reported, our experience with the bijection method and the results of this paper show that it often fails to provide feasible initial parameter values. The lack of a dedicated and reliable initialization method often leads to a trial-and-error procedure to get a feasible initial values, which is not trivial, particularly for modelers with less experience. Furthermore, it increases the likelihood that one of the aforementioned problems will occur.

To address this, we propose using machine learning (ML) to initialize the calibration (Fig. 1). The rationale for this approach is that machine learning models learn from examples, so our ML models learn to incorporate the estimation-in-context principle. Moreover, machine learning methods can draw on related sources of information during calibration in addition to species data. For instance, taxonomy

classifications and habitat information can be linked to species parameter values (Lika et al., 2011b; Marques et al., 2018) and therefore contribute to more accurate calibration.

We develop two types of models: neural networks and nearest-neighbors. The neural networks are inspired by the theoretical framework of the bijection method, taking as input several specific data types and combining them with taxonomic and habitat data. The nearest-neighbor models mimic the *ad hoc* approach of choosing a similar species. The model uses a custom distance metric based on taxonomy and ultimate weight. Both models are designed to produce feasible parameters by directly integrating most feasibility constraints into their structure.

We leverage the available data and parameter estimates in the Add-my-Pet collection for more than 2000 species to train and evaluate the proposed models. We start by verifying how often the methods produce feasible parameter estimates. Then we compare how the predictions of the ML models are to the DEB model parameters. Finally, we run the DEB model calibration with the initial parameter values produced by our models and assess the performance of each initialization method.

The remainder of this paper is structured as follows. Section 2 provides a brief introduction to DEB theory, the DEB model calibration problem, and how to reframe it as a machine learning problem. Section 3 addresses data collection and the techniques and decisions that led to the final dataset used. Section 4 describes the ML models and how to incorporate feasibility constraints into their structure. Section 5 describes the experimental setup used to train and evaluate our models. In Section 6, we present and discuss the results of the evaluation and compare the proposed approaches with the bijection method. Finally, in Section 7, we provide concluding remarks, examine limitations, and outline future work.

## 2. Calibration of dynamic energy budget models

### 2.1. *Dynamic energy budget theory*

Dynamic Energy Budget theory aims to model the metabolism of an organism. DEB theory takes a simplified view of metabolism, but introduces enough complexity to describe how a changing environment, in terms of food and temperature, affects each key metabolic process throughout the life cycle of an organism (Kooijman, 2020; Sousa et al., 2006). DEB models typically consider four state variables, reserve *E*, which accounts for energy reserves, structure *V*, which quantifies growth, maturity *E*<sub>H</sub>, which represents development, and reproduction buffer *E*<sub>R</sub>, the energy reserves used for reproduction. The basic assumptions of DEB theory apply to all organisms, meaning that the difference between the DEB models of two species is mainly in the parameter values, which are presented in Table 1. To avoid unnecessary detail, we briefly describe the standard DEB model and only introduce the equations in Appendix A. A scheme of the model is also shown in Fig. 2.

<footnote n="1">¹ https://www.zotero.org/groups/500643/deb_library/</footnote>

<footnote n="2">² https://www.bio.vu.nl/thb/deb/deblab/add_my.pet/species_list.html, Accessed: Apr. 10, 2025.</footnote>
<figure><img src="image_4.png" /><figcaption>Fig. 2. A scheme of the standard DEB model.</figcaption></figure>

When food is ingested, only part of its mass and energy is assimilated by the organism and converted into reserve *E*. The remainder is wasted and expelled to the environment in the form of feces *P*, nitrogen waste, carbon dioxide, water vapor, and heat. Part of the reserve is stored but most of it is utilized to pay energetic costs for other metabolic processes such as growth, maintenance, development, and reproduction. Growth corresponds to an increase in structure *V*. Structure is not directly impacted by fluctuations in food availability and requires maintenance, contrary to reserve. Development is the increase in maturity *E*<sub>H</sub>, which represents the cumulative energy costs associated with the increasing complexity of the organism throughout ontogeny. Maintenance has a somatic component, which is a function of structure *V*, and a maturity component, which is proportional to maturity *E*<sub>H</sub>.

Thresholds in maturity values define the behavior of the model in the various life stages of the organism. Until birth the organism does not feed, and before puberty it does not reproduce. When it reaches puberty, maturity stops increasing, and the energy is instead used for reproduction, for example, to produce gametes or to provide energy for embryos.

The energy intake of food $\dot{p}_X$ is often expressed as a function of food availability $X$ via the scaled functional response $f$, which can take values between 0 and 1, where 1 means abundant food. Changes in food availability directly affect how much energy is available to the organism and how much is allocated to each process. For example, growth and development occur only if all maintenance costs have been paid. Changes in temperature affect the speed of each metabolic process by modifying the values of parameters that have time units. For more details on how food availability and temperature are modeled, see Appendix A.

Some organisms require changes to the standard DEB model (std) described above in order to properly model key processes that are life-stage-specific. These changes are encompassed into typified DEB models which can be readily applied to most species (Marques et al., 2018). Here, we briefly describe the abj model which is the most common typified DEB model in the Add-my-Pet collection.

The abj model is a one-parameter extension (maturity at metamorphosis $E_H^j$) and is typically applied to fish. These organisms undergo metabolic acceleration from birth until metamorphosis ($E_H^b < E_H^j < E_H^H$), during which growth is exponential (for a constant scaled functional response). Afterwards, the abj model is equal to the standard model, but the parameters {$\dot{p}_{Am}$} and $\dot{v}$ are multiplied by the acceleration factor $s_M = L_j/L_b$, where $L_j$ and $L_b$ are the structural length ($L = V^{1/3}$) at metamorphosis and birth, respectively.

<table><caption>**Table 1**
Core parameters of the standard DEB model.</caption><thead><tr><th>Symbol</th><th>Parameter</th><th>Default value</th><th>Unit</th></tr></thead><tbody><tr><td>{p&#775;<sub>Am</sub>}</td><td>Surface-specific maximum assimilation rate</td><td>—</td><td>J d<sup>-1</sup> cm<sup>-2</sup></td></tr><tr><td>&kappa;<sub>X</sub></td><td>Digestion efficiency</td><td>0.8</td><td>—</td></tr><tr><td>&kappa;<sub>P</sub></td><td>Defecation efficiency</td><td>0.1</td><td>—</td></tr><tr><td>&#775;<sub>v</sub></td><td>Energy conductance</td><td>0.02</td><td>cm d<sup>-1</sup></td></tr><tr><td>&kappa;</td><td>Allocation fraction to soma</td><td>0.8</td><td>—</td></tr><tr><td>[p&#775;<sub>M</sub>]</td><td>Volume-specific somatic maintenance rate</td><td>18</td><td>J d<sup>-1</sup> cm<sup>-3</sup></td></tr><tr><td>[E<sub>G</sub>]</td><td>Specific cost for structure</td><td>2800</td><td>J cm<sup>-3</sup></td></tr><tr><td>k&#775;<sub>j</sub></td><td>Maturity maintenance coefficient</td><td>0.002</td><td>d<sup>-1</sup></td></tr><tr><td>E<sub>H</sub><sup>j</sup></td><td>Maturity at birth</td><td>—</td><td>J</td></tr><tr><td>E<sub>H</sub><sup>H</sup></td><td>Maturity at puberty</td><td>—</td><td>J</td></tr><tr><td>&kappa;<sub>R</sub></td><td>Reproduction efficiency</td><td>0.95</td><td>—</td></tr></tbody></table>

The MATLAB package DEBtool<fnref n="3" /> implements the standard DEB model, extensions such as abj, and simulation code to run the model in different environmental conditions. The package also provides a full pipeline to estimate the parameters of DEB models.

## 2.2. The DEB parameter estimation problem

Formally, the DEB model is a dynamical system of the form (Kearney et al., 2021; Raue et al., 2009; Boussange et al., 2024):

$$ \frac{d}{dt}x(t) = g(x(t), u(t), \theta), \qquad (1) $$

$$ y(t) = h(x(t), \theta, \xi) + \epsilon(t), \qquad (2) $$

where $x(t) = (E, V, E_H, E_R)$ are the state variables, $g$ are the system equations (see Appendix A), $u(t)$ are the inputs to the system (environmental forcings such as temperature and food) and $\theta$ is the set of core parameters (see Table 1).

The state variables $x(t)$ of the DEB model are abstractions and cannot be measured in real life. Instead, observables $y(t)$ are variables that can be measured and as a result tampered by an observation error $\epsilon(t)$. Due to an extensive body of auxiliary theory, $y(t)$ can be almost any measurable quantity of an organism, from length and weight to reproduction rate and oxygen consumption (Kooijman, 2010), all of which are a function of the same parameters of Eqs. (1) and (2).

The relationship between observables $y(t)$ and state variables $x(t)$ is given by the mapping $h$, which also has auxiliary parameters $\xi$ (Lika

<footnote n="3">³ Available at https://github.com/add-my-pet/DEBtool_M.</footnote>
et al., 2011a). Some examples of parameters in ξ are the shape co-efficient, which relates observed length to structural length, and the chemical parameters that appear in the mass balance.

In the Add-my-Pet collection, DEB models are typically calibrated for a species; inter-individual variability is neglected and, when present, contributes to the observation error $e(t)$ (Marques et al., 2019; Oliveira et al., 2024). To obtain the parameters $\Theta^s = \{\theta^s, \xi^s\}$ of species $s$, the parameter estimation problem must be solved: from data $D^s$ on species $s$, find the parameters of the DEB model that lead to the best fit between the data and the DEB model.

The species data $D^s = \{D_1^s, ..., D_{|D|}^s\}$ is a collection of datasets $D_i^s$, where $|D^s|$ is the number of datasets. Each dataset $D_i = (s_{i1}, d_{i1}), ..., (s_{i|D_i|}, d_{i|D_i|})$ is an ordered set of pairs $(s_{ij}, d_{ij})$, where $s_{ij}$ is the independent variable and $d_{ij}$ is the dependent variable. Datasets can have any number of pairs and are categorized as zero-variate if they have a single data point, and uni-variate if they have more than one. Because the state variables are not observable, the dependent variable $d_{ij}$ is sampled from one of the components of $y(t)$. The independent variable $s_{ij}$ can take several forms:

* **Time** – The time since an initial state $x_i(0)$, specific to dataset $D_i^s$. An example is length as a function of time.

* **Observable** – In this case, the independent variable is itself an observable $y(t)$ and therefore a function of the state variables. An example is reproduction rate as a function of weight. The independent variable ‘weight’ is used alongside the forcings $u_i(t)$ and initial condition $x_i(0)$, to determine the state variables, which are then used to calculate the reproduction rate.

* **Event** – Here, the independent variable is categorical, which mathematically can be translated as a condition on the state variables. For example, age at birth is a zero-variate dataset with a single data point where the independent variable is the event ‘birth’, which occurs when the state variable $E_H$ (maturity) equals the parameter $E_H^b$ (maturity at birth).

The data collection $D^s$ can have datasets from various sources and experiments in the literature. This means that the initial condition $x_i(0)$ and inputs $u_i(t)$ of each dataset $D_i^s$ may not be equal for every dataset in the collection. In the case of the initial conditions, since the state variables are not observable, they must either be estimated or assumed. The inputs, however, are not always available because they are difficult to measure in the wild. Although temperature can be obtained from meteorological records, food intake is harder to obtain outside of carefully designed experimental setups. The most common case is for modelers to assume constant food availability and temperature. In particular, for the vast majority of zero-variate datasets, constant temperature and abundant food are assumed.

To assess the fit between the data and the DEB model, a collection of predictions $P$ is generated by simulating the state equations of the DEB model. Each individual prediction $P_i = (s_{i1}, p_{i1}), ..., (s_{i|D_i|}, p_{i|D_i|})$ has the same dimensions as the corresponding dataset $D_i^s$. Each datapoint prediction $p_{ij}$ is computed with a dataset-specific prediction function $q_i$ that takes as input the independent variables $s_{ij}$ and the model parameters $\Theta^s$. In the case where time is the independent variable, the prediction function $q_i$ is one of the components of the mapping $h$. If it is an observable or event, then the prediction function $q_i$ is a composition of components and inversions of components of the mapping $h$.

The set of core parameters $\theta^s$ that are estimated depends on the available data and usually does not include all 11 parameters in Table 1. The same occurs with the set of auxiliary parameters $\xi^s$, where shape coefficients are regularly estimated, but chemical parameters are almost always not. Parameters that are not estimated are set to a default value to avoid structural non-identifiability: when the available data is not sufficiently informative of a particular metabolic process, then parameters related to it cannot be uniquely determined, even in the absence of noise (Raue et al., 2009; Lika et al., 2011a). For example,

the species can have many weight curves, but without additional assumptions, parameters related to reproduction such as the reproduction efficiency $\kappa_R$ cannot be estimated.

We denote the set of parameters that are estimated as $\tilde{\Theta}^s = \{\tilde{\theta}^s, \tilde{\xi}^s\}$, where $\tilde{\theta}^s$ and $\tilde{\xi}^s$ are the set of estimated core and auxiliary parameters, respectively. The loss function is based on the mean squared error between data and predictions, but modified so it can handle the different scales and units in the data and to ensure multiplicative symmetry between data and predictions. It has the following formula:

$$ \underset{\tilde{\Theta}^s}{\text{minimize}} \quad \mathcal{L}(D^s, \mathcal{W}, \Theta^s) = \sum_{i=1}^{|D^s|} \sum_{j=1}^{|D_i^s|} \frac{w_{ij}}{|D_i^s|} \frac{(d_{ij} - p_{ij})^2}{d_i^2 + \bar{p}_i^2}, \quad (3) $$

$$ \text{with } \bar{d}_i = \frac{1}{|D_i^s|} \sum_{j=1}^{|D_i^s|} d_{ij}, \quad \bar{p}_i = \frac{1}{|D_i^s|} \sum_{j=1}^{|D_i^s|} p_{ij} \quad \text{and} \quad p_{ij} = q_i(s_{ij}, \Theta^s), $$

where $\mathcal{W} = \{W_1, ..., W_{|D|}\}$ is a set of weights for each dataset. The weight $w_{ij}$ quantifies the relative reliability of data point $j$ of dataset $D_i$ and is generally set to 1. Dividing it by the number of data points $|D_i|$ results in datasets of the same reliability contributing equally to the loss function regardless of how many data points they have (Marques et al., 2018, 2019).

Besides structural non-identifiability, a key issue in the parameter estimation of DEB models is practical non-identifiability. It occurs when the quality of the data is insufficient, leading to regions of the parameter space where the loss function changes negligibly, but where the behavior of the model can change significantly. The best way to solve this issue is to increase the quantity and/or quality of the estimation data. To address situations where this is not possible, the concept of pseudo-data was introduced (Lika et al., 2011a; Marques et al., 2019). Pseudo-data are a form of regularization and can be interpreted as priors on the values of some parameters or compound parameters (functions of parameters) (Oliveira et al., 2024). These are called pseudo-data because they take on the form of a zero-variate dataset in the loss function, although with a smaller weight (usually 0.1 for pseudo-data vs 1 for a dataset). The values for pseudo-data tend to be equal to the default values in Table 1, but they can also be taxa-specific or based on modeling decisions.

## 2.3. Estimation algorithms and initialization methods

The minimization problem in Eq. (3) is typically solved with optimization algorithms, such as the Nelder-Mead simplex algorithm (Nelder and Mead, 1965) or evolutionary algorithms (Robles et al., 2023; Filgueira et al., 2020). Approaches based on Bayesian inference have also been proposed (Johnson et al., 2013; Boersch-Supan et al., 2017; Boersch-Supan and Johnson, 2019; Palmer et al., 2024).

The Nelder-Mead simplex algorithm is the most widely used method. Nelder-Mead is a direct search optimization method that has been widely applied in multiple scientific areas (Vakil Baghmisheh et al., 2012; Jeslin Drusila Nesamalar et al., 2016). The fact that it does not need derivatives of the loss function makes it a popular choice to solve nonlinear problems, even though the algorithm does not have any convergence guarantees to the optimal value (McKinnon, 1998; Han and Neumann, 2006). In DEB calibration, Nelder-Mead is often restarted if it has failed to converge.

To initialize the estimation, modelers typically choose the estimated parameters of a similar species, where similar means that the species are taxonomically close and of comparable size. This procedure requires finding a similar species in the Add-my-Pet collection; however, a species of the same genus or family may not exist. Currently, there is no automated way to suggest a similar species.

Another option is to use the bijection method from (Lika et al., 2014), a 16-step algorithm with formulas, differential equations and iterative procedures. The bijection method estimates 8 of the core parameters and maturity at metamorphosis $E_H^j$ if it applies to the given
species. The input to the algorithm is a set of 9 easily obtainable zero-variate datasets: weight and age at several life stage events (birth, weaning or metamorphosis, and puberty), maximum weight, lifespan, dry/wet weight ratio, and reproduction rate. If some of these datasets are not available, then a subset of the parameters are not estimated and default values are used. This is the case for the parameters $k_j$, $\kappa_X$, $\kappa_P$ and $\kappa_R$, which are set to default values as the 9 easily obtainable zero-variate datasets do not contain information about these parameters.<fnref n="4" />

In addition to the requirement to minimize the deviation between data and model predictions, the set of parameters must also be feasible. Parameters have domain restrictions, and some parameter sets do not lead to an organism that can start feeding (birth), undergo puberty and reproduce (Marques et al., 2018). Meeting these parameter constraints, which are detailed in Section 4.1.2, is needed for two main reasons: to respect thermodynamic laws and for biological realism. However, in the case of biological realism, they are not sufficient; this should be assessed after calibration through estimation-in-context principles (Kooijman et al., 2024).

If a parameter set is close to infeasibility, then solving the model differential equations might become intractable or even cause the estimation to halt. The initial parameter set must therefore be feasible, otherwise the optimization procedure cannot begin. Any initialization method should then produce (i) feasible parameter sets, (ii) close to the parameters that minimize the loss function, so that the calibration is fast and avoids searching in areas of the parameter space that are close to infeasibility.

## 2.4. Machine learning calibration of DEB models

Let us consider a species s with an associated set of parameters $\Theta^s$. The objective of the calibration problem is to obtain a set of parameters that best describe the species, based on the available collection of data for the species $D^s$. When using a machine learning model, we attempt to approximate a mapping $\phi$ from data to parameters:

$$ \Theta^s = \phi(D^s) + \epsilon, \quad (4) $$

where $\epsilon$ is an irreducible error term that arises not only from problems with data gathering, such as lack of information on food and temperature, gaps in types of data, averaging data from multiple individuals and measurements errors, but also from the fact that the DEB model is a simplified representation of the metabolism of an organism (Marques et al., 2019; Lika et al., 2011a; Oliveira et al., 2024).

To formalize, let us consider a dataset for training a machine learning model with n species. The set of features $\mathcal{X} = \{D^1, \dots, D^s, \dots, D^n\}$ is the collection of data for each species. It is highly heterogeneous, as not all species will have the same types of data. This is due to both gaps in the available data and the different characteristics of each species. For example, species like fish that follow the abj DEB model can have additional datasets on the state at metamorphosis which do not exist for other species.

Besides missing data, each type of data also has different characteristics. Zero-variate data types (e.g. age at puberty) can be seen as tabular data, while uni-variate data types (e.g. weight as a function of time) are sequential data. Moreover, machine learning techniques open up the possibility of including categorical data, such as taxonomy or habitat characteristics, into the parameter estimation procedure.

The set of targets $\mathcal{Y} = \{\Theta^1, \dots, \Theta^s, \dots, \Theta^n\}$ contains the DEB model parameters of each species. Since parameters cannot be measured, these are obtained from solving the minimization problem in Eq. (3) for each species in the dataset. Due to the identifiability problems outlined previously, this set will often contain parameters that were not estimated and therefore have default values. Furthermore, in the

case of species that need DEB model extensions to properly model their life cycle, additional parameters are introduced, which means that the target parameter set is not the same for all species.

To build such a dataset, we can use the Add-my-Pet collection, an open-source database of DEB models. Entries to the collection are submitted independently and vetted by a board of curators (Marques et al., 2018; Kooijman et al., 2024), ensuring their quality and consistency. In addition to the estimation data and parameter estimates, each entry also includes modeling decisions and the code to replicate the estimation, which can be further used to evaluate the quality of parameter estimates produced by ML models. Moreover, each entry contains information on the species, such as habitat, climate, migration, and reproduction characteristics, which can be encoded and integrated into the model.

## 3. Data

The data for most species in the Add-my-Pet collection consists primarily of zero-variate datasets, such as maximum weight and reproduction rate. These datasets are easier to obtain from literature or from other open-source databases such as Avibase (Lepage, 2025) or Fishbase (Froese and Pauly, 2025, 2000). Uni-variate sequential data such as weight or length curves are typically measured within experimental trials, and the experimental data is often not shared or difficult to find (Tedersoo et al., 2021; Gomes et al., 2022; Michener, 2015).

For this reason, we train neural network models that take as input zero-variate data. We focus on the types of zero-variate data that the bijection method utilizes for two reasons. First, they are available for most species in the Add-my-Pet collection. Second, and perhaps most importantly, the mathematical foundation of the bijection method implies that these datasets can theoretically be used to estimate the core parameters of a DEB model.

Leveraging the flexibility of neural networks, we also include categorical data on taxonomy and ecological characteristics. These categorical variables will aid the models in learning the estimation-in-context principle: similar species should have similar parameter values.

### 3.1. Dataset construction

We parsed all species in the Add-my-Pet collection, which at the time of dataset construction (retrieved on 2025-04-10) contained DEB models for 5626 species. From these models, we retrieved only zero-variate data: age, time since birth,<fnref n="5" /> length, and weight at all maturity transitions (birth, metamorphosis, and puberty), as well as lifespan, maximum weight, and maximum length. Data on reproduction such as maximum reproduction rate and total reproduction were also collected. Further numerical columns include the dry/wet weight ratio and the typical body temperature.

In the case that more than one values exists for the same data type, we take the mean value. We also exclude data for males, which often vary from females in the values of parameters {$\hat{p}_{Am}$} and $E_H^p$. For data types that are dependent on body temperature, if multiple values exist, we choose the value that has a reported temperature closest to the typical body temperature of the organism. Additionally, every one of these datasets is normalized to the DEB reference temperature of $T_{ref} = 20^{\circ}$C with a temperature correction factor based on the Arrhenius temperature of the species (see Eq. (A.1)). This is consistent with the procedure employed by the bijection method (Lika et al., 2014; Kooijman, 2010).

In terms of categorical data, for each species we retrieved: the taxonomy, consisting of the taxonomic levels phylum, class, order,

<footnote n="4">⁴ See function get_pars in DEBtool for details.</footnote>

<footnote n="5">⁵ Age in DEB models starts at fertilization, however, many datasets are reported as time since birth (first-feeding) or time since hatching.</footnote>
<figure><img src="image_5.png" /><figcaption>Fig. 3. Percentage of species with a missing value for zero-variate datasets at several life stage events. Value for metamorphosis only include abj species.</figcaption></figure>

family, and genus; the DEB model type (std or abj); and the ecocodes on habitat, climate, food, and migration. Ecocodes are a classification of species traits of the Add-my-Pet collection, based on existing ecological classification systems. These ecocodes are extremely detailed, so to avoid an excessive number of categorical variables, we summarize them into broader categories that we believe significantly influence parameter values. These broader classifications are detailed in Appendix B.

## 3.2. Data imputation and cleaning

Fig. 3 shows the percentage of species with a missing value for zero-variate datasets. Although some datasets exist for almost all species, age at birth and puberty, and length or weight at puberty are missing for more than one third of the species in the collection. Values at metamorphosis are also missing for more than 90% of all abj species, which represent 45% of the species in the Add-my-Pet collection. To address this, we imputed some missing values using theoretical relationships from DEB theory. This theory-based imputation is performed independently for each record since it only uses the data of each species.

Under constant food, which is assumed for the collected zero-variate datasets, weight is proportional to cubed length. If a species has length and weight at the same maturity transition, we compute this proportionality constant and then use it to impute the missing weight values if the corresponding length value is available.

For bird species (class Aves), we noticed that most of them did not have values for weight at puberty. Following from a stylized fact that when birds reach puberty, their weight is close to the maximum, we impute this missing value with 95% of the ultimate weight of the species.

For species where the reproduction rate was missing, but the total reproduction, the lifespan and the age at puberty were not, we computed the reproduction rate as the total reproduction divided by an approximation of the reproductive period (lifespan minus age at puberty).

In the presence of gestation time, data types that report the time since birth data at a maturity transition were converted to age data types.

To ensure that all records had sensible data, we verified that after imputation both ages and weights at each maturity transition always increased from birth to metamorphosis, puberty, and adult. Species that did not meet this criteria were discarded from the dataset.

## 3.3. Final dataset characteristics

The goal of our machine learning models is to predict the core parameters of the DEB model of the species. Our models predict 7 of the 11 core parameters presented in Table 1: {$\dot{p}_{Am}$}, $\kappa$, $\dot{v}$, {$\dot{p}_{M}$}, $E_{H}^{b}$, $E_{H}^{p}$, and $k_{J}$. For abj species, the models also predict the maturity at metamorphosis $E_{H}^{l}$ and the acceleration factor $s_{M}$. The acceleration factor is not a parameter of the model, but instead a prediction of the DEB model. The prediction of the acceleration factor is necessary to ensure the feasibility of the parameter set, as will be detailed in Section 4.1.2.

We choose not to predict the parameters digestion efficiency $\kappa_{X}$, defecation efficiency $\kappa_{P}$, and reproduction efficiency $\kappa_{R}$ since the zero-variate data types collected are not informative enough to accurately estimate their values. Moreover, they are rarely estimated, which means their values are fixed and equal to the default values for the majority of species.

The specific cost for structure [$E_{G}$] is also not predicted by our models. This parameter is defined by the following formula:

$$[E_G] = \frac{\bar{\mu}_V d_V}{\kappa_G w_V}, \qquad (5)$$

where $\bar{\mu}_V$ is the chemical potential of structure V, $d_V$ the density of structure, $\kappa_G$ is the growth efficiency, $w_V$ is the molecular weight of structure. If we assume, similarly to the bijection method, that the densities of reserve $d_E$ and structure $d_V$ (in dry weight per volume) are equal and that the wet densities of reserve $d_{wE}$ and structure $d_{wV}$ (in wet weight per volume) are also equal and equal to 1 g cm$^{-3}$, we obtain that the dry-to-wet weight ratio $\delta_W$, that can be measured, is equal to $\delta_W = d_V/d_{wV} = d_E/d_{wE}$. The specific cost of structure is then predicted by the following equation:

$$[E_G] = \frac{\bar{\mu}_V \delta_W}{\kappa_G w_V d_{wV}}, \qquad (6)$$

with $d_{wV} = 1 g cm^{-3}$, and $\delta_W$ often chosen based on taxonomy. We observed that using this theoretical equation with default values for the three parameters ($\bar{\mu}_V = 550 000 J C^{-mol^{-1}}$, $\kappa_G = 0.8$, $w_V = 23.9 g C^{-mol^{-1}}$) and the dry-wet ratio $\delta_W$ for each species results in a $R^2 = 0.999$. This effectively renders its prediction by an ML model unnecessary and instead the theoretical equation can be applied.

We believe that this occurs due to the large weight in the pseudo-data of growth efficiency $\kappa_G$. This pseudo-data value is standard for all species, but its weight contribution in the loss function is frequently over half of the sum of weights. Minimizing the deviation from the default value of 0.8 is therefore a priority of the optimization procedure, which leads to it being mostly a function of the dry-wet ratio $\delta_W$.

When collecting the parameter values, we observed that most species were not in a minimum of the loss function. To ensure that the
<figure><img src="image_6.png" /><figcaption>Fig. 4. Structure of the neural network models.</figcaption></figure>

target parameter values align with the minimum of the loss function, we reestimated all species in the dataset. We applied repeated runs of the Nelder-Mead algorithm, each with 500 steps and alternating the direction of simplex size (either +0.02 or -0.02 of the current parameter values). We stopped the estimation when a run of Nelder-Mead converged, when the difference in loss function after two runs was less than $10^{-4}$, or when a maximum time limit of 15 h was reached. We refer to the parameters obtained after this procedure as target parameters.

For predictors, initially we considered all zero-variate datasets that the bijection method takes as input. However, this led to only 102 full records, even after data imputation, which we considered too few to train a machine learning model. To address this, we relaxed the requirement for species to have values for the data types age at puberty and wet weight at metamorphosis. This increased the size of the dataset to 2353 species, of which 780 (≈ 33%) have an abj DEB model.

# 4. Machine learning initialization methods

## 4.1. Neural network models

### 4.1.1. Model structure

The neural network model will map a set of features $\mathcal{X} \in \mathbb{R}^{n \times p}$ to a set of values $\mathcal{Y} \in \mathbb{R}^{n \times o}$, where $n = 2353$ is the number of species, $p = 34$ is the number of features, and $o = 10$ is the number of targets for each species.

The set of features includes 8 numerical variables: age at birth $a_b$, wet weight at birth $W_b$, wet weight at puberty $W_p$, lifespan $a_m$, ultimate wet weight $W_\infty$, maximum reproduction rate $R_\infty$, dry-wet weight ratio $\delta_W$, and typical body temperature $T_{typical}$.

The remaining 26 features are categorical and encoded as binary variables. These features are of three types: taxonomic information, habitat and species characteristics given by ecocodes and whether the typified DEB model of the species is abj or not. To avoid an excessive number of dummy variables, we restrict the taxonomy information of each species to the most frequent (> 5%) taxa of the class taxonomic level: Aves, Actinopterygii, Reptilia, Chondrichthyes, Amphibia, Mammalia, Bivalvia, Branchiopoda, Malacostraca. Any species that does not belong to one of these taxa is categorized as 'other'. These inputs are available for all species in the dataset, so the set of features is consistent and homogeneous.

The neural network does not directly output the core parameter values of a species. Instead, it predicts 10 outputs, which are functions of the core parameters, and then converts them to core parameter values. The model outputs are converted into 7 core parameter values or 8 if the species has an abj DEB model. A schematic of the structure of the neural network is presented in Fig. 4.

This encoding serves several purposes. First, it separates predictions for parameters that depend on the typified DEB model of the species. When a species does not have a parameter, then that value is substituted with an appropriate placeholder, and its weight in the loss function is set to zero so that it does not influence training. This addresses the issue described in 2.4 that the target parameters are not the same for all species.

Second, it allows us to predict only some of the core parameters, if there is not enough information to accurately estimate them all. In this case, we imputed the prediction for the fixed parameters and set their weight to zero in the loss function.

Third, and most importantly, it is used to impose constraints on model outputs such that the predicted parameter set respects most feasibility requirements, even when some parameters are fixed. Details on how these constraints are imposed are given in Section 4.1.2.

We represent all numerical inputs and parameter outputs in a transformed space. Specifically, we first apply a logarithmic transformation and then standardize each variable. The logarithmic transformation compacts the data and reduces the skew of the distributions, which facilitates training (Maier et al., 2023). Standardization then places all variables on a comparable scale.

Due to the large variety of species in the dataset, most input data types and parameter outputs span several orders of magnitude. If raw scales were used, the neural network training would prioritize reducing the deviation for species with large parameter values and disregard species with small parameter values. As a consequence, parameters with smaller values would be fitted less accurately, which would hinder the generalizability of the neural networks for these species. This effect is further exacerbated because the loss is averaged over all DEB parameters, so parameters with typically large values, such as {$\hat{p}_{Am}$} or {$\hat{p}_M$}, would receive a disproportionate weight during optimization.

For these reasons, training is carried out entirely in the log-standardized space. To obtain the core parameter values in original scale, we compute the output of the neural network and invert the transformations (see post-processing block of Fig. 4).

The neural network has several hidden layers with ReLU (rectified linear units) activation functions, which are initialized with values drawn from a Kaiming normal distribution (He et al., 2015). Each hidden layer is followed by a dropout layer. The number of hidden layers, the size of each hidden layer, and the dropout rate are tuned hyperparameters. The model is trained with mean squared error loss and optimized with the AdamW optimizer (Kingma and Ba, 2017; Loshchilov and Hutter, 2019). We tuned both the learning rate and the weight decay of the optimizer.
### 4.1.2. Implementing feasibility constraints

DEB parameter sets must respect the following constraints in order to be feasible:

1. All parameters must be positive

2. The allocation fraction to soma must be smaller than 1: $\kappa < 1$

3. All efficiencies must be smaller than 1: $\kappa_X, \kappa_P, \kappa_G, \kappa_R < 1$

4. Maturity thresholds must increase: $E_H^b < E_H^j < E_H^p$

5. The organism must reach birth

6. The organism must reach metamorphosis

7. The organism must reach puberty

Constraints 1–3 reflect thermodynamic principles whilst constraints 4–7 are needed to ensure biological realism. Constraint number 6 applies only to species that are modeled with an abj model.

Given that we apply a logarithm transformation to the outputs before training the model, the predictions of the model will also be in logarithmic scale. Since the logarithm is not defined for negative values, the model always predicts positive parameter values, and therefore constraint 1 is always respected.

To enforce constraint 2, we must restrict the prediction for $\kappa$ to be smaller than 1. However, the model instead predicts $1 - \kappa$ and enforces it to be smaller than 1. The distribution of $\kappa$ in our dataset is very right-skewed as most values are close to 1. This transformation, together with the logarithm transformation, renders it closer to a normal distribution, facilitating training.

The efficiency parameters are set to default values for the reasons detailed in 3.3, which means constraint 3 is always respected and no additional modifications are needed.

To ensure that maturity thresholds always increase, our models predict the maturity threshold $E_H^p$ and ratios of maturity thresholds:

$$s_H^{b/j} = \frac{E_H^b}{E_H^j}, \quad s_H^{j/p} = \frac{E_H^j}{E_H^p}, \quad s_H^{b/p} = \frac{E_H^b}{E_H^p}, \qquad (7)$$

where b refers to birth, j to metamorphosis, and p to puberty. The maturity ratios $s_H^{b/j}$ and $s_H^{j/p}$ are only predicted for abj species, while the maturity ratio $s_H^{b/p}$ is only predicted for non-abj species. By restricting these ratios to values between 0 and 1, constraint 4 is always respected. The estimates for the maturity thresholds $E_H^b$ and $E_H^j$ are then computed by inverting the maturity ratios in Eq. (7).

To check whether birth, and metamorphosis when relevant, can be reached, a differential equation must be solved. For this reason, we do not implement this constraint directly into the model structure. Nevertheless, given that for every example in the dataset birth and metamorphosis are always reached, we make the hypothesis that the model will learn to generate parameter sets that respect these constraints. We assess the validity of this hypothesis in Section 6.1.

Constraint 7, on the other hand, does not require solving a differential equation. For the simulation to reach puberty, the predicted parameter set must respect:

$$s_p = \frac{k_J E_H^p [\hat{p}_M]^2}{\kappa^2(1-\kappa)[\hat{p}_{Am}]^3} < s_M^3 J^3, \qquad (8)$$

where f is the scaled functional response, which measures the availability of food. It can take values between 0 and 1, where 1 means abundant food and is the standard food level for zero-variate data. To impose this constraint, the model predicts the ratio $s_p/s_M^3$ and enforces it to be smaller than 1: $s_p/s_M^3 < f^3 = 1^3$. This requires that the model also predicts $s_M = L_j/L_b$, which must be larger than 1, since $L_j > L_b$ always occurs after birth. To enforce this, the model instead predicts $1/s_M = L_b/L_j$ and we restrict its output to be smaller than 1.

For species that are not modeled with an abj model, the constraint is similar, except that $s_M = 1$ because metamorphosis does not occur. For these species, the prediction for $s_M$ is disregarded and substituted by the value 1.

The prediction of $s_p/s_M^3$ replaces the prediction of $\{\hat{p}_{Am}\}$, which is computed by inverting Eq. (8) and substituting the predictions of the model.

The constraints that the outputs $1 - \kappa$, $s_H^{b/j}$, $s_H^{j/p}$, $s_H^{b/p}$, $s_p/s_M^3$, and $1/s_M$ must be smaller than 1 are directly imposed in model scale, that is, after the logarithm transformation and standardization. Directly imposing the output constraints in model scale circumvents the need to invert the output scaling to ensure the constraints are respected. Due to the large ranges of each parameter, inverting the transformation can lead to excessively large gradients, which, in turn, can cause training instabilities.

Consider a model output $k$ that takes values between 0 and 1, and the corresponding prediction by the model $\hat{y}_k$. To enforce the constraint, the output before inverting the scaling $\hat{y}_k$ must respect:

$$0 < \hat{y}_k < 1 \Rightarrow -\infty < \log \hat{y}_k < 0 \Rightarrow \hat{y}_k = \frac{\log \hat{y}_k - \mu_{\log k}}{\sigma_{\log k}} < -\frac{\mu_{\log k}}{\sigma_{\log k}}, \quad (9)$$

where $\mu_{\log k}$ is the mean of $\log k$ in the training set, and $\sigma_{\log k}$ is the standard deviation. We impose the constraint in Eq. (9) both as a hard constraint and as a soft constraint, resulting in two models: a hard-constrained neural network, which we name HC-InitNet, and a soft-constrained neural network, which we refer to as SC-InitNet.

To implement it as a hard constraint, we apply a logsigmoid activation function to the output neuron and then add the constant shift $-\mu_{\log k}/\sigma_{\log k}$. The logsigmoid function maps $\mathbb{R}$ to $(-\infty, 0)$, which ensures that the standardized output satisfies (9).

To implement it as a soft constraint, we add additional terms to the loss function for each parameter that must be constrained to be smaller than 1. The loss $\mathcal{L}$ used to train this model is then given by:

$$\mathcal{L} = \mathcal{L}_{MSE} + \lambda_k \mathcal{L}_k + \lambda_{s_H^{b/j}} \mathcal{L}_{s_H^{b/j}} + \lambda_{s_H^{j/p}} \mathcal{L}_{s_H^{j/p}} + \lambda_{s_H^{b/p}} \mathcal{L}_{s_H^{b/p}} + \lambda_{s_p} \mathcal{L}_{s_p}, \quad (10)$$

where each $\lambda_k$ is the coefficient of the term in the loss function and is a tuned hyperparameter. Each term has the following form:

$$\mathcal{L}_k = \text{ReLU}\left(\frac{\hat{y}_k}{-\left(\frac{-\mu_{\log k}}{\sigma_{\log k}}\right)^2}\right), \qquad (11)$$

This term does not contribute to the loss function unless the output after unscaling is greater than 1, thus violating the constraint.

### 4.2. Mimicking the biologist with a nearest-neighbor model

The standard ad hoc procedure for initialization involves selecting parameters from a previously estimated species that is taxonomically and morphologically similar. We automate this process by developing 1-nearest-neighbor (1-NN) models that identify the most similar species based on a combined distance measure of taxonomy and body size. Nearest-neighbor models store the parameter values of the species in the training set. During prediction, the models select the parameters of the closest species, i.e., the species in the training set with the lowest distance (Fig. 5).

Since the taxonomic models select parameters from species already included in the Add-my-Pet collection – parameters which are feasible by construction – they always predict feasible parameter sets. However, if some core parameters are not estimated and are instead fixed at values provided by the modeler, the guarantee of feasibility no longer holds. This occurs because fixed parameters may differ from the nearest neighbor’s values, potentially leading to a predicted parameter set that does not respect all necessary constraints. To mitigate this risk, we introduce a penalty based on the logarithmic difference between the fixed parameter values and those of the candidate species. This penalty helps ensure that the predicted parameter set remains within feasible bounds, even after imputing fixed parameter values.

We refer to the standard taxonomic 1-nearest-neighbor model as Taxo1NN. We also develop a variant of the Taxo1NN model that adjusts the parameters of the closest species with body size scaling relationships, which we call SR-Taxo1NN. The scaling relationships emerge from DEB theory and describe how some parameter values and traits change (or do not change) when body size changes (see Chap.8 in Kooijman (2010) for more details).
<figure><img src="image_7.png" /><figcaption>Fig. 5. Structure of the taxonomic 1-nearest-neighbor models.</figcaption></figure>

### 4.2.1. Distance computation

We define the distance $d(s, r)$ between two species $s$ and $r$ as the sum of three terms, a taxonomic distance $d_T(s, r)$, an ultimate weight distance $d_{W_\infty}(s, r)$, and a fixed parameter distance $d_{\tilde{\theta}}$:

$$d(s, r) = d_T(s, r) + \lambda_W d_{W_\infty}(s, r) + \lambda_{\tilde{\theta}} d_{\tilde{\theta}}(s, r), \quad (12)$$

where $\lambda_W$ is a hyperparameter that controls the importance given to the ultimate weight and $\lambda_{\tilde{\theta}}$ is a fixed parameter penalty, which we set to a large value.

To compute the taxonomic distance, we first encode the taxonomy classifications into a vector $\mathcal{T}^s = (\mathcal{T}_1^s, \dots, \mathcal{T}_l^s, \dots, \mathcal{T}_{|\mathcal{T}|}^s)$, where each entry corresponds to the taxon of species $s$ for each taxonomic level $l$. We considered $|\mathcal{T}| = 5$ taxonomic levels: genus, family, order, class, and phylum.

Given this encoding, the taxonomic distance $d_T(s, r)$ is given by:

$$d_T(s, r) = \sum_{l=1}^{|\mathcal{T}|} 1 (\mathcal{T}_l^s \neq \mathcal{T}_l^r), \quad (13)$$

where $1(x)$ is the indicator function, equal to 1 if condition $x$ is true and 0 otherwise. This distance measures the number of different taxonomic levels between the two species.

The ultimate weight distance $d_{W_\infty}(s, r)$ is given by the difference of the logarithm of the ultimate weights:

$$d_{W_\infty}(s, r) = |\log W_\infty^s - \log W_\infty^r| = \left| \log \frac{W_\infty^s}{W_\infty^r} \right|, \quad (14)$$

where $W_\infty^s$ and $W_\infty^r$ are the ultimate weights of species $s$ and $r$, respectively.

Finally, when for some species a parameter is not estimated, we add a penalty based on the logarithm distance to the fixed value:

$$d_{\tilde{\theta}}(s, r) = \sum_{k=1}^{o} 1(k \notin \tilde{\theta}^s) \left| \log y_k^s - \log y_k^r \right| \quad (15)$$

where $y_k^s$ and $y_k^r$ are the values of parameter $k$ of species $s$ and $r$, respectively. The indicator function $1(k \notin \tilde{\theta}^s)$ takes the value 1 if the parameter is fixed and 0 otherwise. The predictions for fixed parameters are imputed with the target value and do not count toward the error of the model.

Similar to the neural network models, the taxonomic model also predicts $s_M$ and $s_p/s_M^3$, with {$\hat{p}_{Am}$} being obtained by inverting Eq. (8). This ensures that the predicted parameter set can reach puberty. However, if some parameters are not estimated, the predicted initial parameters is not found in the reference dataset, and therefore reaching birth and metamorphosis is no longer guaranteed. The large weight of the fixed parameter distance minimizes the chance that this occurs by prioritizing species whose parameters are close to the fixed parameters.

### 4.2.2. Adjusting predictions with scaling relationships

Consider a species $s$ and a species $\hat{r}$ that the SR-Taxo1NN model determined to be the closest to species $s$. The prediction for species $s$ will be the parameters of species $\hat{r}$, $\theta^r$. These species, however, will have different body sizes. DEB theory posits that some parameter values will change with the size of the organism simply as a result of physical principles (Lika et al., 2011a). These scaling relationships are not sufficient to predict the parameter values because selective pressures cause variations in this pattern. Nevertheless, they can be used to adjust the predictions and are typically used by modelers when choosing the initial parameter values. Moreover, adjusting parameters with scaling relationships does not violate any of the feasibility constraints.

The scaling of the parameter values is done with a zoom factor $z = L_m/L_m^{ref}$, where $L_m^{ref} = 1$ cm is the reference maximum structural length. Of the core parameters that the model predicts, the surface-specific maximum assimilation rate {$\hat{p}_{Am}$} scales with the zoom factor and the maturity levels scale with the cube of the zoom factor. Therefore, to adjust parameter values with the scaling relationships, we need to find the ratio between zoom factors of species $s$ and $\hat{r}$. This can be obtained from the ultimate weights using the proportionality between weight and structural length at constant food:

$$\frac{z^s}{z^{\hat{r}}} = \left( \frac{W_{\infty}^s}{W_{\infty}^{\hat{r}}} \right)^{\frac{1}{3}} \quad (16)$$

Then, the predictions of the taxonomic model for the parameters that vary with body size are given by:

$$\hat{y}_{E_H}^s = y_{E_H}^{\hat{r}} \left( \frac{z^s}{z^{\hat{r}}} \right)^3, \quad \hat{y}_{E_H}^s = y_{E_H}^{\hat{r}} \left( \frac{z^s}{z^{\hat{r}}} \right)^3, \quad \hat{y}_{E_H}^s = y_{E_H}^{\hat{r}} \left( \frac{z^s}{z^{\hat{r}}} \right)^3 \quad (17)$$

Because {$\hat{p}_{Am}$} is computed from $s_p/s_M^3$, the scaling relationship for this parameter does not need to be applied. When we scale the prediction for $E_H^p$, it indirectly applies the scaling relationship for {$\hat{p}_{Am}$}.

## 5. Experimental setup

All models were coded in Python 3.10.14. The neural network models, HC-InitNet and SC-InitNet, were implemented in PyTorch 2.5.1. The taxonomic 1-nearest-neighbor models, Taxo1NN and SR-Taxo1NN, were implemented in the framework of scikit-learn, version 1.5.1. To run the DEB model estimation, we used the DEBtool MATLAB package and MATLAB version R2024b. All experiments were performed on a computer with an Intel® Core™ i7-14700F 2.10 GHz CPU.

The data used to train the models was obtained from the Addmy-Pet collection. The estimation files for each species after reestimation and the final version of the data have been uploaded to Zenodo (Oliveira and Marques, 2025). The implementation of each
model as well as the code to reproduce the results of this paper are made available at the repository: https://github.com/diogo-f-oliveira/DEB_ML_Init. For instructions on how to use the proposed models in DEBtool, please refer to Section “Code availability and usage” at the end of the article.

## 5.1. Training and hyperparameter tuning

We set aside 15% of the dataset for testing the performance of all models. For the neural networks, we split the rest of the data to get a train-validation-test data split of 70%-15%-15%. For the Taxonomic 1-NN models, we instead perform 10-fold cross validation to select hyperparameters. We maintained the same proportion of abj and non-abj species in all splits, to ensure the same number of samples for parameters that are model-specific.

The hyperparameters of the neural networks that were tuned are: the number of hidden layers, the size of each hidden layer (equal for all), the dropout rate, the batch size, the learning rate, the strength of the weight decay, and the number of epochs. For the SC-InitNet model, we also tuned the loss function coefficients.

For the taxonomic nearest neighbor models, we only tune the coefficient of the ultimate weight distance $\lambda_W$. The fixed parameter penalty $\lambda_{\tilde{\theta}}$ is set to a large value of 100 to avoid that during hyperparameter tuning, prediction accuracy is prioritized over feasibility.

We used the implementation of HyperOpt (Bergstra et al., 2013) in Ray Tune (Liaw et al., 2018) as the tuning algorithm. Full details on the possible and final values for each hyperparameter are given in the Supplementary Material.

For both proposed approaches, the best configuration is chosen on the basis of the prediction error of the core parameter values instead of the prediction error of model outputs. The objective of our model is to predict core parameter values, and the transformations we apply are intended to improve their prediction. Therefore, model performance should be evaluated based on the prediction error calculated in the original scale of the parameters. Considering that parameters are measured on different scales, the evaluation in parameter scale must be performed using a relative error measure. We choose the logarithm of the accuracy ratio (Tofallis, 2015):

$$ \ln Q_k = \frac{1}{n_k} \sum_{s=1}^{n} \mathbf{1}(k \in \tilde{\Theta}^s) \left| \ln \left( \frac{y_k^s}{y_k^s} \right) \right|, \quad n_k = \sum_{s=1}^{n} \mathbf{1}(k \in \tilde{\Theta}^s), \quad (18) $$

where $y_k^s$ is the target value for a parameter k of species s, $\hat{y}_k^s$ is the corresponding model prediction, and $n_k$ is the number of species that have an estimate for k. The logarithm of the accuracy ratio is a symmetric measure of relative error, where both underestimation and overestimation are penalized equally. We choose the configuration of hyperparameters that leads to the smallest average log-accuracy ratio, averaged over all parameters.

## 5.2. Evaluation criteria

We evaluate the proposed approaches in three different criteria: feasibility, prediction performance, and initialization performance.

To assess feasibility, we run the DEBtool functions to verify the feasibility of the initial parameter set predicted by each model. These functions depend on the typified DEB model of the species. In addition to returning whether the parameter set is feasible, they also return the reason for infeasibility, that is, which of the 7 constraints listed in Section 4.1.2 is not respected. We imposed a time limit of 5 min to run the feasibility function. A longer execution time means that the differential equation numerical solver has likely halted, which suggests that the parameters are not feasible. We analyze the performance of each proposed method and compare them with the bijection method, the only automated initialization method currently available for researchers. These results are presented in Section 6.1.

Prediction performance measures how close the predicted initial parameters by each ML model are to the target parameter values. Since the models were trained to predict the target parameters based on the input data, the prediction accuracy criterion measures the extent to which the models performed this task successfully. This comparison is performed only on the test set to ensure that the models have not memorized the parameters during training and, in fact, learned to generalize.

For the reasons mentioned in the previous section, this comparison will be performed in terms of relative error. In addition to the logarithm of the accuracy ratio used during hyperparameter tuning, we also report the symmetric mean absolute percentage error (SMAPE) (Flores, 1986), which has the following formula:

$$ \text{SMAPE}_k = \frac{1}{n_k} \sum_{s=1}^{n} \mathbf{1}(k \in \tilde{\Theta}^s) \frac{|\hat{y}_k^s - y_k^s|}{|\hat{y}_k^s| + |y_k^s|} \times 100. \quad (19) $$

SMAPE has a minimum value of 0%, corresponding to a perfect prediction, and a maximum value of 100%, and, like the logarithm of the accuracy ratio, it is also a multiplicative symmetric measure of relative error. We report the values of both metrics for each ML model on the test set in Section 6.2.

To assess the performance of the proposed initialization methods, we run the DEB estimation procedure using the parameter predictions of each method. This is possible because the estimation files for each species in the Add-my-Pet collection are open source. For the estimation algorithm, we choose the standard Nelder-Mead algorithm, with the same settings as for the reestimation described in Section 3.3. We set a time limit of 1 h but allowed the current run of Nelder-Mead to finish if it started before the time limit ended. If the execution time exceeds 1.2 h, then we forcefully stop the estimation as it is a sign that the estimation is stuck in the differential equation solver.

We evaluated both the proposed ML initialization methods and the bijection method on the species of the test set. This is necessary to ensure a fair comparison between the ML models and the bijection method. If the initialization performance was assessed on the training set, the ML models would have an unfair advantage: the neural networks were trained to predict the target parameters and the taxonomic 1-nearest neighbor models have effectively memorized them by design.

For metrics, we consider the success rate of the estimation and enumerate the reasons why the estimation failed. We also compare both the initial and final calibration loss with the loss computed with the target parameters. Since each species can have a different number of datasets and a different weight for each dataset, we normalize the calibration loss with the sum of the weights in the estimation:

$$ \bar{\mathcal{L}}(D^s, \mathcal{W}^s, \Theta^s) = \mathcal{L}(D^s, \mathcal{W}^s, \Theta^s) / \sum_{i=1}^{|D^s|} \sum_{j=1}^{|D^s|} \frac{w_{ij}^s}{|D^s_i|} \quad (20) $$

The initialization performance results are presented and discussed in Section 6.3.

# 6. Results and discussion

## 6.1. Feasibility analysis

Table 2 shows the feasibility of the predicted initial parameter sets by the bijection method and each proposed ML method.

The bijection method produces feasible parameter sets for only 40.5% of species. It often fails completely (57.0%), throwing an error and not producing any initial parameter estimates. Even when it succeeds, it is not guaranteed that the parameter set is feasible.

We attribute this high failure rate to the rigidity of the method. Since the bijection method is a direct mapping between data and parameters, it does not consider that predictions may differ from the data. These differences, which can be caused by faulty measurements, individual variability, or assumptions that are not correct, lead to
<table><caption>**Table 2**
Feasibility results for each initialization method.</caption><thead><tr><th>Method</th><th>Feasibility outcome</th><th>Training set<br>(n = 2000)</th><th>Test set<br>(n = 353)</th><th>Success rate</th></tr></thead><tbody><tr><td rowspan="7">Bijection</td><td>Feasible</td><td>810</td><td>142</td><td rowspan="7">40.5%</td></tr><tr><td>Method failed</td><td>1142</td><td>200</td></tr><tr><td>Non-positive parameters</td><td>36</td><td>6</td></tr><tr><td>Non-increasing maturity levels</td><td>–</td><td>1</td></tr><tr><td>Unable to reach birth</td><td>6</td><td>3</td></tr><tr><td>Unable to reach puberty</td><td>3</td><td>1</td></tr><tr><td>Time limit reached</td><td>3</td><td>–</td></tr><tr><td rowspan="3">HC-InitNet</td><td>Feasible</td><td>1997</td><td>353</td><td rowspan="3">99.9%</td></tr><tr><td>Unable to reach birth</td><td>1</td><td>–</td></tr><tr><td>Unable to reach puberty</td><td>2</td><td>–</td></tr><tr><td rowspan="4">SC-InitNet</td><td>Feasible</td><td>1992</td><td>352</td><td rowspan="4">99.6%</td></tr><tr><td>Non-increasing maturity levels</td><td>–</td><td>1</td></tr><tr><td>Unable to reach birth</td><td>1</td><td>–</td></tr><tr><td>Unable to reach puberty</td><td>7</td><td>–</td></tr><tr><td>Taxo1NN</td><td>Feasible</td><td>2000</td><td>353</td><td>100.0%</td></tr><tr><td>SR-Taxo1NN</td><td>Feasible</td><td>2000</td><td>353</td><td>100.0%</td></tr></tbody></table>

**HC-InitNet** - Hard constrained neural network. **SC-InitNet** - Soft constrained neural network. **Taxo1NN** - Taxonomic 1-nearest-neighbors.

**SR-Taxo1NN** - Taxonomic 1-nearest-neighbors with scaling relationships.

deviations from the theoretical mapping and therefore errors in the estimates of each parameter. The errors then compound until the algorithm reaches a point that is incompatible with the next steps of the algorithm, resulting in execution errors. Even when the method does not fail, it can still produce nonsensical parameter estimates, in the form of negative or unreasonably large values, which are often the reason for the infeasibility.

The neural network models consistently produce feasible solutions, except for a very limited number of cases. Both neural network models failed to respect the birth constraint only one single time. As this constraint is not directly implemented, this result supports our hypothesis that, by learning from examples of feasible parameter sets, the model can learn to meet the constraint at birth and at metamorphosis.

The hard constrained neural network (HC-InitNet) did not meet the puberty constraint for two species, even though this requirement should be ensured by the network structure. The reason for this is that both species have a food level *f* (scaled functional response) below 1, while the constraint is implemented assuming the maximum value of *f* = 1. Since in these cases the food level is estimated together with the parameter values, the typical initial value for *f* would be 1. Although the network could be trained with a stricter constraint by assuming a lower value of *f*, this would come at the expense of prediction accuracy. Furthermore, if the user desires an initial *f* < 1, the puberty constraint can be met by adjusting the prediction for {\*p\*<sub>Am</sub>} by dividing by the cube of the desired scaled functional response.

As expected, the soft-constrained neural network (SC-InitNet) produces infeasible parameter sets more often than the HC-InitNet. Nevertheless, it produces infeasible solutions less than 1% of the time. Most infeasibility cases are due to the puberty constraint not being met. The SC-InitNet also fails to meet this constraint for the same two species with lower food level, but the remaining 5 species have a standard food level *f* = 1. All 7 species are birds, which in general have values of *s*<sub>p</sub> close to 1 in the Add-my-Pet collection, making it easier to break the constraint when it is not strictly enforced.

The taxonomic models always produce feasible parameter sets. Because these models are nearest neighbors where a single neighbor is used, the predictions for the training set are the parameters of the species in the training set. As such, these parameters are guaranteed to be feasible. In the test set, feasibility is no longer guaranteed for the reasons described in Section 4.2.1. Nevertheless, the fixed parameter distance term was successful at minimizing the chance that the predicted initial parameter set is infeasible.

## 6.2. Prediction performance of proposed machine learning models

Table 3 and Table 4 show the performance of the proposed ML approaches on the test set, measured with the logarithm of the accuracy ratio (ln Q) and the symmetric mean absolute percentage error (SMAPE), respectively. For completeness, we include the predictions of the bijection method for the 147 species where the method did not fail and where all predicted parameters were positive. Although the comparison with the bijection cannot be made strictly, the results are still indicative of the performance of the method.

Overall, the methods are able, on average, to predict parameter values within one order of magnitude of the target value. For all parameters except $E_H^f$, the best performing model has an average error within a factor of 2 of the target value. The low accuracy on this parameter is expected, particularly for the neural networks, since no data on metamorphosis is provided to the model. Nevertheless, they still achieve better performance than choosing the $E_H^f$ of the most similar species like the taxonomic 1-NN models do.

As evidenced by the large standard deviation of the errors, the error distribution is heavily skewed toward 0, with a tail of large errors. These large errors are the most problematic, since it increases both the calibration time and the likelihood that the search will stop at a local minimum or end in an execution error.

Both neural network approaches outperform the taxonomic 1-NN approaches as well as the bijection method, with the hard-constrained neural network performing best in 5 out of 8 parameters, according to the ln Q metric, and in 6 out of 8 parameters, according to the SMAPE metric. Fig. 6 plots the predictions of model HC-InitNet for all species, color coded by data split. The performance of the model is consistent across all magnitudes and is similar on all data splits, indicating that the model has not overfit.

For the parameters $\hat{v}_s$ [$\hat{p}_M$] and $\hat{k}_J$, a horizontal line can be observed around the default value for these parameters. This is a likely consequence of the pseudo-data values used during the calibration, which influenced the final value. Due to low identifiability, minimizing the deviation between the parameters and the pseudo-data does not come with the trade-off of worse fit between the DEB model and the data, resulting in the target parameters being significantly influenced by the pseudo-data.

The model did not consider pseudo-data during training; it instead learned the estimation-in-context principle by training with the calibrations of hundreds of species. The predictions of the neural network show that the principle is in disagreement with the pseudo-data values for these species, indicating that taxa-specific values may be a better approach in estimations with low identifiability.
<caption>**Table 3**</caption>

Logarithm of the accuracy ratio for each parameter (mean ± standard deviation). Lower values are better. The value in bold is the best result for the parameter. *Values for the bijection method only include only species where the method was successful and all parameters were positive.

<table><thead><tr><th>Parameter</th><th>Bijection*</th><th>HC-InitNet</th><th>SC-InitNet</th><th>Taxo1NN</th><th>SR-Taxo1NN</th></tr></thead><tbody><tr><td>{p&#770;<sub>Am</sub>}</td><td>1.270 ± 2.861</td><td>0.672 ± 0.697</td><td><strong>0.671 ± 0.675</strong></td><td>1.111 ± 1.126</td><td>1.080 ± 1.158</td></tr><tr><td>κ</td><td>0.156 ± 0.290</td><td><strong>0.097 ± 0.147</strong></td><td>0.101 ± 0.152</td><td>0.155 ± 0.209</td><td>0.158 ± 0.213</td></tr><tr><td>ψ̇</td><td>0.818 ± 1.507</td><td><strong>0.373 ± 0.431</strong></td><td>0.377 ± 0.439</td><td>0.489 ± 0.585</td><td>0.505 ± 0.598</td></tr><tr><td>[p&#770;<sub>M</sub>]</td><td>1.481 ± 4.295</td><td><strong>0.829 ± 0.839</strong></td><td>0.844 ± 0.827</td><td>1.205 ± 1.367</td><td>1.251 ± 1.448</td></tr><tr><td>E<sup>b</sup><sub>H</sub></td><td>1.142 ± 1.246</td><td>0.975 ± 0.918</td><td>0.949 ± 0.926</td><td><strong>0.866 ± 1.700</strong></td><td>1.667 ± 1.819</td></tr><tr><td>E<sup>h</sup><sub>H</sub></td><td>2.773 ± 2.517</td><td>2.055 ± 1.729</td><td><strong>2.021 ± 1.694</strong></td><td>2.837 ± 2.435</td><td>2.738 ± 2.371</td></tr><tr><td>E<sup>h</sup><sub>H</sub></td><td>0.911 ± 0.922</td><td><strong>0.87 ± 0.825</strong></td><td>0.893 ± 0.836</td><td>1.444 ± 1.549</td><td>1.387 ± 1.287</td></tr><tr><td>k<sub>j</sub></td><td>1.808 ± 3.659</td><td><strong>0.850 ± 1.395</strong></td><td>0.858 ± 1.342</td><td>1.121 ± 1.219</td><td>1.138 ± 1.211</td></tr></tbody></table>

**HC-InitNet** - Hard constrained neural network. **SC-InitNet** - Soft constrained neural network. **Taxo1NN** - Taxonomic 1-nearest-neighbors. **SR-Taxo1NN** - Taxonomic 1-nearest-neighbors with scaling relationships.

<caption>**Table 4**</caption>

Symmetric mean absolute percentage error (SMAPE) for each parameter (mean ± standard deviation). Lower values are better, with a maximum error of 100%. The value in bold is the best result for the parameter. Values for the bijection method only include only species where the method was successful and all parameters were positive.

<table><thead><tr><th>Parameter</th><th>Bijection*</th><th>HC-InitNet</th><th>SC-InitNet</th><th>Taxo1NN</th><th>SR-Taxo1NN</th></tr></thead><tbody><tr><td>{p&#770;<sub>Am</sub>}</td><td>35.60 ± 34.89</td><td><strong>28.70 ± 28.13</strong></td><td>28.81 ± 28.24</td><td>40.73 ± 39.92</td><td>39.05 ± 38.27</td></tr><tr><td>κ</td><td>7.39 ± 7.25</td><td><strong>4.78 ± 4.69</strong></td><td>4.99 ± 4.89</td><td>7.57 ± 7.42</td><td>7.72 ± 7.56</td></tr><tr><td>ψ̇</td><td>28.18 ± 27.62</td><td><strong>17.25 ± 16.91</strong></td><td>17.39 ± 17.05</td><td>21.55 ± 21.12</td><td>22.15 ± 21.71</td></tr><tr><td>[p&#770;<sub>M</sub>]</td><td><strong>32.84 ± 32.19</strong></td><td>33.74 ± 33.07</td><td>34.51 ± 33.82</td><td>40.55 ± 39.74</td><td>40.83 ± 40.01</td></tr><tr><td>E<sup>b</sup><sub>H</sub></td><td>40.37 ± 39.57</td><td>39.04 ± 38.26</td><td>38.03 ± 37.28</td><td><strong>23.84 ± 23.37</strong></td><td>49.77 ± 48.78</td></tr><tr><td>E<sup>h</sup><sub>H</sub></td><td>66.57 ± 65.24</td><td><strong>60.28 ± 59.08</strong></td><td>60.33 ± 59.12</td><td>68.08 ± 66.72</td><td>66.75 ± 65.42</td></tr><tr><td>E<sup>h</sup><sub>H</sub></td><td>36.14 ± 35.42</td><td><strong>36.12 ± 35.40</strong></td><td>36.61 ± 35.88</td><td>46.87 ± 45.94</td><td>47.97 ± 47.02</td></tr><tr><td>k<sub>j</sub></td><td>46.58 ± 45.65</td><td><strong>32.47 ± 31.83</strong></td><td>32.92 ± 32.27</td><td>39.87 ± 39.08</td><td>40.45 ± 39.65</td></tr></tbody></table>

**HC-InitNet** - Hard constrained neural network. **SC-InitNet** - Soft constrained neural network. **Taxo1NN** - Taxonomic 1-nearest-neighbors. **SR-Taxo1NN** - Taxonomic 1-nearest-neighbors with scaling relationships.

<figure><img src="image_8.png" /><figcaption>**Fig. 6.** Predictions for the hard constrained neural network (HC-InitNet) for all species in the dataset.</figcaption></figure>

The discrepancies in the model rankings are due to the SMAPE metric penalizing larger errors less than the ln Q metric. This can be clearly observed for the parameter $[ṗ_M]$, where the bijection method is the worst method according to ln Q, but the best method according to SMAPE. The bijection method can produce very accurate estimates for $[ṗ_M]$, since it can use data on the age at puberty $a_p$ whenever it is available, which greatly increases the identifiability of the parameter (Lika et al., 2014). This data type is not an input to the neural network models for data availability reasons; nevertheless, reasonable accuracy is still achieved without it. Conversely, the bijection method can sometimes produce extremely large estimates for $[ṗ_M]$, which considerably inflate the SMAPE. We attribute these cases to numerical instabilities that do not cause an execution error.

Despite these problems, the accuracy of the bijection method is on par with that of the taxonomic 1-NN models. Except for the parameter $E_H^b$, where the Taxo1NN model is overall the best, the performance of both taxonomic approaches is similar. The scaling relationships do not significantly improve the prediction of parameter values and, in the case of $E_H^b$, the predictions get considerably worse. This result hints that environmental pressures are likely to have a greater effect on parameter values than physical laws, as discussed in, e.g., Kooijman (2010), pp. 296–298.

Although patterns in parameter values have been found within and between taxa (Lika et al., 2011b; Marques et al., 2018), the higher accuracy of the neural network models shows that taxonomy may not be the main driver of parameter values, implying that species-specific adaptations are likely to have a substantial effect.
<table><caption>**Table 5**
Performance of initialization methods for species in the test set. Values in bold indicate the best performance for the metric.</caption><thead><tr><th>Metric</th><th>Bijection</th><th>HC-InitNet</th><th>SC-InitNet</th><th>Taxo1NN</th><th>SR-Taxo1NN</th></tr></thead><tbody><tr><td>Successful runs</td><td>112</td><td>307</td><td>316</td><td>272</td><td>283</td></tr><tr><td>Method failed</td><td>200</td><td>0</td><td>0</td><td>0</td><td>0</td></tr><tr><td>Infeasible solutions</td><td>11</td><td>0</td><td>1</td><td>0</td><td>0</td></tr><tr><td>Did not pass custom filters</td><td>1</td><td>7</td><td>4</td><td>8</td><td>10</td></tr><tr><td>Execution errors</td><td>10</td><td>20</td><td>18</td><td>49</td><td>41</td></tr><tr><td>Time limit exceeded</td><td>19</td><td>19</td><td>14</td><td>24</td><td>19</td></tr><tr><td>Did not reach dataset loss</td><td>46</td><td>103</td><td>111</td><td>98</td><td>103</td></tr><tr><td>Reached dataset loss</td><td>47</td><td>99</td><td>106</td><td>79</td><td>90</td></tr><tr><td>Improved on dataset loss</td><td>19</td><td>105</td><td>99</td><td>95</td><td>90</td></tr><tr><td colspan="6">Initial minus dataset normalized calibration loss</td></tr><tr><td>Mean</td><td>0.0708</td><td>0.0615</td><td><strong>0.0551</strong></td><td>0.1331</td><td>0.0891</td></tr><tr><td>Standard deviation</td><td>0.1674</td><td>0.0782</td><td><strong>0.0507</strong></td><td>0.4210</td><td>0.0923</td></tr><tr><td>Minimum</td><td><strong>0.0004</strong></td><td>0.0010</td><td>0.0018</td><td>0.0009</td><td>0.0012</td></tr><tr><td>Maximum</td><td>1.8367</td><td>0.9010</td><td><strong>0.4177</strong></td><td>7.4053</td><td>0.8745</td></tr><tr><td colspan="6">Final minus dataset normalized calibration loss</td></tr><tr><td>Mean</td><td>0.0028</td><td><strong>0.0008</strong></td><td>0.0009</td><td>0.0028</td><td>0.0031</td></tr><tr><td>Standard deviation</td><td>0.0230</td><td><strong>0.0083</strong></td><td>0.0087</td><td>0.0248</td><td>0.0236</td></tr><tr><td>Minimum</td><td>-0.0030</td><td>-0.0041</td><td>-0.0041</td><td><strong>-0.0059</strong></td><td>-0.0058</td></tr><tr><td>Maximum</td><td>0.2431</td><td><strong>0.1297</strong></td><td>0.1401</td><td>0.3897</td><td>0.3712</td></tr></tbody></table>

**HC-InitNet** - Hard constrained neural network. **SC-InitNet** - Soft constrained neural network. **Taxo1NN** - Taxonomic 1-nearest-neighbors.

**SR-Taxo1NN** - Taxonomic 1-nearest-neighbors with scaling relationships.

Another factor that may be influencing parameter values is the data available for estimation. The majority of species in Add-my-Pet are calibrated mostly with zero-variate data, which may bias parameter values. The accuracy of the bijection method and, in particular, of the neural network models supports this conclusion, since these methods predict parameters based solely on zero-variate datasets and not based on the full data collection available for parameter estimation.

As the scope and size of the Add-my-Pet collection continues to grow, the effect of scaling relationships, environmental pressures and calibration data can and should be studied in greater statistical detail. This will lead to improved calibration procedures and a greater ecological understanding through the analysis of DEB parameter patterns.

## 6.3. Initialization performance

Table 5 shows the results of running the estimation procedure for each species in the test set when starting from the initial set of parameters provided by each method. We include the numbers for when the models fail or produce infeasible parameters to show the complete overview of the performance of each method.

The estimation runs that were initialized with the proposed methods were successful for most species in the test set. For some runs, the initial parameter sets did not pass custom infeasibility filters sets defined by modelers. These custom filters reflect specific characteristics of the species and cannot be easily integrated into the initialization methods.

Execution errors were uncommon but occurred more often for the taxonomic models. These errors usually occur due to infeasible parameter sets that are not always caught by the code. The infeasibility of the parameters leads to an inability to solve the necessary differential equations, thereby causing an execution error. In addition, infeasible parameter sets can also halt the estimation when the differential equations cannot be solved. The row 'Time limit exceeded' shows the number of species where this occurred for each method.

Failing to complete the estimation is a cause of a suboptimal starting point. If the initial parameter set is closer to the true parameters of the species, then the likelihood of evaluating an infeasible solution decreases. This can be observed in the results of both neural network models — because their parameter accuracy is better, fewer errors are encountered.

In terms of the difference between the normalized loss computed with the initial parameters and the loss computed with the parameters in the dataset, the neural networks methods perform best. This is expected since they have the best parameter accuracy out of all

<figure><img src="image_9.png" /><figcaption>**Fig. 7.** Head-to-head comparison between initialization methods, considering only successful runs for the methods being compared. Annotations in each cell correspond to the Win-Draw-Loss record of the method in the row compared to the method in the column. Cells are colored according to the win rate of the method in the row, where draws count for half.</figcaption></figure>

methods tested. However, despite the hard constrained neural network having the best parameter accuracy, it has larger mean loss and higher dispersion than the soft constrained neural network.

When considering the full calibration, the neural networks are again the best performing method. The mean difference between the normalized loss computed with the calibrated parameters and the target parameters is 0.0008 for the HC-InitNet model and 0.0009 for SC-InitNet, which is roughly three times smaller than the mean values for the other methods.

Although the taxonomic nearest neighbor models achieved the best improvement on calibration loss, this performance is not consistent, as evidenced by the mean, standard deviation, and maximum value of the final calibration loss. In fact, in terms of calibration loss the initialization performance of the taxonomic 1-NN methods is similar to that of the bijection method.

**Fig. 7** shows the head-to-head comparison between each method, accounting only for species where both methods were successful. A
win means that the method in the row achieved a lower loss than the method in the column. If the difference between the losses is less than $10^{-4}$, which is the loss difference stopping criteria used throughout the paper, we declare a tie between the methods.

Many draws occurred, meaning that both methods obtained similar final loss. Neither model completely outperformed another, as all methods obtained wins and losses in the head-to-head comparison. In terms of win rate, the HC-InitNet outperforms all methods, achieving lower loss for more species than all other methods. The SC-InitNet was outperformed by the Taxo1NN model, despite having greater initial accuracy. The bijection method performed the worst, losing versus all other methods.

Although the neural network models were not capable of producing feasible initial parameter estimates for all evaluated species, they still achieved a higher success rate after considering the whole estimation process. On top of this, the neural networks achieve on average a better fit between data and DEB model than the taxonomic 1-NN models. The latter models replicate the typical initialization procedure of DEB practitioners. Our results show that when we consider the complete calibration, this procedure is less reliable and less accurate than employing our proposed neural networks.

Overall, the results demonstrate that the initial parameters significantly influence the estimation. A better starting point leads to an improved fit between the species data and the DEB model. All methods were able to improve on the loss obtained with the dataset parameters for many species, even though these were obtained after continuing the estimation from the parameters in the Add-my-Pet collection. Since all calibrations were performed with Nelder-Mead, a local algorithm, the results also hint that global optimization algorithms, such as particle swarm or genetic algorithms, may more consistently find the best fitting parameters. Given the complexity of the DEB calibration problem, we believe that the calibration procedure can be significantly improved through better algorithms and that this should be a future research direction for the field.

## 7. Conclusions and future work

Parameter estimation is a key step in the application of Dynamic Energy Budget theory, and accurate parameter estimates are crucial to guarantee useful model predictions. Although the process is largely automated with DEBtool, one part has not been suitably automated: where to start the calibration. Most estimations are initialized *ad hoc* from the parameters of a similar species. The only automated initialization is the bijection method of Lika et al. (2014), a 16-step map from data to parameters; however, our results show it often fails to produce a feasible initial parameter set, highlighting the need for new methods.

We addressed this issue with machine learning, developing two types of models: neural networks and nearest-neighbor models. The neural networks are inspired by the bijection method, predicting species parameters from a few commonly available datasets, while the nearest-neighbor models mimic selecting parameters of a similar species. Both models are built around the fundamental requirement of initial parameters: they must respect all parameter constraints that ensure the DEB model is biologically realistic.

We used more than 2000 DEB models from the Add-my-Pet collection to train and evaluate our models against the bijection method. All proposed models produce feasible parameter estimates for more than 99% of the species tested. The neural networks achieve the best prediction accuracy, with average errors within a factor of two of the target parameter values for all but one parameter. This translates into better initialization performance, yielding a lower final calibration loss and a higher success rate. The hard-constrained neural network, which strictly enforces parameter constraints, is the best method overall and outperforms all other methods in a head-to-head comparison.

Our results show that the starting point has a considerable impact on the estimation. Existing options, such as the bijection method or choosing parameters from a similar species, can cause early convergence to

<caption>**Table A.6**</caption>

Equations for state variables in the DEB model and powers associated to each process.

<table><thead><tr><th>State Variable</th><th>Dynamics Equation</th><th>Name</th></tr></thead><tbody><tr><td>E [J]</td><td>$\frac{dE}{dt} = \dot{p}_A - \dot{p}_C$</td><td>Reserve</td></tr><tr><td>V [cm<sup>3</sup>]</td><td>$\frac{dV}{dt} = \frac{\dot{p}_J}{[E_H]}$</td><td>Structure</td></tr><tr><td>E<sub>H</sub> [J]</td><td>$\frac{dE_H}{dt} = \begin{cases} \dot{p}_R & E_H < E_H^p \\ 0 & \text{otherwise} \end{cases}$</td><td>Maturity</td></tr><tr><td>E<sub>R</sub> [J]</td><td>$\frac{dE_R}{dt} = \begin{cases} 0 & E_H < E_H^p \\ \kappa_R \dot{p}_R & \text{otherwise} \end{cases}$</td><td>Reproduction Buffer</td></tr><tr><td>Powers [J d<sup>-1</sup>]</td><td>Equation</td><td>Name</td></tr><tr><td>$\dot{p}_A$</td><td>$\dot{p}_A = \kappa_x \dot{p}_X$</td><td>Assimilation</td></tr><tr><td>$\dot{p}_P$</td><td>$\dot{p}_P = \kappa_p \dot{p}_X$</td><td>Feces</td></tr><tr><td>$\dot{p}_C$</td><td>$\dot{p}_C = c_T E \frac{[E_H]^{v-1} + [\dot{p}_M]}{\kappa E^{v+1} + [E_H]}$</td><td>Mobilization</td></tr><tr><td>$\dot{p}_S$</td><td>$\dot{p}_S = c_T [\dot{p}_M] / V$</td><td>Somatic maintenance</td></tr><tr><td>$\dot{p}_G$</td><td>$\dot{p}_G = \kappa \dot{p}_C - \dot{p}_S$</td><td>Growth</td></tr><tr><td>$\dot{p}_J$</td><td>$\dot{p}_J = c_T k_J E_H$</td><td>Maturity maintenance</td></tr><tr><td>$\dot{p}_R$</td><td>$\dot{p}_R = (1 - \kappa) \dot{p}_C - \dot{p}_J$</td><td>Maturation/Reproduction</td></tr></tbody></table>

a local minimum or even a failure to converge. Our contribution is therefore twofold: we reduce trial-and-error during parameter estimation, thereby simplifying the calibration procedure, and improve the goodness-of-fit by providing a better initial point. We have made our code open source and collaborated with the DEBtool team to include our models as an initialization option, so that they can be easily used by the broader DEB community and facilitate ecological research.

The proposed models suggest initial values for species with an abj or s-model, which represent more than 90% of the Add-my-Pet collection. A logical next step is to extend the approach to other typified DEB models (e.g., h-models for insects). Model-specific feasibility requirements should be translated, whenever possible, to differentiable constraints, so that they can be integrated into machine learning pipelines. However, more species need to be calibrated so that enough data is available to train ML models.

The proposed models do not consider all available species data, which we see as an important direction for future work. Recurrent deep learning architectures that learn from sequence data could incorporate information from uni-variate data, whilst graph neural networks could better leverage the hierarchical taxonomic structure of the tree of life.

Another promising research line would be to use the DEB model prediction error to directly train a neural network. The network would predict the DEB parameters based on the available data, which would then be used to simulate the DEB model. The error between species data and DEB model predictions would then be backpropagated to update the weights of the network. Models with this architecture could potentially learn relationships between different data types, leading to more accurate parameter estimates, particularly for calibrations with limited data.

## CRediT authorship contribution statement

**Diogo F. Oliveira:** Writing – original draft, Visualization, Software, Methodology, Investigation, Formal analysis, Data curation, Conceptualization. **Gonçalo M. Marques:** Writing – review & editing, Validation, Supervision, Resources, Methodology, Funding acquisition, Data curation, Conceptualization. **Filipe M.P. Santos:** Writing – review & editing, Methodology. **Laure Pecquerie:** Writing – review & editing, Validation. **João M.C. Sousa:** Writing – review & editing, Supervision, Methodology, Funding acquisition, Conceptualization. **Tiago Dominatos:** Writing – review & editing, Supervision, Resources, Methodology, Funding acquisition, Conceptualization.
## Code availability and usage

The code to download and process the data, build and train the models, and reproduce the results of this paper is available on the GitHub repository at https://github.com/diogo-f-oliveira/DEB_ML_Init. Please refer to the README file in the repository for details on how it is structured.

The MATLAB port of the neural network model and its integration into DEBtool can be accessed in the directory DEBtool_M/lib/pet/DEBInitNet.

**Usage in DEBtool.** Make sure to download the latest available DEBtool version. To run the InitNet initialization, the code in the directory above must be added to the MATLAB path alongside the other DEBtool directories. Its usage is similar to other initialization methods in DEBtool. Simply choose option '3' as the pars_init_method: estim_options('pars_init_method', 3);

## Declaration of Generative AI and AI-assisted technologies in the writing process

During the preparation of this work the authors used ChatGPT in order to improve readability and language. After using this tool, the authors reviewed and edited the content as needed and take full responsibility for the content of the publication.

## 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.

## Acknowledgments

We would like to thank the anonymous reviewers for their constructive comments that contributed to improve this paper.

This work was supported by the Portuguese Foundation for Science and Technology (FCT), through IDMEC, under LAETA, project UIDB/50022/2020, by FCT/MCTES (PIDDAC) through projects UIDB/50009/2025, UIDP/50009/2025 and LA/P/0083/2020. This work was also supported by Fundo Europeu de Desenvolvimento Regional through project AVALON (COMPETE2030-FEDER-02288900). The work of Diogo F. Oliveira was supported by the PhD Studentship 2022.14030.BDANA from FCT. The work of Filipe M.P. Santos was supported by the PhD Studentship 2022.12077.BDANA from FCT.

## Appendix A. Standard DEB model equations

The dynamic equations for the standard DEB model are presented in Table A.6. The state equation for each state variable is a function of powers that describe how the energy intake of an organism is allocated to each metabolic process. The model takes as input a temperature *T*, with units K, and an ingestion power $\dot{p}_X$, with units J/d.

The variable $c_T$ is a temperature correction factor that describes how rates change with temperature *T* (in K):

$$c_T(T) = \exp\left(\frac{T_A}{T_{ref}} - \frac{T_A}{T}\right), \qquad (A.1)$$

where $T_A$ is the Arrhenius temperature (in K) and $T_{ref} = 20°C = 293.15$ K is a reference temperature. Formulations that model temperature tolerance ranges also exist, but are less common due to lack of species-specific data (Kooijman, 2010).

Food intake is based on food availability, measured by food density in the environment *X*. The scaled function response of the organism *f* follows a Michaelis–Menten or Holling type II function,

$$f = \frac{X}{X + K}, \qquad (A.2)$$

<table><caption>**Table B.7**
Simplified ecocode classifications used as categorical variables in the dataset.</caption><thead><tr><td>Ecocode</td><td>Categories</td></tr></thead><tbody><tr><td>Climate</td><td>Tropical, Dry, Temperate, Continental, Polar</td></tr><tr><td>Habitat</td><td>Terrestrial, Fresh Water, Inland Salt Water, Marine</td></tr><tr><td>Migration/Torpor</td><td>Torpor, No Torpor</td></tr><tr><td>Food</td><td>Carnivore, Herbivore, Planktivore, Omnivore, Other</td></tr></tbody></table>

where *K* is the half-saturation coefficient. The ingestion power can be written as:

$$\dot{p}_X = c_T \{\dot{p}_{Xm}\} f^{2/3}, \qquad (A.3)$$

where {$\dot{p}_{Xm}$} is the surface-specific maximum ingestion rate. The assimilation power $\dot{p}_A$ then becomes:

$$\dot{p}_A = c_T \{\dot{p}_{Am}\} f^{2/3}, \quad \{\dot{p}_{Am}\} = \kappa_X \{\dot{p}_{Xm}\} \qquad (A.4)$$

To model reproduction, the energy in $E_R$ must be emptied according to some reproduction rule that decides when this energy is used to generate offspring. The simplest approach is to assume a constant reproduction rate $\dot{R}$:

$$\dot{R} = \dot{p}_R / E_0, \qquad (A.5)$$

where $E_0$ is the initial reserve of an egg (Kooijman, 1986).

## Appendix B. Simplified ecological variables

The Add-my-Pet collection stores several ecological indicators of the species. In total, there are 8 categories: climate, ecozone, habitat, embryo environment, migration/torpor, food, gender, and reproduction. Full details are available at https://www.bio.vu.nl/thb/deb/deblab/add_myahl/AmPeco.html.

To avoid an excessive number of categories, we focus only on four categories that we believe have the most influence on parameter values. Furthermore, we reduce the number of options in each category to the most common options in our dataset. The possible values are reported in Table B.7.

## Appendix C. Supplementary data

Supplementary material related to this article can be found online at https://doi.org/10.1016/j.ecoinf.2026.103624.

## Data availability

We have added a section to the paper explaining how to access the data as well as the code. Data was shared through Zenodo and the code was shared through GitHub.

## References

Bergstra, J., Yamins, D., Cox, D., 2013. Making a Science of Model Search: Hyper-parameter Optimization in Hundreds of Dimensions for Vision Architectures. In: Proceedings of the 30th International Conference on Machine Learning. PMLR, pp. 115–123, (ISSN: 1938-7228).

Boersch-Supan, P.H., Johnson, L.R., 2019. Two case studies detailing Bayesian parameter inference for dynamic energy budget models. Ecosystem based management and the biosphere: a new phase in DEB research, J. Sea Res. Ecosystem based management and the biosphere: a new phase in DEB research, 143.57–69,

Boersch-Supan, P.H., Ryan, S.J., Johnson, L.R., 2017. deBInfer: Bayesian inference for dynamical models of biological systems in R. Methods Ecol. Evol. 8 (4), 511–518.

Boussange, V., Aceituno, P.V., Schäfer, F., Pellissier, L., 2024. Partitioning time series to improve process-based models with machine learning. bioRxiv:2022.07.25.501365 Section: New Results.

Desforges, J.-P.W., Sonne, C., Dietz, R., 2017. Using energy budgets to combine ecology and toxicology in a mammalian sentinel species. Sci. Rep. 7 (1), 46267, Publisher: Nature Publishing Group.
Edelstein-Keshet, L., 2005. Mathematical Models in Biology. In: Classics in Applied Mathematics, Society for Industrial and Applied Mathematics.

Elsawah, S., Filatova, T., Jakeman, A.J., Kettner, A.J., Zellner, M.L., Athanasiadis, I.N., Hamilton, S.H., Axtell, R.L., Brown, D.G., Gilligan, J.M., Janssen, M.A., Robinson, D.T., Rozenberg, J., Ullah, I.T., Lade, S.J., 2020. Eight grand challenges in socio-environmental systems modeling. Socio-Environmental Syst. Model. 2, 16226-16226.

Filgueira, R., Chica, M., Palacios, J.J., Strohmeier, T., Lavaud, R., Agüera, A., Damas, S., Strand, Ø., 2020. Embracing multimodal optimization to enhance Dynamic Energy Budget parameterization. Ecol. Model. 431, 109139.

Flores, B.E., 1986. A pragmatic view of accuracy measurement in forecasting. Omega 14 (2), 93-98.

Friedrich, T., Wagner, M., 2015. Seeding the initial population of multi-objective evolutionary algorithms: A computational study. Appl. Soft Comput. 33, 223-230.

Froese, R., Pauly, D., 2000. Fish Base 2000: concepts, design and data sources. ICLARM Contrib. (Philippines)

Froese, R., Pauly, D. (Eds.), 2025. FishBase. World Wide Web electronic publication, www.fishbase.org. (Accessed April 2025).

Gergs, A., Hager, J., Bruns, E., Preuss, T.G., 2021. Disentangling Mechanisms Behind Chronic Lethality through Toxicokinetic-Toxicodynamic Modeling. Environ. Toxicol. Chem. 40 (6), 1706-1712.

Giacoletti, A., Bosch-Belmar, M., Di Bona, G., Mangano, M.C., Stechele, B., Sarà, G., 2024. DEBEcoMod: A dynamic energy budget r tool to predict life-history traits of marine organisms across time and space. Ecol. Informatics 84, 102897.

Gomes, D.G.E., Pottier, P., Crystal-Ornelas, R., Hudgins, E.J., Foroughirad, V., Sánchez-Reyes, L.L., Turba, R., Martinez, P.A., Moreau, D., Bertram, M.G., Smout, C.A., Gaynor, K.M., 2022. Why don't we share data and code? Perceived barriers and benefits to public archiving practices. Proc. R. Soc. B: Biological Sci. 289 (1987), 20221113.

Grimm, V., Railsback, S.F., 2013. Individual-based Modeling and Ecology. Princeton University Press.

Gunawardena, J., 2014. Models in biology: 'accurate descriptions of our pathetic thinking'. BMC Biol. 12 (1), 29.

Haberle, I., Marn, N., Geček, S., Klanjšček, T., 2020. Dynamic energy budget of endemic and critically endangered bivalve *Pinna nobilis*: A mechanistic model for informed conservation. Ecol. Model. 434, 109207.

Han, L., Neumann, M., 2006. Effect of dimensionality on the Nelder-Mead simplex method. Optim. Methods Softw. 21 (1), 1-16.

He, K., Zhang, X., Ren, S., Sun, J., 2015. Delving Deep into Rectifiers: Surpassing Human-Level Performance on ImageNet Classification. In: 2015 IEEE International Conference on Computer Vision. ICCV, (ISSN: 2380-7504) pp. 1026-1034.

Jager, T., Goussen, B., Gergs, A., 2023. Using the standard DEB animal model for toxicokinetic-toxicodynamic analysis. Ecol. Model. 475, 110187.

Jager, T., Heugens, E.H.W., Kooijman, S.A.L.M., 2006. Making Sense of Ecotoxicological Test Results: Towards Application of Process-based Models. Ecotoxicology 15 (3), 305-314.

Jakeman, A.J., Elsawah, S., Wang, H.-H., Hamilton, S.H., Melsen, L., Grimm, V., 2024. Towards normalizing good practice across the whole modeling cycle: its instrumentation and future research topics. Socio-Environmental Syst. Model. 6, 18755-18755.

Jeslin Drusila Nesamalar, J., Venkatesh, P., Charles Raja, S., 2016. Managing multi-line power congestion by using Hybrid Nelder-Mead - Fuzzy Adaptive Particle Swarm Optimization (HNM-FAPSO). Appl. Soft Comput. 43, 222-234.

Johnson, L.R., Pecquerie, L., Nisbet, R.M., 2013. Bayesian inference for bioenergetic models. Ecology 94 (4), 882-894.

Kearney, M.R., Jusup, M., McGeoch, M.A., Kooijman, S.A.L.M., Chown, S.L., 2021. Where do functional traits come from? The role of theory and models. Funct. Ecol. 35 (7), 1385-1396.

Kim, Y., Son, J., Lee, Y.-S., Lee, M., Hong, J., Cho, K., 2019. Integration of an individual-oriented model into a system dynamics model: An application to a multi-species system. Environ. Model. Softw. 112.

Kingma, D.P., Ba, J., 2017. Adam: A Method for Stochastic Optimization. arXiv: 1412.6980 [cs].

Kooijman, S.A.L.M., 1986. What the hen can tell about her egg; egg development on the basis of budgets. J. Math. Biol. 23, 163-185.

Kooijman, B., 2010. Dynamic Energy Budget Theory for Metabolic Organisation, third ed. Cambridge University Press, Cambridge.

Kooijman, S.A.L.M., 2020. The standard dynamic energy budget model has no plausible alternatives. Ecol. Model. 428, 109106.

Kooijman, S.A.L.M., Kearney, M.R., Marn, N., Sousa, T., Domingos, T., Lavaud, R., Récapet, C., Klanjšček, T., Yeuw, T.T., Marques, G.M., Pecquerie, L., Lika, K., 2024. From formulae, via models to theories: Dynamic Energy Budget theory illustrates requirements. Ecol. Model. 497, 110869.

Lavaud, R., Filgueira, R., Augustine, S., 2021. The role of Dynamic Energy Budgets in conservation physiology. Conserv. Physiol. 9 (1), coab083.

Lepage, D., 2025. Avibase - The World Bird Database. https://avibase.bsc-eoc.org/ (Accessed April 2025),

Li, Q., Liu, S.-Y., Yang, X.-S., 2020. Influence of initialization on the performance of metaheuristic optimizers. Appl. Soft Comput. 91, 106193.

Liaw, R., Liang, E., Nishihara, R., Moritz, P., Gonzalez, J.E., Stoica, I., 2018. Tune: A Research Platform for Distributed Model Selection and Training. arXiv:1807.05118 [cs].

Lika, K., Augustine, S., Pecquerie, L., Kooijman, S.A.L.M., 2014. The bijection from data to parameter space with the standard DEB model quantifies the supply-demand spectrum. J. Theoret. Biol. 354, 35-47.

Lika, K., Kearney, M.R., Freitas, V., van der Veer, H.W., van der Meer, J., Wijsman, J.W.M., Pecquerie, L., Kooijman, S.A.L.M., 2011a. The "covariation method" for estimating the parameters of the standard Dynamic Energy Budget model I: Philosophy and approach. The AquaDEB project (phase II): what we've learned from applying the Dynamic Energy Budget theory on aquatic organisms, J. Sea Res. The AquaDEB project (phase II): what we've learned from applying the Dynamic Energy Budget theory on aquatic organisms, 66 (4).270-277,

Lika, K., Kearney, M.R., Kooijman, S.A.L.M., 2011b. The "covariation method" for estimating the parameters of the standard Dynamic Energy Budget model II: Properties and preliminary patterns. The AquaDEB project (phase II): what we've learned from applying the Dynamic Energy Budget theory on aquatic organisms, J. Sea Res. The AquaDEB project (phase II): what we've learned from applying the Dynamic Energy Budget theory on aquatic organisms, 66 (4).278-288,

Loshchilov, I., Hutter, F., 2019. Decoupled Weight Decay Regularization. arXiv:1711. 05101 [cs].

Maier, H.R., Galelli, S., Razavi, S., Castelletti, A., Rizzoli, A., Athanasiadis, I.N., Sánchez-Marrè, M., Acutis, M., Wu, W., Humphrey, G.B., 2023. Exploding the myths: An introduction to artificial neural networks for prediction and forecasting. Environ. Model. Softw. 167, 105776.

Marques, G.M., Augustine, S., Lika, K., Pecquerie, L., Domingos, T., Kooijman, S.A.L.M., 2018. The AmP project: Comparing species on the basis of dynamic energy budget parameters. PLoS Comput. Biol. 14 (5), e1006100.

Marques, G.M., Lika, K., Augustine, S., Pecquerie, L., Kooijman, S.A.L.M., 2019. Fitting multiple models to multiple data sets. Ecosystem based management and the biosphere: a new phase in DEB research. J. Sea Res. Ecosystem based management and the biosphere: a new phase in DEB research, 143.48-56,

Marques, G.M., Teixeira, C.M.G.L., Sousa, T., Morais, T.G., Teixeira, R.F.M., Domingos, T., 2020. Minimizing direct greenhouse gas emissions in livestock production: The need for a metabolic theory. Ecol. Model. 434, 109259.

Martins, I., Guerra, A., Vale, C.G., Xavier, C., Martins, I., Pinheiro, M., Neuparth, T., Xavier, J.R., Duarte, P., Santos, M.M., Colaço, A., 2024. Developing a dynamic energy budget model to project potential effects of deep-sea mining plumes on the Atlantic deep-sea mussel, *Bathymodiolus azoricus*. Ecol. Informatics 83, 102803.

McKinnon, K.I.M., 1998. Convergence of the Nelder-Mead Simplex Method to a Nonstationary Point. SIAM J. Optim. 9 (1), 148-158.

Michener, W.K., 2015. Ecological data sharing. Ecol. Informatics 29, 33-44.

Molnár, P.K., Derocher, A.E., Thiemann, G.W., Lewis, M.A., 2010. Predicting survival, reproduction and abundance of polar bears under climate change. Conservation planning within emerging global climate and economic realities, Biol. Cons. Conservation planning within emerging global climate and economic realities, 143 (7).1612-1622,

Nelder, J.A., Mead, R., 1965. A simplex method for function minimization. Comput. J. 7 (4), 308-313.

Oliveira, D.F., Marques, G.M., 2025. Dataset of "Reliable machine learning initialization methods for the calibration of Dynamic Energy Budget models" by Oliveira et al. (2025). Zenodo. https://doi.org/10.5281/zenodo.1752404.

Oliveira, D.F., Marques, G.M., Carolino, N., Pais, J., Sousa, J.M.C., Domingos, T., 2024. A multi-tier methodology for the estimation of individual-specific parameters of DEB models. Ecol. Model. 494, 110779.

Palmer, M., Moro-Martínez, I., Tomàs-Ferrer, J., Grau, A., López-Belluga, M.D., Herlin, M., Stavrakisidis-Zachou, O., Campos-Candela, A., 2024. Assessing between-individual variability in bioenergetics modelling: Opportunities, challenges, and potential applications. Ecol. Model. 498, 110848.

Pouvreau, S., Bourles, Y., Lefebvre, S., Gangnery, A., Alunno-Bruscia, M., 2006. Application of a dynamic energy budget model to the Pacific oyster, *Crassostrea gigas*, reared under various environmental conditions. Dynamic Energy Budgets in Bivalves, J. Sea Res. Dynamic Energy Budgets in Bivalves, 56 (2).156-167,

Raue, A., Kreutz, C., Maiwald, T., Bachmann, J., Schilling, M., Klingmüller, U., Timmer, J., 2009. Structural and practical identifiability analysis of partially observed dynamical models by exploiting the profile likelihood. Bioinformatics 25 (15), 1923-1929.

Robles, J.F., Chica, M., Filgueira, R., Agüera, A., Damas, S., 2023. MultiCalib4DEB: A toolbox exploiting multimodal optimisation in Dynamic Energy Budget parameters calibration. arXiv:2301.07548 [cs] version: 1.

Sousa, T., Domingos, T., Kooijman, S., 2008. From empirical patterns to theory: a formal metabolic theory of life. Phil. Trans. R. Soc. B 363 (1502), 2453-2464.

Sousa, T., Mota, R., Domingos, T., Kooijman, S.A.L.M., 2006. Thermodynamics of organisms in the context of dynamic energy budget theory. Phys. Rev. E, Stat. Nonlinear, Soft Matter Phys. 74 (5 Pt 1), 051901.
Tan, T.Y., Miraldo, M.C., Fontes, R.F.C., Vannucchi, F.S., 2022. Assessing bivalve growth using bio-energetic models. Ecol. Model. 473, 110069.

Tedersoo, L., Küngas, R., Oras, E., Köster, K., Eenmaa, H., Leijen, A., Pedaste, M., Raju, M., Astapova, A., Lukner, H., Kogermann, K., Sepp, T., 2021. Data sharing practices and data availability upon request differ across scientific disciplines. Sci. Data 8 (1), 192.

Tofallis, C., 2015. A better measure of relative prediction accuracy for model selection and model estimation. J. Oper. Res. Soc. 66 (8), 1352–1362.

Vakil Baghmisheh, M.T., Peimani, M., Sadeghi, M.H., Ettefagh, M.M., Tabrizi, A.F., 2012. A hybrid particle swarm–Nelder–Mead optimization method for crack detection in cantilever beams. Appl. Soft Comput. 12 (8), 2217–2226.