Context dependency of maintenance communities of invasive parasites under climate change: a case study of mussels and intestinal copepods in the Wadden Sea

TL;DR

This study explores how climate change affects invasive and native parasites in mussels, showing that rising temperatures can alter parasite persistence and host population dynamics in the Wadden Sea.

Contribution

The study introduces a mathematical model to show how temperature changes affect parasite-host dynamics and maintenance communities in a real-world marine system.

Findings

Temperature increase reduces mussel populations below the critical size for M. intestinalis persistence.

M. orientalis can persist without oysters under higher temperatures.

M. orientalis does not further reduce host populations when co-occurring with M. intestinalis.

Abstract

Climate change can impact the persistence of native and invasive parasites and their effects on hosts. Given the complexity of interactions in natural systems, models based on parasite–host systems can be helpful to explore long-term impacts. We investigate how two intestinal parasitic copepods impact host populations, and how the predicted temperature increase by year 2100 may affect the persistence and impacts of the parasites. We study Mytilicola intestinalis (a specialist established in blue mussels, Mytilus edulis) and Mytilicola orientalis (a recent invader infecting mussels and Pacific oysters, Magallana gigas) in the Wadden Sea. The parasites are non-lethal but can influence host maturation and fecundity. Using a mathematical model parametrized with empirical, field and literature data, we explore how temperature increase affects parasite basic reproduction numbers and the long-term population trends of parasites and mussels. Temperature increase reduces mussel populations below the critical community size for M. intestinalis persistence, while allowing M. orientalis to persist without oysters. M. orientalis does not have a negative effect on the host population in additional to that of M. intestinalis when both are present. We show that environmental change can have qualitatively different effects on related parasites by changing the role of the shared host as a maintenance population.

Figures

Structure of the models of the life cycles of (a) mussels, (b) M. intestinalis and (c) M. orientalis. Solid lines indicate maturation of hosts and input of free-living stages, dashed lines indicate reproduction of mussels and parasite stage transitions, dotted lines indicate loss of mussels or parasites.

Model-derived steady states of the parasite populations at 10∘C as functions of Λ, the rate of immigration of M. orientalis. Magenta (PU) and red (PW): mean numbers of M. orientalis in immature and mature mussels, respectively; cyan (PQ) and blue (PS): mean numbers of M. intestinalis in immature and mature mussels, respectively. Broken lines show M. orientalis values in the absence of M. intestinalis. Filled circles on the vertical axis show M. intestinalis values in the absence of M. orientalis. Parameter values as in tables 1–3 with (a) ξ=0.0025d−1, (b) ξ=0.0030d−1.

As in figure 2, but here the values of the steady states of the parasite populations are calculated when ξ=0.0014d−1 again as functions of Λ, the rate of immigration of M. orientalis. Magenta (PU) and red (PW): mean numbers of M. orientalis in immature and mature mussels, respectively. Parameter values as in tables 1–3 with (a) temperature 10∘C, (b) 15∘C.

On the vertical axis are the values for the mussel steady states (I+M, total immature and mature) at indicated water temperatures and maturity rates ξ, relative to the baseline value of an uninfected population at 10∘C and ξ=0.003d−1. Values are `0' (green): infection free; `in' (blue): with M. intestinalis only present in the population; `or0' (orange): with M. orientalis only and Λ=0; `or' (red): with M. orientalis only and Λ=0.25ℓ.d−1; `io' (purple): with both parasites present in the population and Λ=0.25ℓ.d−1. A striped entry indicates Rin<1, Ror<1, or both as appropriate, meaning that the parasite is not present in the steady state and thus does not affect the mussel population.

Graphic summary of our results on the maintenance of M. intestinalis and M. orientalis at the two temperature scenarios. The terminology for maintenance is based on Haydon et al. [8], see also [1]. We depict the situation in which M. orientalis larval immigration from oysters is always present and fixed at a low rate (Λ=0.03 d−1), and the mussel maturation rate is moderate to high (0.0025d−1<ξ<0.003d−1). We focus on the following situations: whether one or both parasite species are present and maintained and how this differs by temperature (10 and 15∘C). Mussels are a maintenance host for M. intestinalis at 10∘C, but not at 15∘C, either alone or in the presence of M. orientalis (first and third columns). The maintenance community of M. orientalis depends on temperature and whether M. intestinalis is present. Mussels and oysters are a maintenance community for M. orientalis at both temperatures when M. intestinalis is absent (second column). When M. intestinalis is present (fourth column), at 10∘C larval immigration from oysters (arrow) is needed to sustain M. orientalis in mussels (hence without oysters, M. orientalis is not maintained in the system). At 15∘C, mussels and oysters are a maintenance community for M. orientalis. In addition, the model results suggest that as a single species, mussels are a maintenance host for M. orientalis at 15∘C (i.e. without the influx of larvae from an oyster population).

Tables

Table 1.. Mussel parameter values used in the model.

