Single or multiple spawning? Comparison of breeding strategies of freshwater Unionidae mussels under stochastic environmental conditions

The life cycle of unionids is characterized by a obligatory period of larval parasitism on a fish host, any disturbance of which might cause a large-scale decline in the mussel population. Because the probability of fish infestation is so important, we modelled what would happen to a population (in terms of population growth and probability of extinction), if the same number of glochidia were released in one or more separate spawning events, by a hypothetical mussel population living in conditions differing in the probability of fish infestation (a “neutral” scenario, a “good” one and three variants of a “bad” one). The single brood strategy was the best in the “good” scenario. However, when the frequency of unfavourable stochastic events increased (“bad” scenarios), all strategies led to population decline, the single brood strategy being the worst. In “good” and moderately “bad” conditions the double brood strategy performed better than the other multiple brood strategies, but as infestation conditions deteriorated, a greater number of spawning events ensured a slower population decline and longer persistence. Our model can facilitate a better understanding of this problem and set up a framework for further tests in other unionid species and their environmental conditions.


Introduction
Important components of aquatic ecosystems (Vaughn & Hakenkamp, 2001;Gutiérrez et al., 2003), freshwater mussels (Unionida) are among the most threatened animal groups on the planet (Lopes-Lima et al., 2017;. The life cycle of unionids is characterized by a unique obligatory period of larval parasitism on a fish host (Kat, 1984;Haag, 2012;Modesto et al., 2018), which is critical for successful reproduction (Strayer, 2008;Brodie et al., 2014); this aspect of their biology remains poorly studied (Ferreira-Rodríguez et al., 2019). The first question that comes to mind is whether the enormous numbers of larvae produced by mussels (from 2,000 to more than 10,000,000 glochidia; Bauer, 1987;Haag, 2013) are released during one event or over longer periods. In the former case, a hypothetical mussel significantly increases the density of larvae in the water and, all other things being equal, the probability of infestation will increase. In the latter case, although the multiplication of breeding events reduces the density of larvae released at a given moment, it spreads the risk of complete breeding failure caused by environmental stochasticity (Goodman, 1984). Continuous release is unlikely, because the density of the larvae would then be extremely low.
Most unionid species are reported to produce only a single brood per year, although the production of multiple broods has been described in North America for the genera Utterbackia, Glebula, Popenaias (Haag, 2012 and citations within) and Elliptio (one to three broods; Price & Eads, 2011; one to five broods in Elliptio steinstansana ;Haag, 2012). In Europe, multiple broods has been described for Unio crassus (up to five broods; Hochwald, 2001; up to seven broods; Zając & Zając, 2018), U. pictorum and U. tumidus (up to two broods in both species ;Piechocki, 1999).
In this study we present a model of a hypothetical unionid species differing in its glochidia-releasing strategies, which faces stochastic variation in the likelihood of fish host presence/infestation. On the basis of this model, we simulated 100 consecutive breeding seasons of this hypothetical mussel population, capable of generating from one to five broods in each season, in unpredictable environmental conditions (affecting infestation probability). This study aims to evaluate how the number of spawning events affects population parameters important from the point of view of population persistence, i.e. its size, extinction rate, coefficient of variation of mean population size and the mean time to extinction in an environment varying in the probability of an event disturbing mussel reproduction.

Model description
We constructed a population model based on data relating to Anodonta cygnea, collected by Zając (2001) in Lake Zalew Pińczowski, Poland (for a detailed description of the study area, see Zając et al., 2016).
Let us suppose that the modelled mussel population consists of I 2 Z 1 age classes. Let g be the proportion of glochidia which successfully attach to the host fish, r the glochidial survival rate to the end of the first year, s i the survival rate of individuals at age i, A i,j the number of individuals at age i during season j and F i,j the number of females at age i during season j. The number of glochidia G i released by each female in age class i is given by the equation where a and b are fixed parameters. The number of recruits R j by the end of season j is given by The number of adult individuals N j by the end of season j is given by The model's validation and parameterization are given in the Supplementary Materials.

Single versus multiple brood strategy simulations
Let us assume that S spawning events take place during season j, where S 2 {1,2,…,5}. Let us further assume that during each season, one of three events influencing the proportion of glochidia attached to the host fish (g) occurs during the whole season: (1) a ''neutral'' one, when the disturbing event does not affect the mussels' breeding; (2) a ''good'' one, when the disturbing event increases the proportion of attached glochidia (e.g. river not in spate) and (3) a ''bad'' one, when the disturbing event decreases the proportion of attached glochidia (e.g. river in spate, suitable host fish individuals absent). Let us also assume that n spawning events are disturbed in a given season, where n B S. Thus, (1) in the case of a ''neutral'' event, g is set at its default value, (2) in the case of a ''good'' event, g is multiplied by 1 ? n S (e.g. if 2 of 5 spawning events are affected by a ''good'' event, g is multiplied by 1.4), or (3) in the case of a bad event, g is multiplied by SÀn S (e.g. if 2 of 5 spawning events are affected by a ''bad'' event, g is multiplied by 0.6).
Let P b , P n and P g be the respective probabilities of occurrence of a ''bad'', ''neutral'' and ''good''event during season j. We simulated a total of five scenarios: (1) a ''neutral'' scenario, where P b = 1/3, P n = 1/3, P g = 1/3; (2) a ''good'' scenario, where P b = 1/4, P n = 1/4, P g = 1/2 and three alternative ''bad'' scenarios: (3a), where P b = 1/2, P n = 1/4, P g = 1/4; (3b), where P b = 3/5, P n = 1/5, P g = 1/5 and (3c), where P b = 2/3, P n = 1/6, P g = 1/6. Simulations were run for 100 consecutive seasons. In each scenario, we modelled five strategies differing in the number of spawning events: the Single Brood Strategy (SBS, 1 spawning event during each season) and four variants of the Multiple Broods Strategy (MBS, from 2 to 5 spawning events during each season). Both events influencing glochidial attachment and the number of spawning events affected by the stochastic event are randomly selected at the beginning of each simulated season. We obtained 1000 realizations of each scenario for each modelled strategy. Below, we present an example of the random selection of an event influencing glochidial infestation and the number of spawning events affected.
To show up the differences between strategies, we compared the final population size (at t = 100) for each strategy in each scenario. Also, for each strategy in each scenario, we obtained the extinction rate (E r )-the percentage of simulation realizations resulting in the extinction of the modelled population, the mean annual growth rate of the population (k) and the mean time to extinction (E t ), calculated using only realizations resulting in the extinction of the modelled population. We used coefficients of variation of mean population size (CV) to compare the relative variability in the number of individuals between strategies in each scenario.

Example
We exemplify the method of randomly selecting an event influencing glochidial infestation and the number of spawning events affected by the event using the MBS strategy with five spawning events in a ''bad'' scenario (3c). Since in this scenario P b = 2/3, P n = 1/ 6, P g = 1/6 and the number of spawning events n equals 5, at the beginning of each season we randomly draw two numbers x (which will determine the stochastic event) and y (which will determine the number of affected spawning events), where x 2 {1,2,…,6} and y 2 {1,2,…,5}. The selection of a stochastic event occurring during the season is based on the value of x. If x B 4, then a ''bad'' event occurs, if x = 5, then a ''neutral'' event occurs and if x = 6, then a ''good'' event occurs. Let us assume that in 3 consecutive draws the randomly selected values are x 1 = 3, y 1 = 4, x 2 = 1, y 2 = 5 and x 3 = 6, y 3 = 2. This means that during the first simulated season a ''bad'' event occurred, which affected 4 spawning events, during the second season a ''bad'' event also occurred, affecting all five spawning events, while during the third season, a ''good'' event occurred, which affected 2 spawning events.

Statistical analysis
The statistical analyses were conducted using Statistica 13. Because the assumption of homogeneity of variance was not fulfilled, the differences in final population size, as well as the differences in mean annual population growth rates and the differences in the mean time to extinction between the simulated strategies in each scenario were tested using the Kruskal-Wallis test. We tested the differences between extinction rates using the test between proportions. We tested the influence of the probability of occurrence of a ''bad'' event and the number of spawning events on the extinction rate (logit transformed) using a General Linear Model (dependent variable: extinction rate (transformed); continuous predictor: probability of a ''bad'' event occurring; ordinal predictor: number of spawning events).

Single versus multiple broods strategy simulations
After a period of initial fluctuations, all the strategies followed a similar trajectory until 60 seasons, when outcomes of the strategies began to differ from each other to a modest extent. The differences in final population size between the strategies were significant (Table 1; Fig. 1A; Kruskal-Wallis test; H = 28.1; P \ 0.0001), whereas the differences in the mean annual population growth rates between the strategies were not significant (Table 1; Fig. 2A; Kruskal-Wallis test; H = 2.49; P = 0.65). The coefficients of variation were higher in SBS than in MBS (Table 1). None of the realizations of this scenario resulted in the extinction of the modelled population (Table 1), and it was clear that, in the absence of extremely rare coincidences, the population would almost certainly not decrease or become extinct.

''Good'' scenario
In a good environment all the strategies differed distinctly from one another (Fig. 1B). The differences in final population size were significant (Table 1; Fig. 1B; Kruskal-Wallis test; H = 1895.0; P \ 0.0001).The final result of the simulations was inversely proportional to the number of broods: the fewer the broods, the larger the final population size (Fig. 1B). A similar pattern was found in the mean annual population growth rate and the coefficient of variation. The differences in mean annual population growth rates between strategies were significant N t=100 mean final population size after 100 seasons (at t = 100), SD Nt=100 standard deviation of mean final population size, k mean annual population growth rate (expressed as a percentage), CV coefficient of variation of mean population size, E r extinction rate (the percentage of realizations of simulations which resulted in the extinction of the population), E t mean time to extinction (in seasons) calculated for the realizations of simulations which resulted in the extinction of the population (Table 1; Fig. 2A; Kruskal-Wallis test; H = 190.3; P \ 0.0001). The coefficient of variation was the highest in SBS and the lowest in MBS(5), although the difference was small (Table 1). However, even small differences in population growth rate translated into large differences in population size over 100 seasons (Fig. 1B). None of the realizations of this scenario resulted in the extinction of the modelled population (Table 1), and it was clear that, in the absence of extremely rare coincidences, the population would almost certainly not decrease or become extinct.

''Bad'' scenarios
In the ''bad'' scenarios, all the populations declined; nevertheless, the population adopting any variant of MBS reached a much larger population size than the population adopting SBS. SBS decreased the fastest, whereas MBS decreased at a much slower rate ( Fig. 1CDE). An interesting aspect of the bad scenario 3a was that strategy MBS(2) decreased more slowly than all the other MBS variants; in scenarios 3b and 3c, the differences between them disappeared. In all the bad scenarios and in all strategies, the population size decreased more rapidly with increasing probability of bad events (the slowest in 3a, accelerating in 3b and   (2), where SBS declined more rapidly than MBS(2) (Figs. 1CDE, 2A). Nonetheless, it is worth noting that: (1) the mean annual growth rate did not change much along with the increasing number of broods, (2) strategies MBS(3-5) became slightly worse than MBS(2) in scenario 3a; they did not change in scenario 3b, but increased in scenario 3c together with the number of broods, where MBS(2) lost its predominance and MBS(5) became the best. Even so, these changes were minimal in comparison to the change from SBS to MBS(2).
The differences in the mean time to extinction between strategies were significant in each modelled scenario (Table 1 Table 1).
Coefficients of variation were higher in SBS compared to MBS (Table 1). Also, the extinction rate was significantly higher in SBS than in any variant of MBS in every simulated bad scenario (Table 1; test between proportions; 3a: P \ 0.0001; 3b: P \ 0.0001; 3c: P \ 0.0001).
The extinction rate depended on the probability of occurrence of a ''bad'' event. If P b took a value \ 0.3, then E r = 0 in each strategy in each scenario; but if P b [ 0.8, then E r = 1 in each strategy in each scenario (Fig. 2C). GLM analysis performed for values of P b ranging from 0.3 to 0.8 showed that both the number of spawning events (N Se ) and the probability of a ''bad'' event(P b ) had a significant influence on the extinction rate ( Fig. 2C; GLM; P b , F = 181.9, P \ 0.0001; N Se , F = 13.9, P = 0.0009). The greater the probability of a ''bad'' event occurring, the higher the extinction rate, but the more spawning events during the season, the lower the extinction rate (GLM; P b , b = 0.90; N Se , b = -0.25).

Discussion
The natural way of developing a theory is to validate it on the basis of falsifiable predictions, which either confirm the model or suggest that it should be rejected or modified. One can predict that, other things being equal, species adopting SBS should perform much better in good conditions than any of the other strategies; MBS2 is second-best. However, SBS and MBS(2) should also disappear much faster than MBS(3-5) in a ''bad'' environment. It turns out that in practical conservation actions, it should be easier and faster to restore SBS species, whereas the restoration of MBS species may be slower.
In general, species living in ''bad'' environmental conditions (in terms of host fish infestation) should be characterized by MBS, whereas species living in ''good'' conditions should adopt SBS. Nonetheless, since MBS guarantees a less variable population size, the two-brood strategy should be the most effective trade-off between SBS and MBS, and this number of broods should be the most common among species that are capable of producing more than one brood (e.g. Utterbackia, Glebula, Popenaias, Elliptio, Unio). Some problems in testing the predictions of this model may be related to the boundaries between SBS and MBS, which do not have to be clear-cut. Species can exhibit a more ''mixed'' strategy than those assumed in this model, like investing more in one spawning event, or saving some energy for later, smaller broods.
By having the capacity to generate from one to seven broods, it seems that Unio species may be a good object for studying the factors regulating the occurrence of MBS. A real practical problem arising when we wish to confirm these predictions may derive from the fact that SBSvsMBS strategies are poorly understood in freshwater mussels, because we know hardly anything about their reproductive effort. Hochwald (2001) points out that because temperature affects the mussels' metabolism, and hence, also their body growth constant and life span, it is quite likely that the number of spawning repetitions is a trait that varies passively in response to environmental factors (e.g. food, temperature). Also, Haag (2012) suggests that glochidial production is determined mainly by physical and energetic constraints. The possibility of multiple clutches in Unionidae mussels is still a matter of discussion. Haag & Staton (2003) suggest that there is no evidence of multiple clutches in Unionidae. On the other hand, Piechocki (1999) and Hochwald (2001) show that Unio species have the ability to produce up to 5 broods per year, whereas Zając & Zając (2018) report 5 broods in U. crassus per season in the River San and up to as many as 7 broods per season in the River Biała. Other researchers suggest that multiple peaks in the numbers of glochidia released are attributable to variations in physical conditions stimulating glochidia release events rather than the production of multiple broods (Bruenderman & Neves, 1993;Hove & Neves, 1994) and that only a small percentage of the reproductive complement is released during such events (Haag & Warren, 2000).
One important question that remains to be answered is why an organism should not breed as often as resources allow? We agree that it should, but the consequences will be the same as in the case of some genetically fixed strategies. If the mussel under consideration were an ''income breeder'' (sensu Stearns, 1992), having no large energy reserves and only accumulating everyday low level energy income during short periods of time to produce many but small broods, the situation of this species would be similar to that of an MBS species. In contrast, a ''capital breeder'', having stored energy, e.g. from the previous season, would be able to invest it all in one big spawning event, achieving a high density of larvae and a high infestation rate. The capital breeder would be much more effective in good conditions, but in the face of unpredictable infestation opportunities, the income breeder (analogous to MBS) would perform better, not only because it was better adapted to poorer feeding conditions, but because it could deal better with unpredictable infestation by spreading the risk of failure over a longer period of time.
In this study we present a general mathematical framework of the consequences of adopting a certain reproductive strategy under stochastic environmental conditions. However, we are aware that life history traits in Unionidae are highly variable within and among species, and our model must be parameterized and corroborated/falsified using real data. In addition, we have assumed that brood size is the same in each spawning event in MBS, and that negative population trends in MBS may be compensated for by the increased number of broods (e.g.the mean number of spawning events in U. crassus in the River Biała is three, but the maximum can reach as many as seven per season; Zając & Zając, 2018). Nevertheless, the total reproductive effort in all strategies in the model must be assumed equal; otherwise, MBS would always result in greater population growth owing to increased reproductive output.
Very little or nothing in known about energy allocation in reproduction and related trade-offs (number vs size of the offspring, parent growth vs reproduction, etc.). Also, we were unable to assess the relative costs (energetic or otherwise) of single versus multiple brood production. Thus, since we lack the required information to make the model more realistic and informative, we were forced to keep it simple. Even so, the results obtained with our model, which formulates the problem explicitly, identifies knowledge gaps and addresses some hither to unidentified questions, shows that MBS multiple brooding seems to be marginally better in consistently bad conditions, when a population is heading towards extinction. Is this negative trend compensated for by a larger total annual glochidial output in MBS compared to SBS? If so, is reproduction based on stored reserves (allowing for SBS) or current income (forcing MBS)? What is the actual level of environmentally induced variation of infestation success and what level of it would give MBS an advantage over SBS? Field studies answering these questions, though highly desirable (Ferreira-Rodríguez et al., 2019), appear very difficult to carry out. The derivation of explicit mathematical models, validated using field data, can at least clarify the interplay of factors influencing the demography of freshwater mussels, in consequence leading to a better formulation of the problems to be solved and a better understanding of the mussels' ecology. provided you give appropriate credit to the original author(s) and the source, provide a link to the Creative Commons license, and indicate if changes were made.