Impact of cyber-invasive species on a large ecological network

As impacts of introduced species cascade through trophic levels, they can cause indirect and counter-intuitive effects. To investigate the impact of invasive species at the network scale, we use a generalized food web model, capable of propagating changes through networks with a series of ecologically realistic criteria. Using data from a small British offshore island, we quantify the impacts of four virtual invasive species (an insectivore, a herbivore, a carnivore and an omnivore whose diet is based on a rat) and explore which clusters of species react in similar ways. We find that the predictions for the impacts of invasive species are ecologically plausible, even in large networks. Species in the same taxonomic group are similarly impacted by a virtual invasive species. However, interesting differences within a given taxonomic group can occur. The results suggest that some native species may be at risk from a wider range of invasives than previously believed. The implications of these results for ecologists and land managers are discussed.

1 Methods for collecting the data from the field and the literature The ecological network on Flat Holm, a 35ha island in the Bristol Channel, UK, was constructed from multiple sources (Tab. 1). A study of flower visitor data in 2008 divided the island into four habitat categories: grassland, scrub, buildings (any area with a man-made substrate) and coastal. One hundred minutes of surveying was carried out each day, the time divided between the four habitat types according to the proportion of island surface they covered. Each habitat type was surveyed by walking through it haphazardly looking for insects on flowers. Any insects found during the survey period were caught, preserved and later identified to species level by taxonomists.
Data on the species which are not pollinators or flowering plant come from the South East Wales Biological Records Centre. Interactions between species listed in the biological records are derived from Cramp [1], Pilorge [2], Harde [3], Snow & Snow [4], Luiselli [5] and Skinner [6] and vertebrateinvertebrate scavenger interactions were inferred from Atkinson (unpublished data). Atkinson defines scavenging species as species or genera collected from meat-baited pitfall traps set in mainland southwest England in summer 2012. Due to proximity, many of the same invertebrate species are found on Flat Holm.
Overall, this network contains 560 species (this is not necessarily an exhaustive list of all species on the island) described by 7 sub-networks: plant-Lepidoptera larvae, plant-other invertebrate herbivore, invertebrate-reptile, prey species-bird, plant-flower visitor, plant-seed feeding bird and vertebrateinvertebrate scavenger.
We can divide the links between species in the ecological network into three categories, according to the net gains and losses of each species in the interaction:

Parametrization of the Flat Holm food web
The generalized model is explored by randomly selecting sets of parameters subject to constraints imposed by ecology. The parameters ρ, σ, β, χ, µ depend upon the nature of the feeding relations between the different species, contained in the adjacency matrix. α depends on the inherent characteristics of the species. The remaining parameters φ, γ, ψ, λ are chosen from ecologically realistic ranges as contained in Ref. [7] (Tab. 1 in the main paper).
• α i is the species turnover rate. It is a time-scale parameter which is the same for species belonging in the same taxonomic group. Assuming this parameter is the reciprocal of the life spans of the species, we approximated them as 1 year for the plants (most of them are annual species), 10 years for the birds and reptiles, 0.5 year for the invertebrates and 0.1 year for the fungus.
• ρ i is the fraction of growth gained by predation. ρ i = 1 for predators and ρ i = 0 for primary producers. Each species gains biomass either by primary production or by predation, but not both.
• σ i is the fraction of biomass loss resulting from predation. The predation pressure for the top-predators (species without predator) is null, then σ i = 0 for them. As we have no more information on the predation pressure for the other species, we assume that 0.5 < σ i < 1.
Values are drawn from a uniform distribution.
• β i,j is the relative loss of species j to its predator i among all of j's predators. We assume that the prey j is equally eaten by its predators, so if i eats j then β i,j = 1/ (number of j's predators) and β i,j = 0 otherwise.
• χ i,j is the relative contribution of the prey j to predator i's diet. We assume that the predator i consumes equally all of its prey, so if i eats j then χ i,j = 1/ (number of i's prey) and χ i,j = 0 otherwise.
• µ i is the elasticity corresponding to the loss of biomass of species i through natural mortality (disease, natural ageing, etc.). In order to stabilize the system, the top-predators have a quadratic mortality rate, i.e. µ i = 2 which means that there is a quadratic relation between their mortality rate and their biomass X i . The mortality rate of the other species though increase linearly with their biomass: µ i = 1.
• φ i is the elasticity corresponding to the gain by primary production. Because of constraints such as space or nutrient availability, the primary production rate is assumed to be sub-linear with respect to the biomass species i. We then chose φ i uniformly between 0 and 1.
• ψ i is the elasticity corresponding to the gain by predation-ship. Individuals of the same species can interact (either by competition or by mutualism) and this introduces non-linearities in the predation rate. ψ i is taken uniformly between 0.5 and 1.5; ψ i = 1 would mean that the predation rate is linear with species i abundance, i.e. there is no direct interference between individuals of the same species.
• γ i is the elasticity corresponding to the gain by predation, but contrary to ψ i , according to the abundance of i's prey (not i). The availability of prey limits predation, γ i describes the feeding strategy of species i. As an example, a Holling type II response is close to a linear relation between the predation rate of i and i's prey abundance (γ i = 1). We allow γ i to be in the range of [0.5, 1.5].
• Λ i,j is the elasticity of prey switching. It describes the non-linearity of the contribution of the prey j to the diet of the predator population i. If the predator population does not chose its prey, Λ i,j = 1. But if the predator adapted a strategy to a specific prey population, Λ i,j = 2. As we do not know here the strategies of the predator populations, Λ i,j = 1 or 2 is picked with equal probability (0.5). We note that there is a positive correlation between the proportion of 2s and the stability of the system. For the full system, we simulated 10 6 Jacobian matrices, 9352 of them lead to a stable steady state, that is ≈ 1%. Fig. 1 shows the distribution of the random parameters of the kept simulations. For the simple system with 5 species (described below), 50000 Jacobian matrices were simulated, with ≈ 90% (44833) of them leading to a stable steady state. Figure 1: Distribution of the random parameters through the simulations which led to a stable steady state (in black). The red lines are the theoretical distributions of each of the random parameters according to the rules described above. The y-axis is the proportion of kept parameters among the drawn ones.

Derivation of the sensitivity and influence of the native species
As a result of the specific structure of the network, different species do not have the same importance within the food web. The Jacobian matrix holds the information relative to each species on how a perturbation of the system at the stable steady state can affect them and how it can propagate trough the network. The strength of the reaction of each species is held by the eigenvalues λ k . The (right) eigenvectors v k , that satisfy Jv k = λ k v k , show the relative amount of each affected species k. Likewise, the left eigenvectors w k (or the right eigenvectors of J T ), that satisfy w k J = λ k w k , show the direction in which the perturbed species k affects the other species through the links of the network. NB: the right and left eigenvalues are the same.
Aufderheide et al. [7] proposed the following definitions for the sensitivity and influence of species in a food web. The sensitivity Se of a species i represents how it is perturbed by disturbances propagating through the network.
with v k the right eigenvector associated to eigenvalue λ k . The sum is weighted by the invert of the eigenvalues which are relaxation time; so a low λ k corresponds to a high disturbance in reaction of the perturbation of the k th mode.
The influence Inf of a species i represents how it affects the other species if perturbed itself. It is computed as with w k the left eigenvector associated to eigenvalue λ k .
Both sensitivity and influence correspond to intrinsic qualities of the considered populations and their values can be computed directly from the Jacobian matrix. They present a general idea of how a perturbation -e.g. an newly arrived invasive species -can affect the abundances of the native species.

Confidence regions
Confidence regions are computed to give a guide to distinguish groups of species which are significantly different. The boundary an iso-contour of a Gaussian bi-dimensional distribution and the region contains 95% of all randomly drawn samples. We assume that the data are independent and normally distributed. The origin of the ellipse is placed at the mean of the observations for x and y directions independently. The alignment of the axis of the ellipse are the eigenvectors of the covariance matrix of the data. The lengths of the axis depend on its eigenvalues a 2 i = sλ i , where 2a i is the length of the major-axis (i = 1) or minor-axis (i = 2), λ 1 and λ 2 represent the sorted eigenvalues of the covariance matrix and the constant of proportionality s is defined as the scale of the ellipse. Using the Chi-square likelihood (P(s < 5.991) = 0.95), the 95% confidence ellipse has the semi-major axis length equal to √ 5.991λ 1 and the semi-minor axis length to √ 5.991λ 2 (Ref. [8]).