at 10∘C
birth rate of hosts	ν=0.03d−1
minimum death rate of immature hosts	μI=4.0×10−3d−1
minimum death rate of mature hosts	μM=2.0×10−3d−1

Table 2.. Mytilicola intestinalis parameter values used in the model.

at 10∘C
transmission rate to immature hosts	βQ=2.2×10−4p.ℓ.−1m.−1d−1
transmission rate to mature hosts	βS=4.6×10−4p.ℓ.−1m.−1d−1
parasite mortality in immature hosts	σQ=7.8×10−3d−1
parasite mortality in mature hosts	σS=7.8×10−3d−1
parasite transmission from immature hosts	λQ=0.38ℓ.p.−1d−1
parasite transmission from mature hosts	λS=0.38ℓ.p.−1d−1
loss rate of free-living stages	ζS=0.29d−1

Table 3.. Mytilicola orientalis parameter values used in the model.

at 10∘C
transmission rate to immature hosts	βU=1.7×10−4p.ℓ.−1m.−1d−1
transmission rate to mature hosts	βW=3.7×10−4p.ℓ.−1m.−1d−1
parasite mortality in immature hosts	σU=7.8×10−3d−1
parasite mortality in mature hosts	σW=7.8×10−3d−1
parasite transmission from immature hosts	λU=0.41ℓ.p.−1d−1
parasite transmission from mature hosts	λW=0.41ℓ.p.−1d−1
loss rate of free-living stages	ζW=0.25d−1

Table 4.. Basic reproduction numbers R0 of the two parasite species at different values of the water temperature and rates of host maturation ξ. Values are Rin : R0 for M. intestinalis; Ror : R0 for M. orientalis; Rio : if Rin>1 then R0 for M. orientalis when M. intestinalis is at steady state, otherwise Rio=Ror.

at 10∘C	Rin	Ror	Rio
ξ=0.0014d−1	0.867	0.860	0.860
ξ=0.0025d−1	1.014	1.007	0.988
ξ=0.0030d−1	1.054	1.048	0.971

Equations

Full text

Introduction

The ability of a parasite to persist in an ecosystem or to successfully invade new ecosystems depends on the availability of suitable host species and will further depend on the specific environmental context rather than being a fixed characteristic of the host–parasite system [1]. Factors influencing the ultimate outcome include environmental changes that may alter the species composition and population densities of host communities as well as the life cycles of the parasites directly [2,3]. Furthermore, since parasites can affect host populations [4–6], the presence of multiple parasite species may also affect the persistence of parasite species sharing a host species [7]. We will use the term persistence, but also refer to maintenance hosts, maintenance communities and reservoirs as defined in [8] within the boundaries of our study system, in the knowledge that the use of this terminology is ecosystem and context specific [1].

One parasite trait that can influence how environmental change affects their persistence is host specificity. Specialist parasites are dependent on the success of their only host species, while generalists can thrive even if one host species is lost altogether [2]. Other relevant parasite traits are their pre-existing environmental adaptations and the genetic diversity of the population, as a higher variation in the genetic pool increases adaptive potential [2,3]. All of these traits can differ between native and introduced species: typical successful invaders are more generalist on the specialist–generalist continuum [9,10] and can have pre-existing environmental adaptations that do not completely align with the conditions of the introduced range, while species often go through genetic bottlenecks on their invasion pathway, resulting in reduced genetic diversity [3,10]. Introduced parasites can also have pre-existing traits that allow them to better adapt to future changing conditions in their new environment. However, the presence of other parasites and hosts in the system can affect the invasion success of novel parasites, both through direct interactions and through their effects on host populations [2,11]. To predict the effect of environmental changes on parasite invasion success, persistence, dynamics and host population interactions, one thus needs to incorporate the complexity of other relevant hosts and parasites in a specific system. Our current understanding of the interplay of environmental shifts such as climate change and parasite invasions is limited and dependent on the parasite–host system investigated [12]. Given the challenging complexity of natural systems, models based on well-studied parasite–host systems can be helpful to explore the long-term impacts of environmental changes on parasite invasions.

