# On the use of a simple physical system analogy to study robustness features in animal sciences

Bastien Sadoul, Olivier Martin, Patrick Prunet, Nicolas N.C. Friggens

► To cite this version:

Bastien Sadoul, Olivier Martin, Patrick Prunet, Nicolas N.C. Friggens. On the use of a simple physical system analogy to study robustness features in animal sciences. PLoS ONE, 2015, 10 (8), pp.1-14. 10.1371/journal.pone.0137333. [https://doi.org/10.1371/journal.pone.0137333](https://doi.org/10.1371/journal.pone.0137333). [hal-01202872](https://hal.science/hal-01202872)

HAL Id: hal-01202872

https://hal.science/hal-01202872

Submitted on 21 Sep 2015

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>
<figure><img src="image_2.png" /></figure>

CrossMark

*click for updates*

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

OPEN ACCESS

**Citation:** Sadoul B, Martin O, Prunet P, Friggens NC (2015) On the Use of a Simple Physical System Analogy to Study Robustness Features in Animal Sciences. PLoS ONE 10(8): e0137333. doi:10.1371/journal.pone.0137333

**Editor:** Aditya Bhushan Pant, Indian Institute of Toxicology Reserach, INDIA

**Received:** April 2, 2015

**Accepted:** July 22, 2015

**Published:** August 31, 2015

**Copyright:** © 2015 Sadoul et al. This is an open access article distributed under the terms of the Creative Commons Attribution License, which permits unrestricted use, distribution, and reproduction in any medium, provided the original author and source are credited.

**Data Availability Statement:** The annotated code of this model was added to the supporting information (S1 Appendix) with an example of use. The code was written for the free software R used in the study.

**Funding:** This work was supported by funds from the Institut National de la Recherche Agronomique. The first author received a fellowship from the Région Bretagne with grant number 30000675. The funders had no role in study design, data collection and analysis, decision to publish, or preparation of the manuscript.

RESEARCH ARTICLE

# On the Use of a Simple Physical System Analogy to Study Robustness Features in Animal Sciences

Bastien Sadoul¹*, Olivier Martin²,³, Patrick Prunet¹, Nicolas C. Friggens²,³

1 INRA, LPGP Fish Physiology and Genomics UR1037, Rennes, France, 2 INRA, Modélisation Systémique Appliquée aux Ruminants UMR0791, Paris, France, 3 AgroParisTech, Modélisation Systémique Appliquée aux Ruminants UMR0791, Paris, France

* bastiensadoul@hotmail.fr

## Abstract

Environmental perturbations can affect the health, welfare, and fitness of animals. Being able to characterize and phenotype adaptive capacity is therefore of growing scientific concern in animal ecology and in animal production sciences. Terms borrowed from physics are commonly used to describe adaptive responses of animals facing an environmental perturbation, but no quantitative characterization of these responses has been made. Modeling the dynamic responses to an acute challenge was used in this study to facilitate the characterization of adaptive capacity and therefore robustness. A simple model based on a spring and damper was developed to simulate the dynamic responses of animals facing an acute challenge. The parameters characterizing the spring and the damper can be interpreted in terms of stiffness and resistance to the change of the system. The model was tested on physiological and behavioral responses of rainbow trout facing an acute confinement challenge. The model has proven to properly fit the different responses measured in this study and to quantitatively describe the different temporal patterns for each statistical individual in the study. It provides therefore a new way to explicitly describe, analyze and compare responses of individuals facing an acute perturbation. This study suggests that such physical models may be usefully applied to characterize robustness in many other biological systems.

## Introduction

The ability of an organism to respond to an environmental perturbation is of growing scientific concern in animal ecology and in animal production sciences [1-3]. Indeed, in the light of global climate change, the number and the suddenness of extreme meteorological events like drought, flood and storms is increasing. In parallel, despite controlled rearing conditions, environmental perturbations are frequent in animal production and can appear in an unpredictable manner. Depending on animal vulnerability, these environmental perturbations can affect the health, welfare, and fitness of the animal [4]. In this context, a significant effort is being devoted
**Competing Interests:** The authors have declared that no competing interests exist.

characterization of robust animals. However, the robustness of biological systems is known to be a complex feature difficult to access and characterize [5,6]. One major component of robustness is the capacity to cope with environmental perturbations, typically via elasticity in sub-elements of the system. Recent work suggests that quantifying robustness needs a multivariate approach and that key information about the ability of animals to face acute perturbations can be derived from the dynamic of responses to a perturbation [7,8]. The above studies have provided useful statistical tools for doing this and gaining some insight into the underlying biological processes. However, they do not per se attempt to describe basic features of the adaptive system. If response-recovery dynamics could be quantified in terms of generic components, this may facilitate a better characterization of adaptive capacity, and thereby robustness.

Physics and biology share many similarities, starting with adjectives and concepts [9]. Terms like plasticity, resistance or elasticity are, indeed, all borrowed from physics and allow the description of some robustness features of biological systems facing environmental perturbations. Thus, a robust biological system is, for some authors, expected to be resistant to acute environmental perturbations, and to show good recovery capacities [3,10,11]. Taking this one step further, it seems likely that the use of an interdisciplinary approach, combining physics and animal science, might help to characterize some generic features of the complex trait of robustness. Thus, we propose in this study that biological systems facing an acute perturbation can be described in terms of resistance and recovery abilities using the simple physics model of a spring and a damper. Similar models were previously developed for describing aspects of locomotion in human (see for review, [12] or for explaining complex hormonal oscillatory mechanisms [13] but never to characterize robustness features in the animal sciences. Therefore, this study aims to (1) investigate the suitability of this model to describe the responses of animals to perturbations, and (2) show a concrete example of how this model can be applied in practice in studies investigating the effect of a perturbation on animals. This second point required data on individuals facing an acute perturbation that were monitored before, during and after the perturbation. Data from a published study investigating the effect of a confinement perturbation on fish were used [8].

# Materials and Methods

## General approach

The model proposed in this study is mathematically formulated in terms of a generic model in which a spring and a damper are set in parallel, a model also called the Kelvin-Voigt model [14]. The deformation of the model is given by the x coordinates of the end-point (Fig 1). This model is composed of 2 main parameters describing the properties of the spring and damper (K and C). An environmental perturbation is represented as a force ($F_{pert}$) that pulls the spring and damper system. We assume that before the challenge, the system is in a non-challenging environment where the force on the model is zero. Because animals might perceive differently a challenge that is physically the same, we envisage the modulation of the force of the perturbation $F_{pert}$ by a coefficient of perception (perc; range 0 to 1).

## Model formalism

The system follows the differential Eq 1.1, determined using the second rule of Newton, where $a$ is the acceleration (i.e. $d^2x/dt$) and $v$ the velocity (i.e. $dx/dt$) of the end-point of coordinate $x$. During the pre-challenge period ($t < \tau_1$), the environment is considered as not affecting the system. The same assumption is made for the period after the challenge ($t \ge \tau_2$). The term $\text{percF}_{pert}$ is therefore considered as negligible during these periods, leading to the Eq 1.2. When the perturbation occurs ($\tau_1 \le t < \tau_2$), the system is affected by a non-negligible environmental
<figure><img src="image_4.png" /><figcaption>**Fig 1. The Kelvin-Voigt model in a non-perturbed environment (1.a) and during a perturbation (1.b).**
The model is composed of a spring and a damper in parallel and characterized respectively by the parameters K and C. During a perturbation a force of perturbation ($F_{pert}$) pulls on the system and the measure of interest, $x(t)$, is increased.</figcaption></figure>

doi:10.1371/journal.pone.0137333.g001

perturbation and Eq 1.1 can be rearranged to Eq 1.3.

$$ \text{perc}F_{\text{pert}} * a(t) = \text{perc}F_{\text{pert}} - Kx(t) - Cv(t) \quad (1.1) $$

$$ v(t) = -\frac{K}{C} x(t) \qquad t < \tau_1 \text{ or } t \ge \tau_2 \quad (1.2) $$

$$ a(t) = \frac{1}{\text{perc}F_{\text{pert}}} ( \text{perc}F_{\text{pert}} - Kx(t) - Cv(t) ) \quad \tau_1 \le t < \tau_2 \quad (1.3) $$

If the challenge was to persist, for $t \to +\infty$, because a and v would tend towards 0, we obtain from Eq 1.3 the asymptotic elastic solution $x_{inf} = \frac{\text{perc}F_{\text{pert}}}{k}$ (Fig 2). During the recovery ($t \ge \tau_2$), we define the value $T = C/K$ as the decay constant, characterizing the recovery capacity of the system (Fig 2) and $x_{max}$ the value of x at the end of the perturbation.

<figure><img src="image_5.png" /><figcaption>**Fig 2. An example of the outcomes of the model and its extrapolated features when continuously and constantly perturbed between $\tau_1$ and $\tau_2$.** We define the value $T = C/K$ as the decay constant, characterizing the recovery capacity of the system. The value $x_{inf}$ corresponds to the asymptotic value of x if the perturbation continues indefinitely, whereas $x_{max}$ is defined as the value of x at the end of a time-limited perturbation.</figcaption></figure>

doi:10.1371/journal.pone.0137333.g002
## Sensitivity analysis

In order to describe the flexibility of the model, an analysis of sensitivity was performed by varying one of the parameters K, C or $F_{pert}$, fixing arbitrarily the two others and running the model. The range of variation of the studied parameter and the values of the fixed parameters were chosen in order to create figures that clearly illustrate the impact of the variation in the studied parameter. We considered a situation where x is followed during 1000 iterations represented on a time scale of 100, i.e. with a time step of 0.1 time units, and where a constant perturbation $F_{pert}$ occurs between time 20 and 60.

First, values of K were chosen between 0.1 to 1, whilst C and $F_{pert}$ were set to 2 and 1 respectively. Secondly, values of C varied between 0.5 to 20, whilst K and $F_{pert}$ were set to 0.1 and 1 respectively. Finally, values of $F_{pert}$ varied between 1 and 10, whilst K and C were set to 0.1 and 2 respectively. For each analysis of sensitivity, the effect of the changing parameter on $x_{max}$, $x_{inf}$ and T was shown in a continuous graphic.

## Application of the model to a study on the effect of a perturbation in rainbow trout

Data of a study previously published [8] were used in this analysis. These are briefly summarized here, for further details please refer to Sadoul et al [8].

Eight aquaria of 1.7L containing each 16 juvenile rainbow trout of similar size were used. Half of the aquaria contained only fish from an isogenic line A22-B57 (called A), the other half contained fish from an isogenic line R23-B57. Cortisol release rate, oxygen consumption and group behavior (group dispersion and group activity) were followed at different time points distributed within 10 hours (Fig 3). A confinement challenge of 4 hours was performed on all aquaria by grouping all the fish of one aquarium in a net. The density was therefore suddenly increased up to 140 kg/m³ creating a stressful situation for the fish [15–17].

The experiment was carried out in two replicates one month apart with fish from the same fecundation. Within each replicate, both lines were used. Therefore, for each line a total of 8 aquaria was used.

## Data preparation

The data preparation steps were based on two hypotheses consisting of (1) considering that the unit of the measures do not have an impact on the parameters of the model (e.g. an increase of 1000 is not more important that an increase of 100 if the initial values are respectively 1000 and 100), and (2) the individuals in the aquarium were, at least at one time point, close to a non-disturbed state (generally just before the challenge or at the end of the recovery period).

Two data preparation steps were conducted to fit to these two hypotheses. In the first step, each measure of each individual was normalized using the mean of all aquarium of the corresponding measure, creating a fold change of the mean measure. The second step consisted in subtracting, for each individual and measure, the minimum of the normalized individual temporal pattern. All the temporal patterns were therefore shifted so that the minimum value was set to 0.

The 10 hours of the experiment were transformed onto a 0 to 100 scale, with the challenge occurring between 20 and 60. This arbitrary scale from 0 to 100 will help comparison with future analyses initially on different time scales.

## Model assumptions and model fitting

The perturbation force $F_{pert}$ was assumed to be positive, creating an extension of the spring. During the perturbation, $x(t)$ will therefore increase and the variations in the different
<figure><img src="image_6.png" /><figcaption>**Fig 3.** Design of the confinement experiment during time. Aquaria are followed during time for behavior, cortisol release rate and oxygen consumption. Each horizontal line represents a 10 minutes film record taken by camera placed above the aquaria. Vertical lines correspond to the water sampling points and the oxygen consumption measurements. At time 0, a 4-hour confinement challenge is performed by enclosing all the fish of one aquarium in a net, creating a density of 140kg/m3.</figcaption></figure>

doi:10.1371/journal.pone.0137333.g003

measures need to be positive. In this study we assume that all aquarium are subjected to the same perturbation (F<sub>pert</sub>), arbitrarily set to 0.1. For simplicity, the differences between aquaria in perception of the perturbation (perc) are assumed to be negligible here, and therefore, the perception is set to 1 for all individuals. In this study, fitting the model consists therefore in optimizing only the parameters K and C of the model. The system was described to be in its non-disturbed state when x = 0. The optimization procedure was performed using the *optim* function on differential Eqs 1.2 and 1.3 solved using the *ode* function in the free software R 2.14.0 (http://cran.r-project.org/). The code of the model and an example of use is given in S1 Appendix for a use with the software R.

The optimization procedure was first run on the mean temporal pattern of each measure, to assess the ability of the model to fit different experimental data.

The model was then optimized on each individual temporal pattern for each measure, providing a K<sub>i,m</sub>, and a C<sub>i,m</sub> for individual i and measure m.

## Model outputs analysis

Model fitting performance was calculated for each aquarium and each measure using the modeling efficiency statistic (MEF) as described in the review of Tedeschi [18] and calculated using the following formula, with Y<sub>i</sub> the observed value i, Ymod<sub>i</sub> the corresponding model value and ŶY the observed mean:

$$MEF = 1 - \frac{\sum_{i}(Y_{i} - Ymod_{i})^{2}}{\sum_{i}(Y_{i} - \bar{Y})^{2}}$$

A MEF value of one would indicate a perfect fit. On the contrary, a negative value of MEF suggests that the model values are worse than the observed mean.

For the purpose of the group level statistics, the few results where the model was not able to correctly fit the individual temporal pattern (when C or K is outwith ± 1.5 times the interquartile range) were withdrawn from the analysis. This corresponded to 1 aquarium (line A) for the cortisol release rate, 3 (2 line A, 1 line R) for the oxygen consumption, 0 for the group activity and 3 (1 line A, 2 line R) for the dispersion.

For statistical analyses, Ks and Cs were log-transformed to increase normality. For each measure, the effect of line and replicate on the corresponding K and C were assessed using a linear mixed model fitting each parameter with line, replicate and the interaction as fixed effects, and aquarium as a random effect. Non-significant fixed effects were withdrawn from the model and the model was re-run. The mean estimation of K and C for each line and the p-value for the line effect were extracted from the final model. Correlations between parameters were calculated using simple Pearson's correlation coefficients.
<figure><img src="image_7.png" /><figcaption>**Fig 4.** Sensitivity analysis of the model where the parameter K varies and a continuous perturbation occurs between time 20 and 60. For this analysis, perc was set to 1, Fpert to 1 and C to 2. K varied between 0.1 and 1. We define the value T = C/K as the decay constant, characterizing the recovery capacity of the system, $x_{max}$ the value of x at the end of the perturbation and $x_{inf}$ the value of x if the perturbation continues indefinitely.</figcaption></figure>

doi:10.1371/journal.pone.0137333.g004

All statistical analyses were run using the free software R 2.14.0 (http://cran.r-project.org/).

# Results

## Model flexibility

The sensitivity analysis showed that varying C and K in the model resulted in a strong diversity of dynamics (Figs 4 and 5). As expected, the parameter K has a strong effect on the amplitude of the response (Fig 4). The relationship between the maximum of the amplitude and K is approximately hyperbolic, with low values of K strongly increasing the amplitude of the response. Similarly, the time of recovery T is approximately hyperbolic to parameter K, with a strong decrease when K increases. Globally, parameter C has a strong impact on the shape of the response since it influences the speed of the deformation (Fig 5). This parameter has no impact on the asymptotic x ($x_{inf}$), however it impacts the $x_{max}$ value if the time of perturbation is not long enough to enable the system to reach $x_{inf}$. When the time of perturbation is low, the value $x_{max}$ and C are negatively correlated. The relationship between the time of recovery and C is proportional. The variable $F_{pert}$ has a positive linear impact on the amplitude of response $x_{max}$ and $x_{inf}$ (Fig 6). However, it has no impact on the recovery. The effect of K on the response

<figure><img src="image_8.png" /><figcaption>**Fig 5.** Sensitivity analysis of the model where the parameter C varies and a continuous perturbation occurs between time 20 and 60. For this analysis, perc was set to 1, K to 0.1 and Fpert to 1. C varies between 0.5 and 20. We define the value T = C/K as the decay constant, characterizing the recovery capacity of the system, $x_{max}$ the value of x at the end of the perturbation and $x_{inf}$ the value of x if the perturbation continues indefinitely.</figcaption></figure>

doi:10.1371/journal.pone.0137333.g005
<figure><img src="image_9.png" /><figcaption>**Fig 6.** Sensitivity analysis of the model where the parameter *Fpert* varies and a continuous perturbation occurs between time 20 and 60. For this analysis, *perc* was set to 1, *K* to 0.1 and *C* to 2. *F<sub>pert</sub>* varies between 1 and 10. We define the value *T* = *C*/*K* as the decay constant, characterizing the recovery capacity of the system, *x<sub>max</sub>* the value of *x* at the end of the perturbation and *x<sub>inf</sub>* the value of *x* if the perturbation continues.</figcaption></figure>

doi:10.1371/journal.pone.0137333.g006

is highly dependent on *C* (Fig 7). Globally, *K* has a stronger impact on the amplitude of the response when *C* is low. The parameter *K* defines therefore the asymptote value that a given perturbation would produce if it continued for long enough. We call *K* the “deformation potential” or the stiffness of the system. On the contrary, *C* has a strong effect on the shape of the response. The *C* effect is less impacted by the value of *K*, and can be defined as the “capacity to impede the deformation” or the resistance of the system to the deformation.

## Adjustment to mean responses in a real study

After data preparation (see dedicated section in *Materials and Methods*), the mean temporal pattern for each measure was calculated and presented, as points, in Fig 8. Globally, fish subjected to a confinement challenge display an increase of cortisol release rate, oxygen consumption and group behavior (activity or dispersion). The model was adjusted to the four mean responses to evaluate its ability to fit different physiological and behavioral responses to a

<figure><img src="image_10.png" /><figcaption>**Fig 7.** Simulated curves showing the interplay between *K* and *C* in a sensitivity analysis of the model where a continuous perturbation occurs between time 20 and 60. *K* and *C* take each three different (low, mid, high) values.</figcaption></figure>

doi:10.1371/journal.pone.0137333.g007
<figure><img src="image_11.png" /><figcaption>**Fig 8.** Fitting the model to cortisol release rate, oxygen consumption, group activity and group dispersion in rainbow trouts facing and recovering from a confinement challenge (grey rectangle).</figcaption></figure>

Data (observed values ± SE) are expressed as fold changes. The model fit is shown for each measure by the solid line. Stiffness (K) and resistance to change (C) parameters of the model were fitted on the mean of the cortisol release rate (K = 0.028, C = 0.384), group activity (K = 0.143, C = 1.137), oxygen consumption (K = 0.314, C = 3.052) and the group dispersion (K = 0.542, C = 9.809) measures.

doi:10.1371/journal.pone.0137333.g008

perturbation. The adjustments were, despite the amplitude differences, very satisfactory (con-
tinuous lines in Fig 8) with high model efficiency statistics (MEF = 0.98, 0.88, 0.90, 0.67 for the
cortisol release rate, the group activity, the oxygen consumption and the group dispersion,
respectively). The highest model values were obtained at the end of the confinement period.
The 4 measures show very divergent amplitudes of response to the challenge, translated by
contrasted K values. Cortisol release rate increased 4-fold relative to baseline values (x_max =
4.4). In contrast, group dispersion increased only 1.3 fold (x_max = 1.31). These differences in
amplitude result in a stiffness 11 times (K Disp / K Cort) stronger for group dispersion than for cor-
tisol (Table 1).

Adjustment to individual responses

The induced stress responses were adjusted using the model for each variable and each
individual aquarium, producing one K and C for each. The results were summarized in Table 1
presenting the average K and C for each line and each measure, and the significance of the line
effect on the parameters. The cortisol release data show a significant difference in K, with line
R being two times stiffer than line A. Similarly, line R had an almost significantly (p-value =
0.051) greater C value than line A. In consequence, the two lines show very similar recovery
rates (T) for the cortisol release rate. For the oxygen consumption, line R tends (p-value = 0.08)
to show a stronger K than line A. This is translated by the significantly greater xinf value for
line A. In contrast, line A showed a stronger stiffness in the group dispersion measure (p-
value<0.1) with a significantly lower xinf. For group activity, the two lines do not differ for the
<table><caption>Table 1. Model estimated (K and C) and derived (T and xinf) parameters, for two isogenic lines subjected to a confinement challenge.</caption><thead><tr><th rowspan="2"></th><th colspan="4">Cortisol release rate</th><th colspan="4">Oxygen consumption</th><th colspan="4">Group dispersion</th><th colspan="4">Group activity</th></tr><tr><th>K</th><th>C</th><th>T</th><th>x<sub>inf</sub></th><th>K</th><th>C</th><th>T</th><th>x<sub>inf</sub></th><th>K</th><th>C</th><th>T</th><th>x<sub>inf</sub></th><th>K</th><th>C</th><th>T</th><th>x<sub>inf</sub></th></tr></thead><tbody><tr><td>Line A</td><td>0.018</td><td>0.302</td><td>17.501</td><td>6.628</td><td>0.080</td><td>1.074</td><td>17.825</td><td>1.312</td><td>0.453</td><td>3.871</td><td>18.814</td><td>0.26</td><td>0.120</td><td>0.740</td><td>9.275</td><td>1.022</td></tr><tr><td>Line R</td><td>0.039</td><td>0.537</td><td>14.799</td><td>2.808</td><td>0.179</td><td>1.449</td><td>11.889</td><td>0.724</td><td>0.245</td><td>5.791</td><td>41.657</td><td>0.45</td><td>0.126</td><td>2.428</td><td>25.539</td><td>0.920</td></tr><tr><td>p-value</td><td>0.003</td><td>0.051</td><td>0.260</td><td>&lt;0.0001</td><td>0.089</td><td>0.526</td><td>0.270</td><td>0.016</td><td>0.058</td><td>0.609</td><td>0.174</td><td>0.028</td><td>0.8398</td><td>0.075</td><td>0.027</td><td>0.733</td></tr></tbody></table>

The significance of the line effect is given by the p-value extracted from a linear mixed model where K or C is explained by line and, if significant, replicate and interaction, including aquarium as a random effect.

doi:10.1371/journal.pone.0137333.t001

parameter K but a significant difference in the time of recovery (T) was found, due to the C parameter.

## Between parameter correlations

The correlations between the different K and C for each measure are shown in Table 2. Within measure K and C were positively correlated for all measures (except for group dispersion, p-value = 0.13). In addition, Ks of the two physiological parameters (water cortisol and oxygen consumption) were positively correlated. A similar trend was found between (p-value = 0.056) the Ks of the two behavioral parameters. Similar results were found for the Cs. The parameter C of the cortisol release rate tended to be negatively correlated with the behavioral Cs.

## Discussion

The purpose of this study was to see if a simple physics model of elasticity could be useful for describing an important diversity of temporal dynamics during and after an environmental perturbation of a biological system. Fitting the model to experimental data on multiple measures of stress in fish subjected to an acute perturbation demonstrated the capacity of the model to characterize the resistance and the ability to recover from a perturbation for each measure.

## The Model parameters

The environmental perturbation, by “pulling” on the system, creates a “distortion” and changes the values of the system measure. The temporal patterns of the four responses (cortisol release,

<table><caption>Table 2. Pearson's correlation coefficients between model estimated parameters.</caption><thead><tr><th rowspan="2"></th><th colspan="2">Cortisol release rate</th><th colspan="2">Group dispersion</th><th colspan="2">Oxygen consumption</th><th colspan="2">Group activity</th></tr><tr><th>K</th><th>C</th><th>K</th><th>C</th><th>K</th><th>C</th><th>K</th><th>C</th></tr></thead><tbody><tr><td rowspan="2">Cortisol release rate</td><td>K</td><td>1</td><td>0.94***</td><td>-0.17</td><td>-0.42</td><td>0.56*</td><td>0.36</td><td>-0.30</td><td>-0.31</td></tr><tr><td>C</td><td></td><td>1.00</td><td>-0.16</td><td>-0.42</td><td>0.51°</td><td>0.43°</td><td>-0.39</td><td>-0.47°</td></tr><tr><td rowspan="2">Group dispersion</td><td>K</td><td></td><td></td><td>0.39</td><td>0.43</td><td>0.07</td><td>0.49°</td><td>0.11</td><td></td></tr><tr><td>C</td><td></td><td></td><td>1.00</td><td>0.13</td><td>-0.22</td><td>0.33</td><td>0.52*</td><td></td></tr><tr><td rowspan="2">Oxygen consumption</td><td>K</td><td></td><td></td><td></td><td>1.00</td><td>0.61*</td><td>0.27</td><td>0.08</td><td></td></tr><tr><td>C</td><td></td><td></td><td></td><td></td><td>1.00</td><td>-0.18</td><td>-0.15</td><td></td></tr><tr><td rowspan="2">Group activity</td><td>K</td><td></td><td></td><td></td><td></td><td></td><td>1.00</td><td>0.54*</td><td></td></tr><tr><td>C</td><td></td><td></td><td></td><td></td><td></td><td></td><td>1.00</td><td></td></tr></tbody></table>

Significant correlations or tendencies (P<0.1) shown in bold (***: P<0.001, *: P<0.05 and °: P<0.1).

doi:10.1371/journal.pone.0137333.t002
O2 consumption, dispersion, and activity) were found to fit the dynamic equation of a spring and a damper. The parameter describing the stiffness of the system (K) provides an explicit measure of the deformation potential of the biological system when facing a given perturbation. The Ks can be compared between responses and individuals to describe which one can be the most modified by the perturbation. On the other hand, the parameter C of the equation, corresponding to the friction of the damper, can be considered as a measure of the capacity of the biological system to resist a change. With a similar K, a low C translates a rapid change of the system response. Rigidity, defined as a small response of the system to the perturbation, can be obtained using two different strategies when facing an acute perturbation, either by showing a strong K, reducing the potential of deformation of the system, or by showing a strong C, impeding the deformation. Therefore, during an acute perturbation, a strong K or a strong C lead to a good rigidity of the system (see Fig 7). These two strategies are probably the results of different underlying biological mechanisms.

The combination of the two parameters also enables the calculation of the recovery capacity of the system (T). This capacity is strong when C is low and K is high. Therefore, by considering the system using a physics model, important features of the system; resistance, rigidity and recovery capacity, can be explicitly described.

In animal science and ecology, resistance is considered to be an important ability of the animal. Resistance allows the animal to minimize changes in life functions that are important for productivity or reproduction when facing acute environmental perturbations [2,3]. Similarly, an animal taking short period of time to recover from a disturbance is of great interest since it can return quickly to its target levels of functioning [10,19,20]. However, the extent to which resistance is advantage or disadvantage for the animal depends on the nature of the perturbation. Given that, in our confinement experiment the perturbation used is not life-threatening, since it does not affect water quality or fish integrity, the individuals showing strong resistance and good recovery capacities are the individuals that avoid triggering unnecessary physiological and behavioral responses. In this experiment, we can therefore define individuals having strong resistance and good recovery capacities as robust animals [3].

## Physiological and behavioral dynamics analyzed using the model

All of the measures obtained from fish subjected to a confinement showed an increase during the perturbation compared to before. These variations were detectable and quantifiable using the present model. The model was found to adequately fit the data at the level of the individual, and therefore seems to be a relevant method for describing the responses of perturbed systems. In contrast to a previous analysis of these data [8], where the different periods of the measurement profiles (increase and decrease) were analyzed separately and sequentially, the use of the present model allowed us to analyze simultaneously all points of the temporal patterns and reduce therefore the effect of unexplained variance of isolated time points. Furthermore, the model provides comparable data in standardized common units facilitating between-challenge interpretation, between-response correlations or between-subject comparisons. Being able to describe temporal patterns in a way that allows inference across measures, or between individuals is often an important challenge in experimental studies [7,21–23].

The present model describes an increase and an early attainment of values close to the maximum response in group behavior during the challenge. Interestingly, although these dispersion and activity behaviors cannot be directly measured during the challenge–because it is a confinement challenge– the model provides a prediction of how they would evolve had the fish been able to express these behaviours. Even if in this experiment there was no way of testing
the validity of these predictions, it remains an interesting property of the model to be able to intuit the motivation to express a behaviour [24].

With respect to the physical measures, oxygen consumption and cortisol release, line R showed stronger resistance (stronger K and stronger or similar C) to the challenge compared to line A. The time to recover from the challenge was similar between the two lines for these physiological measures. These results confirm the findings in a previous analysis of the data [8]. However, the present analysis showed divergent group behavior responses between the two lines that were not detected in the previous analysis. Indeed, line A showed stronger stiffness to the challenge when considering the group dispersion behavior, and better recovery ability when considering the group activity behavior. Thus, these results indicate the existence of divergent coping strategies between the two lines, with one line being more sensitive in its physiological mechanisms and one line being more affected in its behavioral parameters.

## Model implications and potentials

The present model was constructed to be sufficiently generic to be able to describe very different responses. We have shown in this study its ability to fit with physiological and behavioral stress related measures. By using a common model across measures, we can extract from diverse temporal patterns with various scale units, comparable features of resistance and recovery capacities. Given this, it seems likely that such a model could help characterize perturbations of various types, e.g. social stress versus physical stress, and compare them on a same common scale. This model would, for example, have been of particular use for estimating the difference of perturbation intensity ($F_{pert}$ in our model) in welfare studies evaluating stress dynamics (sensitivity and recovery) depending on different intervention methods (see for review [25]).

To illustrate its genericity, the model was tested on data obtained from a different type of challenge, in a different species; a nutritional challenge in goats [26]. In this study, 16 lactating goats were followed every day for dry matter intake (DMI) and milk fat content (MFC). Goats initially fed a standard diet were suddenly subjected to a 100% straw diet during 2 days. The measures were recorded before, during and after the nutritional challenge. The present model was fitted to MFC and DMI measures as described in the Materials and Methods. Fig 9 shows that the model was able to correctly fit the responses to the challenge providing values of resistance and recovery for MFC. However, the fitting of the recovery period for DMI was not totally satisfactory since the temporal pattern of DMI seems to display a rebound phenomenon that the model cannot capture in its present form. This type of phenomenon has been observed in other measures [23,27] and could be captured by adapting the damped spring model to include a mass during the recovery phase. However, this and the possible elaborations of the model would require further work, to explore the conceptual value of such modifications.

In addition to using the model to characterize planned challenges, we can also envisage using the model for extracting and analyzing the perturbation dynamics due to unplanned micro-environmental variations in natural or rearing conditions. Indeed, new technologies based on video analysis or biosensors are rapidly developing and are increasingly being deployed to monitor domestic or wild animals. These can produce large amounts of dynamic data on behavioral and physiological stress mechanisms, which require new tools to analyze and interpret [8,28–30]. Recent work described a method to identify periods of deviance from the basic performance trajectory of cows monitored for their milk yield [31]. Combining this method and our model could potentially enable the identification and characterization of resistance and recovery capacities from dynamic data of monitored animals under natural conditions or in normal rearing systems. This would be an important step towards large scale
<figure><img src="image_12.png" /><figcaption>**Fig 9.** Fitting the model to Dry Matter Intake (DMI) and Milk Fat Content (MFC) in goats before, during and after a 2 days nutritional challenge (grey rectangle). Data (observed values ± SE) are expressed as fold changes. The inverse of DMI is illustrated since DMI is decreasing during the challenge. Stiffness (K) and resistance to change (C) parameters of the model were fitted on the mean of the DMI (K = 0.04, C = 0.06) and MFC (K = 0.72, C = 0.97) measures.</figcaption></figure>

doi:10.1371/journal.pone.0137333.g009

quantification of the adaptive capacity of animals under normal applied conditions. The basic model presented in this study could be extended to be able to accommodate other types of perturbation, e.g. chronic, combined or repeated perturbations. Indeed, by adding some specific features to the model, an updated model could potentially accommodate notions of stress accumulations, sequential stressors, simultaneous stressor, habituation or compensation.

Finally, the model was developed in order to characterize the effect of an acute stress during life history of an animal. However, we assume that this type of model can be used at smaller or broader dimension scales, from cell to species by looking at shorter to longer time periods. Thanks to the genericity of the model, it becomes for example conceivable to apply the model to the effect of climate change perturbations on species fitness traits. A sudden change in environment is expected to lead to strong alterations of fitness traits like physiological performances and reproductive outputs [32,33]. In consequences, populations affected by this perturbation avoid extinction through adaptive plasticity or evolutionary adaptation. These two mechanisms lead to a shift in phenotypes bringing the fitness traits back to the values prior to the perturbation [34]. The rate and the capacity of adaptation throughout generations can be assessed using the model described in this study. Indeed, by fitting the model to longitudinal data on fitness values of species subjected to an environmental change can help extracting values of the sensitivity and the adaptive capacities of the species, two indispensable features to estimate the vulnerability of a species to climate change [35].

## Conclusions

The spring-damper model appears to be a good, simple, solution for modeling responses to an acute perturbation. It provides a new way to explicitly describe, analyze and compare responses of individuals facing an acute perturbation. The model fitted satisfactorily to physiological and behavioral responses to a confinement challenge in fish. Given that it also fitted responses to a nutritional challenge in goats, it seems likely that it can adjust to many types of biological system responses. This suggests that such physical models may be usefully applied to characterize robustness in biological systems.
**Supporting Information**

S1 Appendix. R code for the spring and damper model.

(R)

**Author Contributions**

Conceived and designed the experiments: BS OM PP NCF. Performed the experiments: BS. Analyzed the data: BS NCF. Contributed reagents/materials/analysis tools: PP. Wrote the paper: BS OM NCF.

**References**

1. Korte SM, Koolhaas JM, Wingfield JC, McEwen BS. The Darwinian concept of stress: benefits of allostaticasis and costs of allostatic load and the trade-offs in health and disease. Neurosci Biobehav Rev. 2005; 29: 3–38. doi: 10.1016/j.neubiorev.2004.08.009 PMID: 15652252

2. Strandberg E, Klopčič M, Reents R, Philipsson J, Kuipers A. The role of environmental sensitivity and plasticity in breeding for robustness: lessons from evolutionary genetics. Breed Robustness Cattle M Klopčič R Reents J Philipsson Kuipers Ed EAAP Sci Ser. 2009; 17–33.

3. Wingfield JC. The comparative biology of environmental stress: behavioural endocrinology and variation in ability to cope with novel, changing environments. Anim Behav. 2013; 85: 1127–1133. doi: 10.1016/j.anbehav.2013.02.018

4. Korte SM, Olivier B, Koolhaas JM. A new animal welfare concept based on allostasis. Physiol Behav. 2007; 92: 422–428. doi: 10.1016/j.physbeh.2006.10.018 PMID: 17174361

5. Carlson JM, Doyle J. Complexity and robustness. Proc Natl Acad Sci U S A. 2002; 99 Suppl 1: 2538–2545. doi: 10.1073/pnas.012582499 PMID: 11875207

6. Kitano H. Biological robustness. Nat Rev Genet. 2004; 5: 826–837. doi: 10.1038/nrg1471 PMID: 15520792

7. Højsgaard S, Friggens NC. Quantifying degree of mastitis from common trends in a panel of indicators for mastitis in dairy cows. J Dairy Sci. 2010; 93: 582–592. doi: 10.3168/jds.2009-2445 PMID: 20105530

8. Sadoul B, Leguen I, Colson V, Friggens NC, Prunet P. A multivariate analysis using physiology and behavior to characterize robustness in two isogenic lines of rainbow trout exposed to a confinement stress. 2014;

9. Lesne A. Robustness: confronting lessons from physics and biology. Biol Rev. 2008; 83: 509–532. doi: 10.1111/j.1469-185X.2008.00052.x PMID: 18823391

10. De Goede DM, Gremmen B, Blom-Zandstra M. Robust agriculture: Balancing between vulnerability and stability. NJAS—Wagening J Life Sci. 2013; 64–65: 1–7. doi: 10.1016/j.njas.2012.03.001

11. Sauvant D, Martin O. Robustesse, rusticité, flexibilité, plasticité... les nouveaux critères de qualité des animaux et des systèmes d'élevage: définition systémique et biologique des différents concepts. Prod Anim. 2010; 23: 5–9.

12. Nikooyan AA, Zadpoor AA. Mass-spring-damper modelling of the human body to study running and hopping—an overview. Proc Inst Mech Eng [H]. 2011; 225: 1121–1135.

13. Zeeman ML, Weckesser W, Gokhman D. Resonance in the menstrual cycle: a new model of the LH surge. Reprod Biomed Online. 2003; 7: 295–300. PMID: 14653886

14. Marques SPC, Creus GJ. Computational Viscoelasticity. Springer Science & Business Media; 2012.

15. Auperin B, Baroiller JF, Ricordel MJ, Fostier A, Prunet P. Effect of Confinement Stress on Circulating Levels of Growth Hormone and Two Prolactins in Freshwater-Adapted Tilapia (Oreochromis niloticus). Gen Comp Endocrinol. 1997; 108: 35–44. doi: 10.1006/gcen.1997.6938 PMID: 9378272

16. Pottinger TG, Carrick TR. Modification of the Plasma Cortisol Response to Stress in Rainbow Trout by Selective Breeding. Gen Comp Endocrinol. 1999; 116: 122–132. doi: 10.1006/gcen.1999.7355 PMID: 10525368

17. Quillet E, Krieg F, Dechamp N, Hervet C, Bérard A, Le Roy P, et al. Quantitative trait loci for magnitude of the plasma cortisol response to confinement in rainbow trout. Anim Genet. 2014; 45: 223–234. doi: 10.1111/age.12126 PMID: 24444135

18. Tedeschi LO. Assessment of the adequacy of mathematical models. Agric Syst. 2006; 89: 225–247. doi: 10.1016/j.agsy.2005.11.004

19. Knap PW. Breeding robust pigs. Aust J Exp Agric. 2005; 45: 763–773.
20. Napel J ten, Bianchi F, Bestman M. Utilising intrinsic robustness in agricultural production systems: Inventions for a sustainable development of agriculture. TransForum Agro & Groen; 2006. Available: http://core.kmi.open.ac.uk/display/10929438

21. Fertuck KC, Eckel JE, Gennings C, Zacharewski TR. Identification of temporal patterns of gene expression in the uteri of immature, ovariectomized mice following exposure to ethynylestradiol. Physiol Genomics. 2003; 15: 127-141. doi: 10.1152/physiolgenomics.00058.2003 PMID: 12915738

22. Schuring E, Quadt F, Kovacs EMR, Meullenet J-F, Wiseman S, Mela DJ. A quantitative method for estimating and comparing the duration of human satiety responses: Statistical modeling and application to liquid meal replacers. Appetite. 2012; 59: 601-609. doi: 10.1016/j.appet.2012.07.003 PMID: 22796948

23. Desnoyers M, Giger-Reverdin S, Duvaux-Ponter C, Sauvant D. Modeling of off-feed periods caused by subacute acidosis in intensive lactating ruminants: Application to goats. J Dairy Sci. 2009; 92: 3894-3906. doi: 10.3168/jds.2008-1784 PMID: 19620672

24. Marshall JA, Favreau-Peigne A, Fromhage L, Mcnamara JM, Meah LF, Houston Al. Cross inhibition improves activity selection when switching incurs time costs. Curr Zool. 2015; 61: 242-250.

25. Mellor DJ, Cook CJ, Stafford KJ. Quantifying some responses to pain as a stressor. Biol Anim Stress Basic Princ Implic Welf. 2000; 171-198.

26. Friggens NC, Duvaux-Ponter C, Tessier J-A, Schmidely P. Responses to a nutritional challenge in early and late lactation. 2012 Joint Annual Meeting. Phoenix, United States; 2012. p. 671. Available: https://hal.archives-ouvertes.fr/hal-01019075

27. Walker JJ, Spiga F, Waite E, Zhao Z, Kershaw Y, Terry JR, et al. The Origin of Glucocorticoid Hormone Oscillations. PLoS Biol. 2012; 10: e1001341. doi: 10.1371/journal.pbio.1001341 PMID: 22679394

28. Ellis T, James JD, Stewart C, Scott AP. A non-invasive stress assay based upon measurement of free cortisol released into the water by rainbow trout. J Fish Biol. 2004; 65: 1233-1252. doi: 10.1111/j.0022-1112.2004.00499.x

29. Wilson GS, Gifford R. Biosensors for real-time in vivo measurements. Biosens Bioelectron. 2005; 20: 2388-2403. doi: 10.1016/j.bios.2004.12.003 PMID: 15854814

30. Takase M, Murata M, Hibi K, Huifeng R, Endo H. Development of mediator-type biosensor to wirelessly monitor whole cholesterol concentration in fish. Fish Physiol Biochem. 2014; 40: 385-394. doi: 10.1007/s10695-013-9851-1 PMID: 24037271

31. Codrea MC, Højsgaard S, Friggens NC. Differential smoothing of time-series measurements to identify disturbances in performance and quantify animal response characteristics: An example using milk yield profiles in dairy cows. J Anim Sci. 2011; 89: 3089-3098. doi: 10.2527/jas.2010-3753 PMID: 21478456

32. Pankhurst NW, Munday PL. Effects of climate change on fish reproduction and early life history stages. Mar Freshw Res. 2011; 62: 1015-1026.

33. Visser ME, Noordwijk AJ van, Tinbergen JM, Lessells CM. Warmer springs lead to mistimed reproduction in great tits (Parus major). Proc R Soc B Biol Sci. 1998; 265: 1867-1870. doi: 10.1098/rspb.1998.0514

34. Visser ME. Keeping up with a warming world: assessing the rate of adaptation to climate change. Proc R Soc B Biol Sci. 2008; 275: 649-659. doi: 10.1098/rspb.2007.0997

35. Williams SE, Shoo LP, Isaac JL, Hoffmann AA, Langham G. Towards an Integrated Framework for Assessing the Vulnerability of Species to Climate Change. PLoS Biol. 2008; 6: e325. doi: 10.1371/journal.pbio.0060325