A Theoretical Comparison of Alternative Male Mating Strategies in Cephalopods and Fishes

We used computer simulations of growth, mating and death of cephalopods and fishes to explore the effect of different life-history strategies on the relative prevalence of alternative male mating strategies. Specifically, we investigated the consequences of single or multiple matings per lifetime, mating strategy switching, cannibalism, resource stochasticity, and altruism towards relatives. We found that a combination of single (semelparous) matings, cannibalism and an absence of mating strategy changes in one lifetime led to a more strictly partitioned parameter space, with a reduced region where the two mating strategies co-exist in similar numbers. Explicitly including Hamilton’s rule in simulations of the social system of a Cichlid led to an increase of dominant males, at the expense of both sneakers and dwarf males (“super-sneakers”). Our predictions provide general bounds on the viable ratios of alternative male mating strategies with different life-histories, and under possibly rapidly changing ecological situations.

Authors and Affiliations

1. Department of Mathematics, University of Pittsburgh, Pittsburgh, PA, USA

   Joseph A. Landsittel & G. Bard Ermentrout

2. Neurolinx Research Institute, La Jolla, CA, USA

   Klaus M. Stiefel

3. Silliman University Angelo King Center for Research & Environmental Management, Hibbard Avenue, 6200, Dumaguete City, Philippines

   Klaus M. Stiefel

Corresponding author

Correspondence to Klaus M. Stiefel.

Publisher's Note

Springer Nature remains neutral with regard to jurisdictional claims in published maps and institutional affiliations.

Appendix: Equations and Parameters of the Models

Typically models simulating the evolution of populations of animals are based on dynamical systems. In such a dynamical system a state vector, x, represent the fractions of a population. If the population is divided into n unique states, we can view x as lying in a phase space \({\varvec{x}}\in X\subset {R}^{n}\), and the population can be seen as evolving according to a map, \(\varphi\):

Such a mapping can either be discrete (i.e., \(t\in {\mathbb{Z}}\)) or continuous (i.e., \(t\in {\mathbb{R}}\)). Additionally, it may or may not be smooth over time. In other words, it is conceivable that the population will evolve in different ways depending on separate periods or “seasons”. In our model, we simulate a mating season alternating with a growth season. The evolution of the state during the growth season will depend on how the state evolved during the previous mating season.

The model we used is a discrete, piecewise-smooth dynamical systems. Each timestep represents 1 day, and at each step the population is updated according to a time dependent mapping. With \(\mu\) as a set of parameters, the master equation of our model is:

Separating the growth and mating season in our model is implemented as:

where P and Q are matrices, each depending on a set of parameters, \(\upmu\). Q also depends on the values taken on by the state vector, x, during the previous mating season. The evolution of a state vector according to a transition matrix such as this is generally referred to as a Markov chain.

The parameter range chosen for figures is intended to demonstrate the full portrait of behavior, not guarantee convergence. In general, it is a classic result from probability theory that Markhov chain converges to a unique stationary distribution (equilibrium population) if the transition matrix is irreducible (Allen 2010); in this case, that transition matrix would be P⁶⁰ + Q³⁰⁰. However, there is a time delay in our model (to let prior mating influence future birth), so instead we demonstrate convergence with longtime simulations.

Figures including results iterate this by applying \(P\) 60 times and then Q 300 times for a total of 50 cycles. Figure 10 demonstrates that 50 years is sufficient for producing the same behavior as is present in 100 year and 1000 year simulations.

Comparison of the state of the 2-dimensional parameter sweeps after 10, 50 and 1000 simulated years

Full size image

Base Model: Wrasses

To model the alternative male mating strategies for wrasses, we used a piecewise system of Markov chains (Fig. 1). On any given day, a fish can move from state x to state y with probability P_x→y. A one by seven state vector x[t] is initialized as follows:

The population at t = 0 is assumed to be uniformly distributed across the high and low energy states for the two mating strategies (dominant/sneaker). The population is simulated with two transition matrices:

P is applied to x[t] for the first 60 days (during the mating season), and then Q is applied to x[t] for the following 300 days (during the growth season).

The values of x[t][i] take on real numbers between zero and one representing the proportion of the population that any given state represents. Hence the following holds true for all t:

During the mating season:

During the growth season:

The element in the location (i,j) of a transition matrix is the probability that a subject in state i will transition towards state j in 1 day. The \({i}^{th}\) state is defined as the state corresponding to the \({i}^{th}\) entry of x[t]. Observe that the placement of zeros is consistent with the diagram given by Fig. 1. Making choices for these parameter values requires a consideration of both (1) how the previous mating season played out and (2) to what extent predation and food are readily present in the environment.

To incorporate (1) into the model, we utilize \({Q}_{d\to ls}\) and \({Q}_{d\to ld}\). More accurately, we devised a deterministic method for pulling subsets of the population out of the dead state and into the two low energy states. This re-drawing of individuals from a pool of the dead is a mathematical approach reflects the birth of new fish during the growth season resulting from mating in previous months, devoid of biological interpretation. We used the following method for choosing \({Q}_{d\to ls}\) and \({Q}_{d\to ld}\) at the beginning of the growth season for the \({k}^{th}\) year of iteration:

The likelihood that a new sneaker is born depends explicitly on the relative fraction of matings by sneakers in the previous cycle relative to dominant males. \(\updelta \in \left(\text{0,1}\right)\) reflects volatility in the population, or the proportion of the dead pool distributed towards sneakers and dominants during the subsequent growth season. We choose \(\delta = 0.2\).

The remaining transition probabilities were chosen according to the following relationships:

Here z corresponds to an arbitrary state. More specifically, let \(\alpha\) correspond to food availability and \(\beta\) predation. The transition probabilities are as follows:

Low Energy Sneaker Transitions:

High Energy Sneaker Transitions:

Low Energy Dominant Transitions:

High Energy Dominant Transitions:

For the two-dimensional parameter sweeps (Figs. 4 and 5), \(\alpha\) was swept from 1.20 to 2.00, \(\beta\) from 10⁻⁵ to 0.02. This choice magnifies the curve of change in ESS. These have been re-scaled to a 1–100 scale for the purpose of figures.

The specification of these proportionality constants was done with the following principles in mind:

• The sneaker male is more likely to increase its energy level in a given day than the dominant male.

• Though the dominant male is slower to evolve, it is less likely to subsequently loose energy.

• The sneaker male is more likely to die than the dominant.

• The sneaker male is less likely to transition from the high energy state to the mating state than the dominant male.

Continuing to fill in our matrix, we must consider the values along the diagonals, which describe the probabilities that a subject does not change its state. There is no immediately obvious biological intuition for these values. The columns of each transition matrix need to sum to one. Hence, we will choose each x and y such that this property is necessarily true.

The last parameter to determine is \({Q}_{hs\to ld}\): the frequency for which a male mating strategy change occurs. For this model, we take that to be a constant. The mating strategy change is known to depend on environmental variables in some species, however this dependence is not included here.

In the simulations the updates of the values of the matrix were iterated for 50 years. Note that a simulated year does not correspond directly to a chronological year, since the simulations were set up to reach an equilibrium at the fastest possible rate and aim to reproduce the equilibrium, but not the path to equilibrium of biological evolution. In Fig. 10 we compare the state of the simulations after 10, 50 and 1000 simulated years, and observe that the sharpness of the transition between sneakers and dominant males, but not the qualitative shape of the result changes.

Cephalopods

We model cuttlefish differently from wrasses by instituting four characteristic changes. First, the sneaker to dominant transition is removed as cuttlefish are not observed to change strategies. Next, cuttlefish are observed to be semelparous, so the state of the population undergoes a mass death directly following every mating season. This is executed by resetting the state vector, X[t], to be heavily biased towards the dead state once per year. The model does not reset though, given that the information from the previous mating cycle will still be encoded in the parameters \({Q}_{d\to ls}\) and \({Q}_{d\to ld}.\)

Thirdly, a cannibalism mechanism is introduced. This is executed by letting \({Q}_{ls\to d}\) and \({Q}_{hs\to d}\) be dependent on the total population of cuttlefish playing the dominant strategy at any given time. More specifically, at each time step, t, these two transition probabilities are redefined in the following way:

Additionally, the values of \({y}_{1}\) and \({y}_{2}\) are updated appropriately at each time so that \(\left(\star \right)\) holds. Finally, the growth rates of cuttlefish are higher than those of fishes. Subsequently, we increased the values of \({Q}_{ls\to hs}\) and \({Q}_{ld\to hd}\) by a factor of two.

Figure 3 includes parameter sweeps over varied predation that result from each of these four changes applied individually to the basal wrasse model. The normalized food availability of \(\alpha = 1.52\).

Lamprologus callipterus

Lastly, we shall consider the case of L. callipterus. A third mating strategy, the dwarf male, is introduced as a genetically determined sneaker. We do this by introducing three new dimensions for our state vector and transition matrices:

And during the mating season:

During the growth season:

We also included a link between the dominant and dwarf strategies. Specifically, the dwarf males are assumed to be genetically related to the dominant males. This leads to an incentive for altruistic behavior for the dominant males.

To incorporate this, we introduced Hamilton’s rule by assuming a mean relatedness between any two male members of the breeding population. The dominant Lamprologus individuals receive a benefit for their altruism that is proportional to the total number of dwarf male in the population. This benefit will materialize as an increase to their probability of transitioning to the high energy state. At every time interval we performed the following update:

The remaining transition probabilities are given as follows.

Low Energy Dwarf Transitions:

High Energy Dwarf Transitions:

Given these values, one could fill in the transition matrices P and Q and simulate the population. Figure 7 includes the resulting ESS incorporating parameter sweeps over predator presence and food availability (with and without the altruism bonus). In 7a and 7b, predation is fixed at \(\beta = 0.001\); in 7c and 7d, food availability is fixed at \(\alpha = 1.0.\) In both cases, the other parameter was swept over the same aforementioned ranged as for the wrasses.

Stochastic Considerations

We recognize that food availability and predator presence are not strictly fixed parameters in practice, but instead quantities that fluctuate with time around some initialized state. This is realized in our model but introducing stochasticity in the low to high and high to low energy transitions.

Consider the baseline wrasse model. We can introduce multiplicative noise with a stochastic variable, \(\xi\), which acts as a multiplier on the low to high energy transitions and a divider on the high to low transitions. \(\xi\)> 1 corresponds to a relatively high concentration of food and relatively low predator presence, specifically due to environmental variability. \(\xi\)< 1 is interpreted as the natural inverse. \(\xi\) is chosen via the following mechanism: at the beginning of a new season (mating or growth), reset \(\xi\)= 1. For each following day, update the variable according to the following equation:

That is, a new value of \(\xi\) is chosen from a random normal distribution with mean \(\xi\) and variance \({\sigma }^{2}\). In other words, \(\xi\) takes a random walk during each season before being reset. During a mating season, transition probabilities are updated in the following way each day:

Similarly, during the growth season:

We can then run simulations to observe that the sneaker or dominant strategy may be present in parameter regimes despite being unstable in the deterministic sense.

Springer Nature or its licensor (e.g. a society or other partner) holds exclusive rights to this article under a publishing agreement with the author(s) or other rightsholder(s); author self-archiving of the accepted manuscript version of this article is solely governed by the terms of such publishing agreement and applicable law.

Reprints and permissions