In this study, we investigated the effect of environmental change in the form of increasing temperature on a community consisting of two species of intestinal parasitic copepods with simple life cycles and two partially shared ectothermic host species in a temperate coastal ecosystem. This parasite–host system is particularly suited to investigate climate change effects on parasite invasions as a wealth of observational and experimental data, including on temperature effects, exists from the study area. In addition, the Wadden Sea is subjected to significant warming. The mean annual water temperature was approximately $[eqn]$ before climate change impacts and is predicted to increase to $[eqn]$ – $[eqn]$ by the end of the century [13,14]. The parasitic copepod, Mytilicola intestinalis, is a specialist that only infects native blue mussels (Mytilus edulis) and arrived in the system approximately $[eqn]$ years ago [15], hence it can be considered to be an established species. Mytilicola intestinalis was first detected in the Mediterranean in $[eqn]$ , but the low genetic diversity in its invaded and putative native region may indicate that it has experienced a genetic bottleneck prior to its arrival in the North Sea [15]. A second congener, Mytilicola orientalis, arrived only $[eqn]$ years ago with its original host, the invasive Pacific oyster (Magallana/Crassostrea gigas) from the western Pacific Ocean [15]. The genetic diversity of M. orientalis is high even in the invaded range due to hitchhiking on numerous oyster transports from the native range [15]. It is unclear if oysters are required for the persistence of M. orientalis in the invaded range. What is currently known is that the oysters and parasites arrived together and coexist throughout the invaded range. In the Baltic Sea, where M. orientalis was first detected in $[eqn]$ , this coincided with the first successful spatfall of oysters [16,17], suggesting that oysters need to be present for M. orientalis to establish. Neither parasite species increases host mortality within a temperature range of $[eqn]$ – $[eqn]$ , but M. orientalis causes a slight reduction in mussel shell growth and M. intestinalis a reduction in mussel condition [14,18]. Direct competition between the parasite species does not affect infection rates [19]. An increase in temperature accelerates the development of both parasites and also shortens the lifespan of infectious free-living juvenile stages in most of the temperature range [14]. For M. intestinalis, the decrease in juvenile lifespan occurs across the whole range of temperatures, but for M. orientalis juveniles, the longest lifespan occurs at $[eqn]$ instead of $[eqn]$ before a decline with increasing temperature [14], as a possible sign of inadequate adaptation for current winter temperatures in the Wadden Sea, but simultaneously perhaps aligning better with the predicted future temperature range.

Given the current presence of both parasites in the Wadden Sea, we investigated three questions with the aid of a mathematical model, parametrized with our empirical, field and literature data:

(a) What conditions were required for M. orientalis to invade the system already occupied by M. intestinalis, and are oysters required for the persistence of M. orientalis in the community?(b) Did the addition of the second invasive parasite increase the parasite impact on mussel populations?(c) How is the predicted increase in temperature by year $[eqn]$ likely to affect the maintenance communities of the two parasites and their impact on blue mussel populations?

Methods

The dynamics of the mussel (host) and parasite populations were modelled using a system of ordinary differential equations [20]. In our model, the host population consists of two simplified life stages: immature mussels that grow $[eqn]$ , and mature mussels that reproduce $[eqn]$ . Both stages can get infected with both parasites. Populations of both parasites consist of parasitic adult life stages that live in the intestine of hosts and infectious free-living juvenile stages that hatch from the eggs of the sexually reproducing adult parasites and spread to new hosts via seawater (see figure 1). Steady-state (equilibrium) values of the parasite populations were found using MatLab [21] as functions of the effect on host maturation due to infection with M. orientalis and the rate of immigration of M. orientalis from adjacent oyster populations. The steady states where just one parasite species is present were found using the fzero function; for both species present, fsolve was used. A steady state of M. intestinalis only exists when there is no immigration of M. orientalis. When there is immigration of M. orientalis, then steady states of that parasite and of both parasites together were computed. Stability of the steady states was checked by computing the eigenvalues of the appropriate Jacobian matrices. In addition, the differential equations were solved numerically over time using the ode45 function of MatLab [21], confirming the steady-state values and their stability. The results are presented as mean numbers of M. intestinalis in immature and mature hosts, $[eqn]$ and $[eqn]$ , respectively, and mean numbers of M. orientalis in immature and mature hosts, $[eqn]$ and $[eqn]$ , respectively. Steady states in the absence of the other parasite are indicated with a $[eqn]$ superscript, and in the presence of the other parasite are indicated with a $[eqn]$ superscript.

Structure of the models of the life cycles of (a) mussels, (b) M. intestinalis and (c) M. orientalis. Solid lines indicate maturation of hosts and input of free-living stages, dashed lines indicate reproduction of mussels and parasite stage transitions, dotted lines indicate loss of mussels or parasites.

The model

2.1.

Let the population densities of immature and mature mussels be $[eqn]$ and $[eqn]$ , respectively. Assume that each subpopulation exhibits logistic dynamics, but without competition for space or resources between the immature and mature populations, and that competition for resources within each subpopulation increases mortality. Let the rate at which uninfected mature mussels produce immature offspring be $[eqn]$ , the rate at which uninfected mussels mature be $[eqn]$ and the rates of mortality of immature and mature mussels be $[eqn]$ and $[eqn]$ , respectively. The parameter $[eqn]$ depends on the local predator population, and $[eqn]$ depend on the sea temperature and food availability. The population dynamics of the mussel population are described by the differential equations

$$[eqn]$$