Results for the full system
We note that: • The bird-to-bird links are removed to avoid unstable states in the system. Indeed, if the species A eats the species B and B eats A as well, there is a perpetuum mobile and thus energy cannot be produced.
• To achieve stability, no species can have the same prey and predator links, though this is a widely used approximation (niche model [9]). Thus, species with exactly the same prey and predators are assembled in aggregated trophic groups.
The data on the full system can be found on-line [10] in the form of a label adjacency matrix (see definition below). Species in the same taxonomic group have similar sensitivity and influence (Fig. 2 Left). The structure of the network, rather than specific values of individual parameters, determines its dynamical behaviour. The influence and sensitivity of the species are not correlated (R 2 = 0.0196, F = 5.518 on 1 and 225 d.f., p = 0.02) and there is a greater variation in influence than sensitivity between species. The birds and reptiles are most influential, with reptiles also showing a relatively high sensitivity to perturbation. Quantities for the plant species show a large variation.
Meanwhile, we can see two different behaviours for the invertebrate group (Fig. 2 Right). A deeper analysis suggests that the invertebrates without prey in the network are amongst the least sensitive and influential ones ( figure 2 left). Indeed, the model considers invertebrates without prey as primary producers. As we have only considered simple trophic interactions in the adjacency matrix (and not other interactions such as pollination), some structurally important links are missing from the network.

Results for the simple system
The species of the Flat Holm ecosystem can be classified into five taxonomic groups: plant (P ), invertebrate (I), bird (B), fungus (F ) or reptile (R).
The label adjacency matrix A which gathers the information about the trophic interaction between the species is if there is no interaction between the species i and j. The bird is a top predator, the plant and the fungus are primary producers, the invertebrate feeds on the plant and the reptile feeds on the invertebrate.
Following the rules described in table 1 in the main text and this adjacency matrix, the non-random parameters are As an example, the plant species has a life-span of one year (α P = 1). It is considered as a primary producer (ρ P = 0). It is eaten as much by the bird as by the invertebrate (β ·,P = 1/2) and its mortality rate is a linear function of its abundance (µ P = 1). As another example, the life span of the bird is 10 years (α B = 1/10). This species is a predator (ρ B = 1) which eats equally its four prey: the plant, the invertebrate, the fungus and the reptile (χ B,· = 1/4). It is a top-predator (σ B = 0), then its mortality rate is a quadratic function of its abundance (µ B = 2). We find that the distribution of the influence and sensitivity is multi-modal: the modes correspond to each species of the network (Fig. 3). The multi-modality reflects the structure of the network, whereas the spread of the simulations reflects the variety of systems permitted by the generalized approach. . The regions are expected to be wider than for the full system, because of the lower variation between the means of the species than between the outcome of the simulations.

Bird groups: generalist and specialist species
The birds with the highest sensitivity and lowest influence are positively affected ones by the cyberherbivore (Fig. 4 Top and Middle). These birds have a lower number of trophic links than the birds with high influence and lower sensitivity (Fig. 5). The cyber-carnivore and cyber-insectivore affect all the birds in similar ways (Fig. 2 in the main paper). In Fig. 4 Bottom, the group of birds which is positively impacted by the cyber-omnivore is not the same one which is positively impacted by the cyber-herbivore. These positively affected birds are those which are not part of the diet of the cyber-omnivore.
The influence of each native species of the Flat Holm food web (apart from the invertebrates) is strongly correlated to its number of trophic links in the Flat Holm network in a non-linear way (Fig.  5). This suggest that generalist species may have larger and more widespread effects on the food webs.  and their number of trophic links with the rest of the Flat Holm network. The parameters of the non-linear model (in red) are adjusted with a Gauss-Newton algorithm (p < 10 −16 for both parameters) -the invertebrates are not taken in account in the fit, but are shown on the diagram.