where $[eqn]$ is the mean number of M. intestinalis per mature mussel and $[eqn]$ is the mean number of M. orientalis per immature mussel. The functions $[eqn]$ and $[eqn]$ quantify a negative effect of the within-host parasite load on the maturation rate of immature hosts and offspring production per mature host, respectively. We model these as functions of the mean parasite loads in the two mussel life stages. As infection with M. orientalis delays maturity, we have $[eqn]$ and $[eqn]$ . Similarly, as infection with M. intestinalis reduces fecundity, we need $[eqn]$ and $[eqn]$ . Specific functional forms for $[eqn]$ and $[eqn]$ were determined related to observed parasite load distributions (see §2.2).

Setting $[eqn]$ for now, the steady-state (equilibrium) solutions of equation (2.1) satisfy $[eqn]$ and $[eqn]$ . Let these solutions be $[eqn]$ and $[eqn]$ . Hence, $[eqn]$ solves the cubic equation

$$[eqn]$$

We have $[eqn]$ and $[eqn]$ for $[eqn]$ . Hence, there is a unique solution $[eqn]$ if and only if $[eqn]$ . The parameter values used for modelling the mussel population are summarized in table 1.

Let the total number of M. intestinalis in the immature and mature mussel populations be $[eqn]$ and $[eqn]$ , respectively, and the density of parasite-free-living stages in the water be $[eqn]$ . Hence, the dynamics of the parasite population may be modelled by the differential equations

$$[eqn]$$

The parameter values used for the M. intestinalis population are summarized in table 2.

The next-generation matrix for M. intestinalis is

$$[eqn]$$

where $[eqn]$ and $[eqn]$ (see Diekmann et al. [20] chapter 11). The basic reproduction number ( $[eqn]$ ) for M. intestinalis is the largest eigenvalue of $[eqn]$ ,

$$[eqn]$$

The first term under the square root represents that a parasite in an immature mussel produces free-living stages at the rate $[eqn]$ for time $[eqn]$ , and that free-living stages infect immature mussels at the rate $[eqn]$ for time $[eqn]$ . Hence, it represents the expected number of parasites produced in immature mussels by a single parasite in an immature mussel, without the involvement of mature mussels. Similarly, the second term represents the expected number of parasites produced in mature mussels by a single parasite in a mature mussel, without the involvement of immature mussels. The third term is the proportion of parasites in immature mussels that are still infecting the mussel after it has matured, $[eqn]$ , multiplied by the rate at which parasites in mature mussels produce free-living stages and the time for which they do that, $[eqn]$ , multiplied by the rate at which free-living stages infect immature mussels and the time for which they do that, $[eqn]$ . Hence, it represents the expected number of parasites produced in immature mussels by a single parasite in an immature mussel, following the maturity of the host mussel. Each of these ingredients is added to produce the expected number of secondary parasites due to a primary parasite. The square root commonly arises when characterizing $[eqn]$ for a parasite with a two-stage life cycle [20]. This is because $[eqn]$ describes growth/decline in terms of subsequent generations of infected individuals. If a system has different host stages or species that alternate in generations, such as here and for many host–vector systems, one has to look two generations ahead to characterize growth in the same stage/species. Per-generation growth, hence, comes with a square root. If $[eqn]$ the parasite population cannot persist in the host population, and if $[eqn]$ the parasite can invade and persist [20].

Similarly, let the total number of M. orientalis in the immature and mature mussel populations be $[eqn]$ and $[eqn]$ , respectively, and the density of parasite-free-living stages in the water be $[eqn]$ .

The dynamics of the parasite populations are described by

$$[eqn]$$

The term $[eqn]$ accounts for immigration of M. orientalis from adjacent oyster host populations. As our focus is on the population dynamics of mussels and their parasites, we do not explicitly model M. orientalis in oyster populations. We assume that mussels are infected by M. orientalis free-living stages generated by infected mussels and oysters and present results over a range of $[eqn]$ values reflecting different population densities of oysters and M. orientalis adjacent to mussel populations. The parameter values used for the M. orientalis population are summarized in table 3.

The basic reproduction number for M. orientalis is

$$[eqn]$$

where $[eqn]$ and $[eqn]$ . The value $[eqn]$ is the invasion threshold for M. orientalis in uninfected mussel populations. Where the host population is infected with M. intestinalis at steady state, the threshold is $[eqn]$ . The expression for $[eqn]$ is identical to that for $[eqn]$ , except that $[eqn]$ and $[eqn]$ are replaced with the steady-state population densities of immature and mature mussels, respectively, in the presence of M. intestinalis.

Derivation of parameter values

2.2.

The parametrization of the model was based on a combination of a search of pre-existing literature, the results from our own experiments to quantify the relationship between key parameter values and temperature, and the collection of field data from the southern Wadden Sea in $[eqn]$ – $[eqn]$ . Here, we briefly describe the sources of each parameter value in tables 1–3; for more details, see electronic supplementary material, tables S1–S3. The details of temperature experiments with the two parasites and mussels can be found in [14,18]. The host birth rate $[eqn]$ included the whole process from spawning until reaching a settled juvenile size that can get infected. Therefore, the estimate was derived from reported differences in mussel settlement success in the Wadden Sea based on winter temperatures that affect the presence of mussel predators [22]. Mussel mortality rates $[eqn]$ and $[eqn]$ were based on rates reported in [23], with a 10% increase in mortality at the higher temperature range as reported in [24]. The competition coefficients $[eqn]$ and $[eqn]$ , that increase mortality in our model, included all competitive interactions between hosts. The values for these parameters were approximated based on two field studies [25,26]. To keep the model simple, we decided not to add a threshold level above which competition starts, but to lower the estimate of density-dependent increased mortality derived from [26] by one order of magnitude. Because immature mussels are more affected by predation than mature mussels [27], they are more likely to benefit from the predation protection provided by other mussels. Hence, the competition coefficient of immature mussels was chosen to be half that of mature mussels. The host maturation rate $[eqn]$ was based on results reported in [28]. The effects of temperature on the growth of mussels were obtained from [29]. The four transmission rates from hosts ( $[eqn]$ values) were based on a combination of the number of eggs per parasite [30,31], an estimate that approximately $[eqn]$ of all adult parasites are females co-inhabiting an intestine with a male, and on development rates derived from our experiments [14,32]. The loss rates of free-living stages ( $[eqn]$ values) were also based on results from [14,32]. The transmission rates to hosts ( $[eqn]$ values) could not be derived from field data or experiments, because we could not determine the concentration of free-living stages in water in our system. Therefore, these parameters at the lower (current) temperature range were estimated by balancing the model so that the $[eqn]$ values were consistent with field observations of the order of magnitude of parasite burdens. At the lower temperature range, most parasite loads in the field are in the range $[eqn]$ – $[eqn]$ copepods/host (unpublished field data of $[eqn]$ mussel dissections from Dutch Wadden Sea 2022–2024, electronic supplementary material, figure S1).

It was assumed that the ratio of $[eqn]$ values of immature and mature mussels should be the same as the ratio of their filtration rates [33,34], because these intestinal parasites most likely enter their hosts with food. The change in $[eqn]$ values by temperature was determined from the change in host filtration rates [33,34]. The rates of parasite mortality within hosts ( $[eqn]$ values) were estimated from our experimental data [14,32]. The parasite load distribution from field data was also used to choose the simplest accurate shapes for parasite effect functions, $[eqn]$ and $[eqn]$ , (electronic supplementary material, figure S1, see tables 2 and 3), along with the observation that in our experiment a single parasite caused a host effect with little increase in effect by each additional parasite [18]. There is insufficient knowledge on the biology of the processes involved on which to base a derivation of a suitable function, and we hence chose simple functions that at least mimic part of the intended feedback. These functions were assumed not to change with temperature and to be the same for all scenarios. In our experiment, the effect of $[eqn]$ , relevant for M. orientalis, did not change by temperature, and even though the effect of $[eqn]$ , relevant for M. intestinalis, was clearest at the lower range of temperatures, an earlier experiment had shown a similar effect at a higher temperature [35].

Results

At $[eqn]$ and low values of the maturation rate of mussels, $[eqn]$ , the basic reproduction number of both parasite species is below one, see table 4. At moderate and high values of $[eqn]$ both parasite species have $[eqn]$ , but M. intestinalis outcompetes M. orientalis: $[eqn]$ . This means that when M. intestinalis is present, M. orientalis can only be maintained in the mussel population if there is immigration of parasite free-living stages from adjacent oyster populations. When M. intestinalis is not present, this immigration is not necessary for M. orientalis to persist in the mussel host.

In figure 2, parasite steady states are shown for two values of $[eqn]$ and increasing values of the rate of immigration of M. orientalis, $[eqn]$ . When $[eqn]$ the two parasite species coexist in the mussel populations for $[eqn]$ , but for higher rates of immigration, infection with M. orientalis reduces the size of the mussel population to a level where M. intestinalis is no longer supported. The parasite population dynamics are sensitive to the mussel population densities, going from $[eqn]$ to $[eqn]$ only changes the steady-state values of $[eqn]$ and $[eqn]$ by less than $[eqn]$ (result not shown). In contrast, when $[eqn]$ both parasite species persist over the range of $[eqn]$ considered, see figure 2.

Model-derived steady states of the parasite populations at 10∘C as functions of Λ, the rate of immigration of M. orientalis. Magenta (PU) and red (PW): mean numbers of M. orientalis in immature and mature mussels, respectively; cyan (PQ) and blue (PS): mean numbers of M. intestinalis in immature and mature mussels, respectively. Broken lines show M. orientalis values in the absence of M. intestinalis. Filled circles on the vertical axis show M. intestinalis values in the absence of M. orientalis. Parameter values as in tables 1–3 with (a) ξ=0.0025d−1, (b) ξ=0.0030d−1.

When $[eqn]$ , the lowest value considered, and at both temperatures considered, $[eqn]$ and M. intestinalis cannot persist in the mussel population. Even though $[eqn]$ , M. orientalis can persist at steady state due to parasite immigration, $[eqn]$ . These results are illustrated in figure 3. When the two higher values of $[eqn]$ were considered at $[eqn]$ , we had $[eqn]$ and $[eqn]$ , see table 4. Very high steady-state values of M. orientalis were obtained, with M. intestinalis failing to persist (see electronic supplementary material, figure S2).

As in figure 2, but here the values of the steady states of the parasite populations are calculated when ξ=0.0014d−1 again as functions of Λ, the rate of immigration of M. orientalis. Magenta (PU) and red (PW): mean numbers of M. orientalis in immature and mature mussels, respectively. Parameter values as in tables 1–3 with (a) temperature 10∘C, (b) 15∘C.

Figures 2 and 3 present results at different values of $[eqn]$ , the maturation rate of mussels uninfected with M. orientalis. The effects of maturation rate, parasite infection and water temperature on mussel populations are shown in figure 4. The figure shows that increasing the water temperature to $[eqn]$ reduces the mussel population density to between 50 and $[eqn]$ of its baseline value at $[eqn]$ . At both temperatures, the lowest maturity rate corresponds to a much lower host population density, and infection with either or both parasite species only reduces the mussel population density by between 5 and $[eqn]$ at the two highest maturation rates. When M. intestinalis is present, the addition of M. orientalis does not cause a further reduction in mussel population, and at $[eqn]$ , M. intestinalis is not maintained in the population. Going from $[eqn]$ to $[eqn]$ only changes the range of $[eqn]$ at which M. orientalis is present, not its effect on hosts.

On the vertical axis are the values for the mussel steady states (I+M, total immature and mature) at indicated water temperatures and maturity rates ξ, relative to the baseline value of an uninfected population at 10∘C and ξ=0.003d−1. Values are 0' (green): infection free; in' (blue): with M. intestinalis only present in the population; or0' (orange): with M. orientalis only and Λ=0; or' (red): with M. orientalis only and Λ=0.25ℓ.d−1; `io' (purple): with both parasites present in the population and Λ=0.25ℓ.d−1. A striped entry indicates Rin<1, Ror<1, or both as appropriate, meaning that the parasite is not present in the steady state and thus does not affect the mussel population.

Discussion

In our model, a moderate temperature increase changed the maintenance communities of two closely related parasites in qualitatively different ways. In [1], it was shown that the concepts of maintenance host, maintenance community and reservoir (as proposed in [8]) depend on the ecosystem context of the species that are considered. The same species can be a maintenance host or a reservoir for a parasite in one context, but that capacity can be lost due to changes in one or more non-host species in the ecosystem. In the present analysis, we see that even without considering the wider ecosystem context, the interactions of two parasites with a single host species are important in determining whether that species plays a role in maintenance [1]. Not only direct changes in interacting species, but also abiotic factors that affect the communities change the populations required for parasite persistence. Another factor that affected the ability of a parasite to be maintained in a host species was the presence of another parasite: at the lower temperature, the second parasite required the presence of another host species to act as a reservoir when the first parasite was present, whereas without the first parasite mussels sufficed as a maintenance community (solid red and magenta lines start at zero while broken lines start at non-zero values when immigration of M. orientalis is zero in figure 2). This means that the presence of the first parasite provided resistance against the invasion of the second parasite.

Our model suggests that due to an increase in temperature within the range predicted by year $[eqn]$ and with other factors remaining constant, the established specialist parasite M. intestinalis will lose its ability to persist in the local mussel population, whereas the invading generalist parasite M. orientalis will gain the ability to persist. For that persistence, the original host, the Pacific oyster, is no longer necessary in this system. This transformation is illustrated in figure 5.

Graphic summary of our results on the maintenance of M. intestinalis and M. orientalis at the two temperature scenarios. The terminology for maintenance is based on Haydon et al. [8], see also [1]. We depict the situation in which M. orientalis larval immigration from oysters is always present and fixed at a low rate (Λ=0.03 d−1), and the mussel maturation rate is moderate to high (0.0025d−1<ξ<0.003d−1). We focus on the following situations: whether one or both parasite species are present and maintained and how this differs by temperature (10 and 15∘C). Mussels are a maintenance host for M. intestinalis at 10∘C, but not at 15∘C, either alone or in the presence of M. orientalis (first and third columns). The maintenance community of M. orientalis depends on temperature and whether M. intestinalis is present. Mussels and oysters are a maintenance community for M. orientalis at both temperatures when M. intestinalis is absent (second column). When M. intestinalis is present (fourth column), at 10∘C larval immigration from oysters (arrow) is needed to sustain M. orientalis in mussels (hence without oysters, M. orientalis is not maintained in the system). At 15∘C, mussels and oysters are a maintenance community for M. orientalis. In addition, the model results suggest that as a single species, mussels are a maintenance host for M. orientalis at 15∘C (i.e. without the influx of larvae from an oyster population).

We also showed that adding a second invasive parasite species, one that is harmful to the hosts on the individual level, does not necessarily increase harm on the host population. When M. intestinalis was already present (blue bars in figure 4), the arrival of M. orientalis at the current temperature did not increase harm on the mussel population (purple bars in figure 4), but rather the population density trivially increased (both purple bars are $[eqn]$ higher than the blue bars) when the mussel maturation rate was moderate to high ( $[eqn]$ ), corresponding to $[eqn]$ – $[eqn]$ days from settlement to spawning which is realistic for the Wadden Sea [28]. If M. orientalis was the only parasite in the system, it did decrease the host population (green versus orange and red bars in figure 4), an effect that seems slightly stronger at the elevated temperature. This is reflected in the positive trend in the mussel population observed since the arrival of M. orientalis in the Wadden Sea [36]. This finding fits with the theory that adding biodiversity can balance systems [37], even if the added species are not native. A facilitation between invasive species did occur in the form of oysters enabling the maintenance of M. orientalis at the lower temperature. Still, the effect of this facilitation on the population of native hosts was positive, and oysters could be seen to offer an ecosystem service for mussels by enabling the second parasite to be present. This effect became detrimental only at the higher temperature when M. orientalis became the only parasitic copepod in the system, which was also when oysters were no longer required to maintain their persistence.

This observation adds to the list of benefits mussels gain from Pacific oysters, despite them being an invasive species. Oysters have previously been shown to protect mussels from predation, storms and ice scouring [38,39]. The role of oysters is even more complex, because they also filter out larvae of both themselves and mussels [40] and they compete with mussels for food [38]. In this simplified model, oysters were included only as a source for infectious stages of M. orientalis. It would be interesting to add oysters as a full entity in the model. Our model is already complex, and adding all oyster interactions would have risked losing the insight from the simpler model.

In our model, no direct competition between the two parasites was added, so the negative effect M. orientalis had on M. intestinalis at the lower temperature occurred solely through host effects. In the northern Wadden Sea, M. intestinalis prevalence has decreased since the arrival of oysters with M. orientalis, and that effect has been attributed to M. intestinalis avoiding mussels already infected by M. orientalis, while M. orientalis tends to aggregate in already-infected mussels [19]. The decline in M. intestinalis numbers could also be because of a dilution effect by oysters that are not competent hosts but filter high volumes of water, as they have been previously shown to protect mussels from trematode infections by filtering out the infectious stages [41]. Both the avoidance of hosts already infected with M. orientalis and the dilution effect by oysters may play a role in the decline of M. intestinalis, but our model shows that the minor host effect due to M. orientalis alone suffices to explain a decline in M. intestinalis numbers.

We assumed that oysters suffice for the maintenance of M. orientalis at the pre-climate-change temperature. This implies the assumption that their critical community size is high enough to keep $[eqn]$ . Based on the genetic evidence of multiple introductions of M. orientalis from the native range, where it infects oysters, and the first detection of M. orientalis typically coinciding or taking place after the invasion of oysters [15–17], it seems plausible that already-low initial populations of introduced oysters are able to maintain the parasite. In the current Wadden Sea, oysters dominate most intertidal reefs and are not likely to suffer from the predicted increase in temperature, because the temperature in their original range is already higher than that predicted for the Wadden Sea by $[eqn]$ [39,42].

When considering experimental data on temperature effects on parasite life cycles, the two parasites appear to be strikingly similar, and the number of life cycles they can complete at different temperature scenarios differs little [14]. The main difference in the parameter values of the two parasites was the poorer adaptation of M. orientalis for low temperatures, which is logical considering how water temperature rarely goes below $[eqn]$ in the parasite’s original range [14,42]. In our model, including data on host effects and other host species changed the outcome of temperature increase on the two parasites, indicating that very simplified predictions for one species can be misleading if interacting species in the ecosystem are ignored.

Adaptation to the changing environment through evolution in the next $[eqn]$ years was not included in the model. This was to enable the analysis of the already complex model and because of uncertainty in the speed of warming. The annual mean water temperature in the southern Wadden Sea has stayed above $[eqn]$ for $[eqn]$ – $[eqn]$ , which is $[eqn]$ higher than the annual mean of $[eqn]$ – $[eqn]$ and $[eqn]$ higher than that of $[eqn]$ – $[eqn]$ , making it possible that the existing predictions severely underestimate the speed of warming [43]. Furthermore, the warming might not be spread evenly across the seasons in this region but rather be more pronounced in the summer [44]. These factors can limit the ability of adaptation to keep up, while evolution still increases the uncertainty of our predictions. The adaptive capacity of the two parasite species is likely to differ because of their large differences in genetic diversity; the chances for M. intestinalis to adapt to the fast change in temperature do not appear promising due to it having passed through a genetic bottleneck on its invasion pathway. In contrast, the high genetic diversity and environmental tolerance of Wadden Sea M. orientalis makes it more likely to adapt to the changing environment [15,16]. The maintenance of M. intestinalis may in future instead be facilitated by distributional shifts of a congeneric species, the southern Mediterranean mussel (Mytilus galloprovincialis), which is adapted to higher temperatures, making it possible that the blue mussels of the southern Wadden Sea will be replaced and hybridized with this species rather than merely decline [45]. The drop in the critical community size of M. edulis to below the threshold required for the maintenance of M. intestinalis is less likely to happen with temperature increase if the species is hybridized with or replaced by M. galloprovincialis. The continued presence of M. intestinalis in M. galloprovincialis populations in southern sites of the European Atlantic decades after the first occurrence of M. intestinalis supports this possibility [15]. This underlines the possibility that species increasing their range with or without the help of humans is not necessarily detrimental but can also increase the resilience of ecosystems.

Our model has $[eqn]$ parameters, $[eqn]$ of which change with temperature. All parameter values, except the transmission rates to hosts, $[eqn]$ , were based on best available information for the system, and extensive effort was put into getting them as accurate as possible. Only the final mortality rate of free-living M. orientalis free-living stages, $[eqn]$ at $[eqn]$ ( $[eqn]$ ), differed from the value derived from experimental results. This was for two reasons: when we analysed the model with the value based only on larval lifespan in the experimental data at $[eqn]$ ( $[eqn]$ ), M. orientalis loads per mussel were in the hundreds when $[eqn]$ , which is not realistic for this parasite based on field data [15,19]. This value also resulted in slightly higher effect on host populations at $[eqn]$ , but it did not affect any other outcomes of our model. The second reason was that in real ecosystems, the mortality of free-living stages is determined by their innate lifespan, predation, getting washed away by tides and salinity [16,46], and temperature increase is likely to increase predation [46]. Increasing the mortality rate of M. intestinalis free-living stages at $[eqn]$ to reflect this increase in predation did not affect the model outcome, as already the low estimate resulted in elimination of the parasite population, but increasing the mortality rate of M. orientalis free-living stages increased the realism of the model output. The final $[eqn]$ values still reflect the longer lifespan of M. orientalis free-living stages at $[eqn]$ compared with $[eqn]$ , and compared with M. intestinalis free-living stages at $[eqn]$ , as seen in the experiment [32]. The effect of temperature change on mussel maturation is not straightforward, because even though temperature may boost growth by increasing ectotherm metabolic rates when no food limitation is present [47], food availability is affected by several environmental variables, including temperature and even the location of mussels in a mixed oyster–mussel reef [38,48]. Also, the genetics of the particular mussel population can modify the effect of environmental variables on growth [29]. For this reason, we decided to analyse the model with three levels of mussel maturation at both temperatures, the most appropriate range being determined by the system and situation of interest.

In our conceptual approach, we have ignored many dimensions of heterogeneity in the real system influencing the mussel–oyster–parasite interaction that we have studied in isolation. Hence, our aim is decidedly not to predict the time evolution of population abundance of mussels and their parasites. Our aims and conclusions are more qualitative and potentially more robust than predicting spatio-temporal variation in detail. In our isolated system, everything is kept constant apart from our parameters of interest, and many aspects are not considered at all. For example, we have assumed that the temperatures are constant, whereas in reality they vary with time (over the days, weeks and months) and with space at any given moment (with shallow parts experiencing different regimes and extremes). On top of that, there are many factors that we ignore such as spatial heterogeneity, tidal movement, predators, currents, food dynamics and competitors. One important competitor of the mussel is the oyster. We have also chosen to include the oyster population as a constant and have explored the influence of oyster population through its smaller or larger (constant) production of infectious stages of the M. orientalis. This means that we have ignored the population dynamics of oysters and any other feedback to and from the mussel population (for example, through competition for food, but also through providing shelter for mussels against their predators). In addition, we only consider the potential effect of a long-term rise in temperature, whereas there are other effects of climate change in systems like the Wadden Sea that may be relevant and impact the mussel–oyster–parasite system. One can think of sea level rise, changes in salinity, future invasions and abundance, composition and distributional changes in the ecosystem (see for example [49]).

In our system of interest, a small change in the average environmental temperature sufficed to cause a change in the maintenance communities and host population effects of two invasive parasites. Mussels did not act as maintenance hosts for M. orientalis at the pre-climate-change temperatures of the Wadden Sea because of the parasite’s poor adaptation for cold temperature and the presence of another parasite, M. intestinalis. This adaptation for higher temperatures proved beneficial in the higher temperature range predicted to occur due to climate change by $[eqn]$ , and being a generalist that had invaded the system with its original host enabled M. orientalis to persist also at the lower temperature. The future of the specialist M. intestinalis in this system appears less promising because of its dependence on one host species that is predicted to suffer from the increase in temperature. These findings highlight how small changes in the environment can change parasite persistence and host effects, and the importance of including the environmental change effects of not only the species of interest but also the species it interacts with in predictions.