More stories

  • in

    Rotating and stacking genes can improve crop resistance durability while potentially selecting highly virulent pathogen strains

    Overview of the model
    This model simulates the population and evolutionary dynamics of different pathogen strains, as they interact with different crop resistant cultivars planted in a single field over successive years. We assume one cultivar is planted each year and we consider a field divided into a finite number m of spatial units (representing limited spaces for infections, or potential lesion sites), in which the spatial aspect is implied rather than explicitly represented. For each year during the cropping season, a number of pathogen spores are released from the infested crop residues, it then lands on the crop plants leading to infections (Fig. 1). These infections are apportioned between the different pathogen strains depending on their previous abundance and interactions with the crop cultivar. At the end of the year, during the non-cropping season, these strains are assumed to sexually recombine in the crop residue. The number of spores released and the number of infections are considered as random variables. We denote these both quantities with an uppercase letter (for example N) in general sense, while their particular realisation or draw in the simulation will be noted with a lowercase letter (for example n). The model was developed using the R Language and Environment for Statistical Computing49.
    Figure 3

    Case 1, model predictions of total infection by each pathogen genotype (proportion of total locations infected, left), and the corresponding frequencies of each virulent allele (right) changing over time under different rotation strategies (From top to bottom: (S1) no rotation; (S2) rotation every year; (S3) rotation every 5 years; and (S4) rotation every year with stacked resistance genes). The parameters are at baseline values: the initial frequency of each virulent allele equals (5%), the fitness modifier is set at 0.9, the modifier of increase rate for non-virulent strains equals 0.05, and the initial amount of inoculum represents (10 %) of available locations.

    Full size image

    Figure 4

    Case 2, model predictions of total infection by each pathogen genotype (proportion of total locations infected, left), and the corresponding frequencies of each virulent allele (right) changing over time under different rotation strategies (from top to bottom: (S1) no rotation; (S2) rotation every year; (S3) rotation every 5 years; and (S4) rotation every year with stacked resistance genes). The initial frequency of each virulent allele equals (50%), and other parameters are at baseline values: the fitness modifier is set at 0.9, the modifier of increase rate for non-virulent strains equals 0.05, and the initial amount of inoculum represents (10 %) of available locations.

    Full size image

    We applied the general model described above to a specific situation with four genes of interaction where there are four different resistance genes that may or may not be deployed within each crop cultivar, and four virulence genes that may or may not be present within each pathogen strain. We assume that the presence of each virulence gene reduces the fitness of a strain independently. Specifically, for each strain i, we define the fitness of that strain (delta _i = delta ^{n_{vir,i}}), where (n_{vir,i}) is the number of virulence alleles present in strain i, and (delta) is a fixed model parameter with potential values between zero and one (Table 1).
    Table 1 Baseline model parameter values used for our analysis, with alternative values shown in parentheses.
    Full size table

    We first set the model parameter values to define a baseline situation where there is a relatively small fitness penalty for virulence alleles (i.e. (delta) is very close to 1, where the value 1 means no penalty); the pathogen has a relatively low ability to reproduce if it does not carry effective virulence genes (i.e. low value for (epsilon), in this baseline situation equal to 0.05); the initial virulence allele frequency (Init.freq) is relatively low, reflecting a low historical selection pressure and lastly the initial quantity of pathogen (Init.path) is also low at (10 %) of carrying capacity (Table 1). We then considered and compared four different strategies for rotating resistant crop cultivars:
    S1.
    No rotation, the same cultivar with only one gene of resistance is employed every year;

    S2.
    A cultivar with a single gene of resistance is employed each year, and the gene of resistance in the cultivar is changed every year, giving a 4-year rotation;

    S3.
    A cultivar with a single gene of resistance is employed each year, and the gene of resistance in the cultivar is changed every 5 years, giving a 20-year rotation; and

    S4.
    A cultivar with two genes of resistance (i.e. pyramided resistance) is employed each year, and the genes of resistance in the cultivar are changed every year, giving a 2-year rotation

    We then investigated how different parameterisations of the model would interact with the selected rotation strategies. We develop four cases in addition to the baseline case described above:
    Case 1.
    Baseline scenario (Table 1).

    Case 2.
    Baseline scenario, except for Init.freq which was increased from 0.05 to 0.5.

    Case 3.
    Baseline scenario, except for Init.freq which was increased from 0.05 to 0.5 and (delta) which was decreased from 0.9 to 0.7.

    Case 4.
    Baseline scenario, except for Init.freq which was increased from 0.05 to 0.5, (delta) which was decreased from 0.9 to 0.7 and (epsilon) which was increased from 0.05 to 0.5.

    Case 5.
    Baseline scenario, except for (epsilon) which was increased from 0.05 to 0.5.

    Genetics
    Cultivar and pathogen strain are both defined through their genotype being restricted to a specific set of interaction genes (loci) related to resistance (for the cultivar) or virulence (for the pathogen). Each gene has two versions (alleles): virulence or avirulence allele for the pathogen and resistance or susceptibility for the cultivar. Virulence and resistance are represented with a 1 and avirulence and susceptibility are represented with a 0 (Fig. 2). If we call ({mathcal {I}}) the set of strains and if (nu) genes of interaction are involved, then the total number of strains will be (left| {mathcal {I}}right| = 2^{nu }). During the infection process, after pathogen spores land on the cultivar, an interaction factor (beta (i,c)) defines the relative rate at which strain i can reproduce within a field of cultivar c, for each strain and cultivar combination (Fig. 2). We consider that a strain overcomes the cultivar genotype when the strain has a virulence allele for every resistance allele of the cultivar (Fig. 2), in which case (beta (i,c)=1), indicating maximum reproduction rate. Otherwise, if the strain does not have a virulence allele for every resistance allele of the cultivar, (beta (i,c)= epsilon), where (epsilon) is a model parameter with constant value (0 le epsilon < 1), indicating a less-than-maximum reproduction rate. As such, (epsilon) is the model parameter modifying the growth and reproduction of pathogen strains not carrying multiple virulence alleles (e.g. 0100) and/or an avirulent pathogen strain (e.g. 0000) (Fig. 2). Accordingly, lower (closer to 0) (epsilon) values represent reduced ability to grow and reproduce in pathogen strains with increasing number of avirulence alleles. Moreover, any strain i with one or more virulence genes is also assumed to suffer a fitness penalty (delta _i) depending on the number of genes involved. Together these interaction factors make a cultivar-strain interaction matrix (B = (beta (i,c))). This code and method for modelling resistance and virulence interactions (without fitness penalty) is similar to those in previous studies48,50. Figure 5 Case 3, model predictions of total infection by each pathogen genotype (proportion of total locations infected, left), and the corresponding frequencies of each virulent allele (right) changing over time under different rotation strategies (from top to bottom: (S1) no rotation; (S2) rotation every year; (S3) rotation every 5 years; and (S4) rotation every year with stacked resistance genes). The fitness modifier is set at 0.7, the initial frequency of each virulent allele equals (50%), and other parameters are at baseline values: the modifier of increase rate for non-virulent strains equals 0.05, and the initial amount of inoculum represents (10 %) of available locations. Full size image Figure 6 Case 4, model predictions of total infection by each pathogen genotype (proportion of total locations infected, left), and the corresponding frequencies of each virulent allele (right) changing over time under different rotation strategies (from top to bottom: (S1) no rotation; (S2) rotation every year; (S3) rotation every 5 years; and (S4) rotation every year with stacked resistance genes). The fitness modifier is set at 0.7, the initial frequency of each virulent allele equals (50%), the modifier of increase rate for non-virulent strains equals 0.5, and other parameters are at baseline values: the initial amount of inoculum represents (10 %) of available locations. Full size image Initial genetic structure of pathogen population At the start of each case, we define the initial proportion of each pathogen genotype using the equation: $$begin{aligned} strains.init = Init.freq^{nr}left( 1-Init.freqright) ^{4-nr} end{aligned}$$ (1) where strains.init is the initial proportion of each pathogen genotype; Init.freq is the frequency of the virulent genes as set by each case and nr is the number of virulent genes present in a given pathogen genotype. We then used a random Poisson distribution generator (rpois function from the stats package in R) to obtain the initial number of spores for each pathogen genotype, where the mean of the Poisson distribution is the proportion of a given pathogen genotype multiplied by the pre-determined pathogen load (Init.path, Table 1). Model dynamics The annual dynamics (Fig. 1) can be divided into two main phases: the phase of parasitic activity, representing events occurring through the cropping season, and the phase of dormancy, representing events occurring between the cropping seasons. During the phase of parasitic activity, the pathogen produces spores which are spread both through the air (sexual ascospores) and via water splash (asexual conidia). These spores may then infect leaves and stems of the cultivar, resulting in new lesions of different strains. During the phase of dormancy, the pathogen remains within the infected crop residue and sexual recombination occurs. These processes are modelled with four steps, three for the parasitic phase and one for the dormancy phase. Figure 7 Case 5, model predictions of total infection by each pathogen genotype (proportion of total locations infected, left), and the corresponding frequencies of each virulent allele (right) changing over time under different rotation strategies (from top to bottom: (S1) no rotation; (S2) rotation every year; (S3) rotation every 5 years; and (S4) rotation every year with stacked resistance genes). The modifier of increase rate for non-virulent strains equals 0.5 and other parameters are at baseline values: the initial frequency of each virulent allele equals (5%), the fitness modifier is set at 0.9, and the initial amount of inoculum represents (10 %) of available locations. Full size image Total amount of spores released First, the model generates the amount of pathogen spores of each strain that is released, using the equation: $$begin{aligned} lambda _{{ released},i}(t) = alpha . n_{{ recombined},i}(t-1) end{aligned}$$ (2) where (lambda _{{ released},i}(t)) represents the expected dispersed propagule (spore) pressure in the field due to strain (i in {mathcal {I}}) during the year t, the parameter (alpha) represents the normal rate of growth for the pathogen from 1 year to the next, and (n_{{ recombined},i}(t-1)) represents the number of spatial units or locations infected by the strain i at the end of the previous year and after genetic recombination. The actual quantity of pathogen strain i released in the current year, (n_{{ released},i}(t)) is then simulated as a Poisson random variable: $$begin{aligned} N_{{ released},i}(t) hookrightarrow {mathcal {P}}(alpha . n_{{ recombined},i}(t-1)) end{aligned}$$ (3) The infective pressure (lambda _{{ infected},i}(t)) is then calculated as: $$begin{aligned} lambda _{{ infected},i}(t) = beta (i,c(t)) . delta _i . n_{{ released},i}(t) end{aligned}$$ (4) where (beta (i,c(t))) is the interaction factor between the strain i and the cultivar c(t) i.e. the cultivar grown in year t, and (delta _i) is the fitness penalty for the particular strain i. Total number of infections Second, the model calculates the total number of infected sites, following a binomial distribution: $$begin{aligned} N_{ infected}(t) hookrightarrow {mathcal {B}}left( m, 1 - prod _{i = 1}^{2^{nu }} (1- rho _i(t))^{n_{{ released},i}(t)}right) end{aligned}$$ (5) where (rho _i(t)) is the probability that a particular given location (among the m possible locations in the field) during year t, will have a given spore from strain i fall down on it and cause a lesion, and thus (displaystyle 1 - prod nolimits _{i = 1}^{2^{nu }} (1- rho _i(t))^{n_{{ released},i}(t)}) represents the probability that at least one of the (displaystyle n_{ released}(t) = sum nolimits _{i = 1}^{2^{nu }} n_{{ released},i}(t)) spores produces a lesion. This equation can be justified in more detail as follows: $$begin{aligned}&P({At; least; one; of; the; n_{ released}(t); spores; produces; a; lesion}) \&quad = 1 - P({ No; released; spores; produces; a; lesion}) \&quad = 1 - prod _{i = 1}^{2^{nu }} P({ A; single; released; spore; of ; strain ; i ; doesn't; produce; a; lesion})^{n_{{ released},i}(t)} \&quad = 1 - prod _{i = 1}^{2^{nu }} (1 - P({A; single; released; spore; of ; strain ; i ; produces; a; lesion}))^{n_{{ released},i}(t)} end{aligned}$$ We assume that a spore will fall on any of the m specific locations with the same probability independently of its infection capabilities. The number of locations m is assumed to be the same for all years whatever the cultivar grown and thus, this probability is independent of the time dimension. Next, we assume that the probability that a spore will induce an infection depends on the interaction factor between the crop cultivar genotype and the pathogen strain (beta (i,c(t))) together with the fitness penalty for that strain (delta _i) . These assumptions mean that: $$begin{aligned} rho _i(t) = P({a; spore; fall; down; on; a; given; location; during ; year; t; and; causes; a; lesion}) \ = P({a; spore; falls; on; a; given; place; where; c(t); is; grown })times \ P({ the; spore; causes; a; lesion} mid { the; spore; falls; on; a; place; where; c(t); is; grown}) \ = frac{1}{m} . beta (i,c(t)) . delta _i end{aligned}$$ Number of infections for each strain Third, the number of infections of each strain is derived from the total number of infections depending on genetic interactions between each strain and crop cultivar being employed that year. Specifically, the total number of infections (N_{ infected}(t) = n_{ infected}(t)) is distributed among the different strains using the multinomial distribution: $$begin{aligned} left( N_{{ infected},1}(t), ldots , N_{{ infected},2^{nu }}(t)right) hookrightarrow {mathcal {M}}left( frac{lambda _{{ infected},1}(t)}{lambda _{ infected}(t)},ldots , frac{lambda _{{ infected},2^{nu }}(t)}{lambda _{{ infected}}(t)}, n_{ infected}(t)right) end{aligned}$$ (6) where (displaystyle lambda _{{ infected}} (t) = sum nolimits _{i = 1}^{2^{nu }} lambda _{{ infected},i}(t)). The number of infected sites due to strain i, without no loss of generalities, follows then the binomial distribution (displaystyle {mathcal {B}}left( n_{ infected}(t),frac{lambda _{{ infected},i}(t)}{lambda _{{ infected}}(t)}right)). Genetic recombination The fourth step involves simulating the process of sexual recombination, where new quantities of each strain are generated based on the previous quantities of each strain. At the end of the year t, we calculate the frequencies (f_j(t)) of each virulent version of each gene from the different genotypes of strains in the crop stubble. We let the genotype of any new spore be represented by a random vector (displaystyle G_i(t) = left( G_{i,1}(t), ldots ,G_{i,nu }(t)right)), where each (G_{i,j}(t)) is a Bernoulli random variable (displaystyle {mathcal {B}}(1,f_j(t))). This vector representation of genotype follows the coding illustrated in (Fig. 2). Assuming that strains recombine independently gene by gene, the probability that (G_i(t)) will be a particular genotype (displaystyle g_i(t) = left( g_{i,1}(t), ldots ,g_{i,nu }(t)right)) is given by: $$begin{aligned} p_i(t) = Pleft( G_i(t) = g_i(t)right) &= prod _{j=1}^{nu } P(G_{i,j}(t) = g_{i,j}(t)) nonumber \ &= prod _{j=1}^{nu } f_j(t)^{g_{i,j}(t)} left( 1 - f_j(t)right) ^{left( 1 - g_{i,j}(t)right) } end{aligned}$$ We can also confirm that across all possible genotypes these probabilities sum to one: $$begin{aligned} sum _{i=1}^{2^{nu }} Pleft( G_i(t) = g_i(t)right) = sum _{i=1}^{2^{nu }} prod _{j=1}^{nu } f_j(t)^{g_{i,j}(t)} left( 1 - f_j(t)right) ^{left( 1 - g_{i,j}(t)right) } = 1 end{aligned}$$ (7) If we shorten the notation for (P(G_i(t) = g_i(t))) to be (p_i(t)) then we can define the recombined version of infected numbers of units of each strain with the following multinomial distribution: $$begin{aligned} left( N_{{ recombined},1}(t), ldots , N_{{ recombined},2^{nu }}(t)right) hookrightarrow {mathcal {M}}left( p_1(t),ldots , p_{2^{nu }}(t), n_{ infected}(t)right) end{aligned}$$ (8) Poisson, binomial and multinomial distribution In plant pathology, it is often relevant to model infections by a random variable. Let’s imagine a released spore flying in the air, we can say that this spore will land on a specific leaf and infect it with a given probability p, then it won’t with probability (1 - p) because these are the only two possible events. We can define Y a random variable to model the situation. If we say the event ({Y = 1}) represents the success of the event (landing and infection) and ({Y = 0}) represents the failure, with this definition we say that Y follows a Bernoulli distribution. The values attributed to the variable depending on the events allow the following generalisation: If we consider n spores, each of them realizing an infection on a specific plant area they fell on with the same probability p, then we can associate to each spore a Bernoulli distribution (Y_i) where (i in {1,ldots ,n }). If we are interested in the total number of infections occurring on this leaf, assuming the fact that they will happen independently of each other, we can model this situation by the variable (displaystyle S = sum nolimits _{i = 1}^n{Y_i}), called binomial variable. We can also denote briefly (S hookrightarrow {mathcal {B}}(n,p)), where n represents the number of events and p the probability of success of each event. Moreover, the Bernoulli variable Y is related to binomial distribution in the way that we can write (Y hookrightarrow {mathcal {B}}(1,p))51,52. Usually it is more likely to model such events by a Poisson law rather than binomial law53,54. When the number of events is so big that we can approximate it by infinity, and when the probability of success of each event is very small, close to zero, it is possible to link both Poisson and binomial distribution through their respective expectations. So if (lim nolimits _{begin{array}{c} n nearrow +infty \ p searrow 0 end{array}} n*p = lambda ,) then if we define (X hookrightarrow {mathcal {P}}(lambda )), we have (S xrightarrow {text {distribution}} X). Returning to our example, that means that if we have a ‘close to infinity’ number of spores that could fall onto a given plant and infect it with a very small probability p for each of them and still acting independently, we can model the total number of infections by both S or X. Even if there are millions and millions of spores released, this amount is still small compared to infinity, so using X is still a modelling approximation. The use of binomial or Poisson laws depends on the complexity of the situation. For example, if the modeller wants to simulate a model where he anticipates 15 infections, they can use (X hookrightarrow {mathcal {P}}(15)) or (S hookrightarrow {mathcal {B}}(10000,0.0015)) or (S hookrightarrow {mathcal {B}}(1000000,0.000015)). We consider now a situation where the plant is attacked by a big number of spores, but with different genotypes modifying their ability to infect, some strains being more efficient than others. To model this situation, we can use a vector of variables, each component representing the number of successes due to a specific genotype. We can choose a vector of binomial number or Poisson number. If we consider the case of a threshold in terms of available space to be infected (a maximum number of infections for the plant), such that spores of different strains are competing for those places, we suggest using a vector of random numbers that follows a multinomial law. This distribution derives from the binomial law, although each component is a specified binomial distribution defined from the parameters of the multinomial distribution. But, it is still possible to interpret some of these components via a conditional Poisson distribution. From binomial to multinomial distribution The binomial distribution is a particular case of the multinomial distribution. We consider S a binomial distribution of parameters (n, p) counting the number of success of n independent events where the basic probability of success is p. Let U the random variable be defined by (n-S) the number of failures. In the case where S represents the number of infections, U represents the number of uninfected places. The probability to get k infections is given by: $$begin{aligned} P(S = k) = {n atopwithdelims ()k} p^k (1-p)^{n-k} = {n atopwithdelims ()n-k} p^k (1-p)^{n-k} = P(U = n-k) end{aligned}$$ (9) As a result, the probability of having k success is the same that having (n-k) failures. Then the Eq. (9) shows that U follows a binomial distribution with parameters ((n, 1-p)). We can also say that the couple (S, U) follows a multinomial distribution of parameter ((p, 1-p, n)), that we can denote ((S,U) hookrightarrow {mathcal {M}}(p, 1-p, n)). In a more general way, the analogue of the binomial distribution is the multinomial distribution, where each trial results in exactly one of some fixed finite number k possible successes, with probabilities (p_1), ..., (p_k) (so that (p_ige 0) for i = 1, ..., k and (sum nolimits _{i=1}^k p_i = 1)), and there are n independent trials. Then if the random variables (X_i) indicate the number of times outcome number i is observed over the n trials, the vector (X = (X_1, ldots , X_k)) follows a multinomial distribution with parameters n and p, where (p = (p_1, ldots , p_k)), that we can also write ({mathcal {M}}left( p_1,ldots ,p_n, N = kright))55. From Poisson to multinomial distribution We consider here a total number of successes (meaning in our example a number of spores that fall on a place and infect it) X being the sum of the infections due to (omega) different strains (X_i) ((1le i le omega)). If we consider that each (X_i) follows a Poisson law of parameter (lambda _i) and that they are all independent, then X follows a Poisson law of parameter (displaystyle lambda = sum nolimits _{i=1}^{omega } lambda _i). The distribution of each (X_i) conditionally to (X = n) is a binomial distribution ({mathcal {B}}(n,frac{lambda _i}{lambda })). We can prove it for all variable (X_j), with (j in {1,ldots ,omega }): $$begin{aligned} Pleft( X_j = k left| right. sum _{i=1}^{omega } X_i = n right) &= frac{P left( X_j = k,displaystyle sum nolimits _{begin{array}{c} i = 1 \ ine j end{array}}^{omega } X_i = n-k right) }{P left( displaystyle sum nolimits _{i = 1}^{omega } X_i = n right) } \ &= frac{Pleft( X_j = k right) Pleft( displaystyle sum nolimits _{begin{array}{c} i = 1 \ i ne j end{array}}^{omega } X_i = n-k right) }{P left( displaystyle sum nolimits _{i = 1}^{omega } X_i = n right) } end{aligned}$$ that we obtain using the Bayes formula for conditioning and the use of independence between the (X_i)’s. Then we replace the probabilities by their Poisson values: $$begin{aligned} Pleft( X_j = k left| right. sum _{i = 1}^{omega } X_i = n right) &= frac{e^{-lambda _j}{lambda _j}^k}{k!} frac{e^{- displaystyle sum nolimits _{begin{array}{c} i = 1 \ ine j end{array}}^{omega } lambda _i}{left( displaystyle sum nolimits _{begin{array}{c} i = 1 \ ine j end{array}}^{omega } lambda _i right) }^{n-k}}{(n-k)!} frac{n!}{e^{-displaystyle sum nolimits _{i = 1}^{omega } lambda _i}{left( displaystyle sum nolimits _{i = 1}^n lambda _i right) }^{omega }} \ &= {n atopwithdelims ()k} frac{{lambda _j}^k {left( displaystyle sum nolimits _{begin{array}{c} i = 1 \ ine j end{array}}^{omega } lambda _i right) }^{n-k}}{{left( displaystyle sum nolimits _{i = 1}^{omega } lambda _i right) }^n} = {n atopwithdelims ()k}{left( frac{lambda _j}{displaystyle sum nolimits _{i = 1}^{omega } lambda _i}right) }^k {left( frac{displaystyle sum nolimits _{begin{array}{c} i = 1 \ ine j end{array}}^{omega } lambda _i}{displaystyle sum nolimits _{i = 1}^{omega } lambda _i}right) }^{n-k} end{aligned}$$ Generalizing this result to the random vector of the (displaystyle (X_i)_{1 le ile omega }) for (omega) strains, the distribution of this vector conditionally to the total number X is a multinomial distribution ({mathcal {M}}left( frac{lambda _1}{lambda },ldots ,frac{lambda _n}{lambda }, X = nright))55. Properties of the model Let (X_1),..., (X_{2^{nu }}) independent random variables such that (X_j hookrightarrow {mathcal {P}}(lambda _{{ infected},j}(t))) for all (j in {1,ldots ,2^{nu }}), we have the following results: A. When (m rightarrow infty), (displaystyle N_{ infected}(t) hookrightarrow {mathcal {P}}(sum nolimits _{j=1}^{2^{nu }} lambda _{{ infected},j}(t))), B. For all (j in {1,ldots ,2^{nu }}), (displaystyle N_{ infected, j}(t) xrightarrow {text {distribution}} X_j left| right. sum nolimits _{i = 1}^{2^{nu }} X_i = n), C. With A and B when (m rightarrow infty), it is equivalent to either simulate (N_{ infected}(t)) then the conditional multinomial vector (displaystyle left( N_{{ infected},1}(t), ldots , N_{{ infected},2^{nu }}(t)right)) conditionally to the realisation (n_{ infected}(t)), or to simulate directly the previously defined variables (X_1),..., (X_{2^{nu }}). The number of infected sites due to strain j, without any loss of generalities, follows the binomial distribution (displaystyle {mathcal {B}}left( n_{ infected}(t),frac{lambda _{{ infected},j}(t)}{lambda _{{ infected}}(t)}right)). It is important to notice that it is the same law as a Poisson variable with parameter (lambda _{{ infected},j}(t)) conditionally to the realisation (n_{ infected}(t)) of a Poisson variable with parameter (lambda _{{ infected}}(t)). Referring to formula (10), we can see that when the number of sites available for infection goes towards infinity, meaning that (N_{{ infected}}(t)) behaves like a Poisson law of parameter (sum nolimits _{i = 1}^{2^{nu }} lambda _{{ infected},i}(t)), then the variables (displaystyle left( N_{{ infected},i}(t)right) _{1 le i le 2^{nu }}) behave like independent Poisson law of respective rates (displaystyle left( lambda _{{ infected},i}(t)right) _{1 le i le 2^{nu }}). Proof of the properties of the model A. With the help of the reminder, we just have to prove this result: $$begin{aligned} lim _{m rightarrow infty } Eleft( N_{ infected}(t)right) = sum _{i=1}^{2^{nu }} lambda _{{ infected},i}(t), end{aligned}$$ (10) which could be obtained with the mean value theorem56. It means that if the total number of places available for infections was unlimited, these infections could be regarded as being Poisson distributed, with infection pressure as defined previously. We consider the notation of (5), and to simplify the formula we will note: (rho _i = frac{1}{m} . beta _i) and because the result (10) does not depend on time we reduce the notation such that (10) is equivalent to: $$begin{aligned} lim _{m rightarrow infty } Eleft( N_{ infected}right) = sum _{i=1}^{2^{nu }} lambda _{{ infected}, i}, end{aligned}$$ (11) and then we want to prove that: $$begin{aligned} lim _{m rightarrow infty } m . left( 1 - prod _{i=1}^{2^{nu }} left( 1- frac{beta _i}{m}right) ^{n_{{ released},i}}right) = sum _{i=1}^{2^{nu }} lambda _{{ infected}, i} end{aligned}$$ (12) Replacing m by (frac{1}{x}), with (xne 0), the latest equation equals: $$begin{aligned} lim _{x rightarrow 0} frac{1}{x}. left( 1 - prod _{i=1}^{2^{nu }} (1- xbeta _i)^{n_{{ released},i}}right) = sum _{i=1}^{2^{nu }} lambda _{{ infected}, i} end{aligned}$$ (13) We define (displaystyle f_{beta , n_{released}}(x) = prod _{i=1}^{2^{nu }} f_{i,({beta , n_{released}})}(x) = prod _{i=1}^{2^{nu }} (1 - xbeta _i)^{n_{{ released},i}}). Taking into account the fact that $$begin{aligned} f_{beta , n_{released}}'(x) = left( prod _{i=1}^{2^{nu }} f_{i,({beta , n_{released}})}(x)right) ' = sum _{i=1}^{2^{nu }} left[ f_{i,({beta , n_{released}})}'(x) prod _{begin{array}{c} i=1 \ jne i end{array}}^{2^{nu }}f_{j,({beta , n_{released}})}(x)right] , end{aligned}$$ (14) we apply the mean value theorem (56) to the derivable function (f_{beta , n_{released}}), we got the following result: $$begin{aligned} lim _{x rightarrow 0} frac{left( 1 - f_{beta , n_{released}}(x)right) }{x} &= - lim _{x rightarrow 0} frac{left( f_{beta , n_{released}}(0) - f_{beta , n_{released}}(x)right) }{0 - x} nonumber \ &= -left( f_{beta , n_{released}}'(0)right) = sum _{i=1}^{2^{nu }} beta _i n_{released, i} end{aligned}$$ (15) that finishes the proof of point A. B. The result is immediate knowing the upper reminder concerning the Poisson–Multinomial laws relationship. We just have to take the value of (omega = 2^{nu }). C. When m is close to infinity, (N_{ infected}(t)) follows a Poisson distribution whose parameter (expectation) is a sum of parameters. A property of Poisson distribution is that the law of a sum equals in distribution the sum of independent Poisson variables with the respective terms. So that we can rewrite B: For all (j in {1,ldots ,2^{nu }}), (displaystyle N_{ infected, j}(t) xrightarrow {text {distribution}} X_j left| right. N_{ infected}(t) = n). More

  • in

    Researchers using environmental DNA must engage ethically with Indigenous communities

    1.
    Day, K. et al. Endanger. Species Res. 40, 171–182 (2019).
    Google Scholar 
    2.
    Fernandes, K. et al. Restor. Ecol. 26, 1098–1107 (2018).
    Google Scholar 

    3.
    van der Heyde, M. et al. Mol. Ecol. Resour. 20, 732–745 (2020).
    Google Scholar 

    4.
    Thomsen, P. F. & Willerslev, E. Biol. Conserv. 183, 4–18 (2015).
    Google Scholar 

    5.
    Giguet-Covex, C. et al. Nat. Commun. 5, 3211 (2014).
    PubMed  Google Scholar 

    6.
    Willerslev, E. et al. Science 300, 791–795 (2003).
    CAS  PubMed  Google Scholar 

    7.
    Willerslev, E. et al. Nature 506, 47–51 (2014).
    CAS  PubMed  Google Scholar 

    8.
    Slon, V. et al. Science 356, 605–608 (2017).
    CAS  PubMed  Google Scholar 

    9.
    Claw, K. G. et al. Nat. Commun. 9, 2957 (2018).
    PubMed  PubMed Central  Google Scholar 

    10.
    Garrison, N. A. et al. Annu. Rev. Genomics Hum. Genet. 20, 495–517 (2019).
    CAS  PubMed  Google Scholar 

    11.
    Kowal, E. in Biomapping Indigenous Peoples: Towards an Understanding of the Issues (eds Berthier-Folgar, S. et al.) 329–347 (Rodopi, 2012).

    12.
    Adams, K., Faulkhead, S., Standfield, R. & Atkinson, P. Women Birth 31, 81–88 (2018).
    PubMed  Google Scholar 

    13.
    National Health and Medical Research Council (NHMRC), Australian Research Council (ARC) & Australian Vice-Chancellors’ Committee (AVCC) National Statement on Ethical Conduct in Human Research (2007) – Updated 2015 (National Health and Medical Research Council, 2015).

    14.
    Pawu-Kurlpurlunu, W. J., Holmes, M. & Box, L. A. Ngurra-kurlu: A Way of Working with Warlpiri People DKCRC Report 41 (Desert Knowledge CRC, 2008); https://go.nature.com/3jl4TR4

    15.
    Rose, D. B. Dingo Makes Us Human: Life and Land in an Australian Aboriginal Culture (Cambridge Univ. Press, 1992).

    16.
    Stanner, W. E. H. in White Man Got No Dreaming: Essays 1938-1973 198–248 (Australian National Univ. Press, 1979).

    17.
    Lewis, C. M. Jr, Obregón-Tito, A., Tito, R. Y., Foster, M. W. & Spicer, P. G. Trends Microbiol. 20, 1–4 (2012).
    CAS  PubMed  Google Scholar 

    18.
    Ma, Y., Chen, H., Lan, C. & Ren, J. Protein Cell 9, 404–415 (2018).
    PubMed  PubMed Central  Google Scholar 

    19.
    Kistler, L., Ware, R., Smith, O., Collins, M. & Allaby, R. G. Nucleic Acids Res. 45, 6310–6320 (2017).
    CAS  PubMed  PubMed Central  Google Scholar 

    20.
    Andersen, K. et al. Mol. Ecol. 21, 1966–1979 (2012).
    CAS  PubMed  Google Scholar 

    21.
    Haile, J. et al. Mol. Biol. Evol. 24, 982–989 (2007).
    CAS  PubMed  Google Scholar  More

  • in

    A birdstrike risk assessment model and its application at Ordos Airport, China

    Study area and data survey method
    Ordos Civil Aviation Airport, our study site, is located in the southwestern part of the Inner Mongolia Autonomous Region, China. It is characterized by mid-temperate continental climate and is located on the northeastern edge of the Mu Us Desert. Its main characteristics are a long winter and short summer, but has four distinct seasons. The mean annual temperature is 6.2 ℃, and the mean annual precipitation is 358 mm, mostly concentrated between June and August. The mean annual wind speed is 3.6 m/s17. There are five main land cover types in Ordos civil aviation airport and its surrounding areas: farmland, residential areas, woodland, shrub grassland and wetland.
    We used the line transects method and point counts method to investigate the environment and birds within our study area. It is a commonly used bird survey method, and is also widely used in the survey of birds in and around airports. The line transects were 1000 m × 100 m and the walking speed was 1.5–2.0 km/h. We observed and recorded birds with 10 × 50BA and 30 × 77BA Leica telescopes, SLR digital cameras (Canon 5D Mark III) with telephoto lenses (Canon 100–400 mm). The point count method we used had an observation radius of 200 m. We observed birds with 10 × 50BA Leica binoculars and a 30 × 77BA Leica fixed-mount spotting scope18,19,20,21,22. The flight altitude was estimated using a visual comparison method: an altimeter was used to measure the height of trees and buildings in the observation area and these heights were used to estimate the flight altitude of observed birds. Bird identifications were based on A Field Guide to the Birds of China23.
    The survey areas were divided into three areas: A, B and C. Section A was located within the boundary of the airport. Section A surveys consisted of five shrub grassland transects, with one transect on the runway and one on the apron. Section B was the area within 4 km of the center of the airport (but excluding Section A). There were five woodlands, nine shrub grassland, two farmland and four residential areas. Section C consisted of all areas located with 8 km of the center of the airport, excluding Sections A and B. In Section C we established three woodlands transects, three shrub-grassland transects, two farmland transects and four wetland transects. The species, quantity, distribution, cluster and flight altitude of birds in 39 transects or point counts set up within 8 km of the airport and its surrounding areas were investigated monthly by the method of sample strip or sample point (wetland using sample point method). A total of 468 individual point count or transect surveys were conducted over the study year.
    Birdstrike risk assessment model
    Investigating the bird situation in the airport and surrounding areas is a prerequisite for birdstrike prevention. The establishment of a scientific and standardized risk assessment process for birdstrike prevention (Fig. 2) is helpful for the systematic evaluation of birdstrike risk. This model is based on the ISO 31000 risk management process24—risk identification, risk analysis, risk assessment, risk response, risk recording and reporting, communication and consultation, monitoring and review. A flow chart for bird strike risk assessment was constructed.
    Figure 2

    Flow chart of the airport birdstrike risk assessment process.

    Full size image

    The occurrence of a birdstrike is a matter of probability. The consequences of a birdstrike are a matter of severity, with loss of aircraft or life occurring in extreme cases. Together they combine to determine birdstrike risk, and thus our five risk factors are meant to capture severity and likelihood. The first risk factor is the comparative number, which is important for the simple reason that if a bird species collides with an airplane, a greater number of birds have more serious consequences for an airplane. Among the bird strike events between 2007 and 2014 with the largest record impact energy, half of them involved species in the family Anatidae, and they were all birds with a relatively large comparative number25,26. The second risk factor is bird weight. The greater the weight of a bird, the greater the force generated by an aircraft impact, and the severity of birdstrikes will also increase. Flight altitude is an important factor in the analysis of birdstrike risk12. According to ICAO data, we use 40 m as the critical value of the risk zone. If the average flight height of a bird species is closer to the critical value, the risk of birdstrike will be higher12. Our fourth risk factor is a clustering coefficient, which relates to the living habits of a bird species to move in large groups. If a bird species often gather in large numbers, then the possibility of encountering an aircraft and causing a birdstrike event is greater. This is due to the nature of the collective behavior of birds while flying in flocks of murmurations. Following large, tight formations, birds make fewer independent moving decisions, being forced to constantly react to the movements of their neighbors and having their view partially obscured. They may not have space to avoid oncoming aircraft, or may lack the freedom and alert to choose a successful escape path leading to a higher probability of collision with the aircraft27,28. About 80% of birdstrikes occur during the take-off, climb, approach, and landing phases of flights12,13,29,30, so the distance between bird activity from the flight zone is also an important factor in assessing the probability of birdstrikes. Combining with the above analyses, a risk assessment matrix based on the five factors of bird number, weight, flight altitude, cluster coefficient and range of activity was proposed to assess the risk level of bird species in the airport and its surrounding area within 8 km.
    Risk factor assignment
    1.
    Comparative number = (the number of individual birds/the number of individuals with the most number of birds) × 100.

    2.
    Comparative weight = (estimated weight of all birds of a single species/the largest weight of all birds of any species) × 100.

    3.
    Risk coefficient of flight height:

    Flight height H (m)
    Risk coefficient of flight height
    H  > 100
    0.1
    100 ≥ H  > 50
    0.5
    50 ≥ H  > 30
    1
    30 ≥ H  > 5
    0.5
    5 ≥ H
    0.1

    4.
    Clustering coefficient assignment:

    Number of individuals of a cluster
    Cluster coefficient
    N  > 100
    1
    100  > N ≥ 20
    0.5
    20  > N ≥ 3
    0.2
    3  > N ≥ 1
    0

    5.
    Activity range risk coefficient assignment: according to the bird species observed area, it could be divided into three levels: activities in flight area, activities within 4 km from flight area, activities within 8 km from flight area but not within 4 km. If a bird species has activity in each area, the nearest one to the flight zone will be used as the input for the risk assessment model. The birds distributed in these three regions were assigned 0.9, 0.6 and 0.3 respectively.

    Risk assessment matrix

    $$ {text{Likelihood }} = , left( {{text{cluster coefficient }} + {text{ Risk coefficient of flight height }} + {text{ Activity range risk coefficient}}} right) , times { 1}00 , /{ 3} $$

    $$ {text{Severity }} = , left( {{text{comparative number }} + {text{ comparative weight}}} right) , times { 1}00/{2} $$

    The expert evaluation method is used to determine the numerical range31 (Table 1).
    Table 1 Birdstrike likelihood and severity rating.
    Full size table

    According to the very low, low, moderate, high and very high levels of possibility and severity (Table 1), the level of potential threatening birds are divided into three risk levels: high danger (level 3), moderate danger (level 2), and low danger (level 1). (Table 2).
    Table 2 Airport birdstrike risk assessment matrix.
    Full size table

    Adjust the risk level of individual bird species according to the actual situation of the airport:
    1.
    If the bird is a raptor, increase the risk level by one.

    2.
    The risk level for bird species that are seen crossing a runway or passing through the sky above the runway more than three times should be increased by one.

    Raptors fly fast, and collisions with airplanes can have very serious consequences. Among the birdstrike events with the largest record of birdstrike impact energy from 2007 to 2014, half of them were raptors. However, because their weight is actually low compared to birds like ducks, and their solitary habits, the risk level calculated by this method is often lower than the actual risk, so the risk level of the raptor is increased by one level. Most birdstrikes occur when the aircraft takes off and lands. If the bird’s movement often crosses the runway or the nearby sky, it is more likely to cross an aircraft’s flight trajectory, and therefore is very dangerous for the aircraft. For this reason, when a bird species is seen crossing the runway and flying over the top of the runway three times, the risk level of that species should be increased by one.
    Each airport should adjust their assessments based on locally collected empirical data on strike likelihood and severity as well as ongoing bird monitoring at the airport and its surrounding environment. More

  • in

    Lipid content and stable isotopes of zooplankton during five winters around the northern Antarctic Peninsula

    Survey area
    The U.S. AMLR Program conducted five winter surveys (August and September 2012–2016) around the northern Antarctic Peninsula aboard the U.S. National Science Foundation research vessel/ice breaker (RVIB) Nathaniel B. Palmer (Fig. 1)19. We surveyed a historical grid of 110 fixed stations located 20–40 km apart around the northern Antarctic Peninsula and South Shetland Islands, which was divided into four sampling areas: the Elephant Island Area (EI; 43,865 km2), the South Area (the Bransfield Strait, SA; 24,479 km2), the West Area (the west shelf immediately north of Livingston and King George Islands, WA; 38,524 km2), and the Joinville Island Area (JI; 18,151 km2). In 2016, we also surveyed the Gerlache Strait (GS; 24,479 km2).
    At-sea sampling
    Detailed methods for all at-sea sampling by the U.S. AMLR Program have been previously reported19. At each sampling station, we performed a Conductivity-Temperature-Depth (CTD) cast to 750 m or to within 10 m of the bottom in shallower areas (SBE9/11; Sea-Bird Electronics). The CTD rosette was equipped with 24 10 l Niskin bottles triggered to close on the upcast at 750, 200, 100, 75, 50, 40, 30, 20, 15, and 5 m. We defined the Upper Mixed Layer (UML) depth (m) as the depth at which the density of the water changed by 0.05 kg m−3 relative to the mean density of the upper 10 m of the water column20. UML temperature and salinities were defined to be means over the depth range of the UML. We defined daylight conditions according to three categories. Day (D) was defined as one hour after local sunrise to one hour before local sunset; night (N) was defined as one hour after sunset to one hour before sunrise; and Twilight (T) was defined as one hour before and after sunrise and sunset.
    Chlorophyll-a (hereafter chl-a) was determined from water samples collected between 5 and 200 m21,22. Samples (285 ml) were filtered at a pressure differential of  More

  • in

    Trait convergence and trait divergence in lake phytoplankton reflect community assembly rules

    1.
    Götzenberger, L. et al. Ecological assembly rules in plant communities—Approaches, patterns and prospects. Biol. Rev. 87, 111–127. https://doi.org/10.1111/j.1469-185X.2011.00187.x (2012).
    Article  PubMed  Google Scholar 
    2.
    Díaz, S. et al. Incorporating plant functional diversity effects in ecosystem service assessments. P. Natl. Axad. Sci. USA 104, 20684–20689. https://doi.org/10.1073/pnas.0704716104 (2007).
    ADS  Article  Google Scholar 

    3.
    McGill, B. J., Enquist, B. J., Weiher, E. & Westoby, M. Rebuilding community ecology from functional traits. Trends Ecol. Evol. 21, 178–185. https://doi.org/10.1016/j.tree.2006.02.002 (2006).
    Article  PubMed  Google Scholar 

    4.
    Nobel, I.R. & Slatyer, R.O. Post-fire succession of plants in Mediterranean ecosystems in: Proceedings of the Symposium on the Environmental Consequences of Fire and Fuel Management in Mediterranean Ecosystems (eds. Mooney, H.A. & Conrad, C.E.) 27–36 (California Palo Alto, 1977).

    5.
    Cornwell, W. K. & Ackerly, D. D. Community assembly and shifts in plant trait distributions across an environmental gradient in coastal California. Ecol. Monogr. 79, 109–126. https://doi.org/10.1890/07-1134.1 (2009).
    Article  Google Scholar 

    6.
    Kraft, N. J. B., Valencia, R. & Ackerly, D. D. Functional traits and niche-based tree community assembly in an amazonian forest. Science 322, 580–582. https://doi.org/10.1126/science.1160662 (2008).
    ADS  CAS  Article  PubMed  Google Scholar 

    7.
    MacArthur, R. & Levins, R. The limiting similarity, convergence, and divergence of coexisting species. Am. Nat. 101, 377–385. https://doi.org/10.1086/282505 (1967).
    Article  Google Scholar 

    8.
    Pásztor, L., Botta-Dukát, Z., Magyar, G., Czárán, T. & Meszéna, G. Theory-Based Ecology: A Darwinian Approach (Oxford University Press, Oxford, 2016).
    Google Scholar 

    9.
    Diamond, J.M. Assembly of Species Communities in Ecology and Evolution of Communities (eds. Cody, M.L. & Diamond, J.M.) 342–444 (Belknap Press, 1975).

    10.
    Litchman, E., Klausmeier, C. A., Schofield, O. M. & Falkowski, P. G. The role of functional traits and trade-offs in structuring phytoplankton communities: Scaling from cellular to ecosystem level. Ecol. Lett. 10, 1170–1181. https://doi.org/10.1111/j.1461-0248.2007.01117.x (2007).
    Article  PubMed  Google Scholar 

    11.
    Reynolds, C. S., Huszár, V., Kruk, C., Naselli-Flores, L. & Melo, S. Towards a functional classification of the freshwater phytoplankton. J. Plankton Res. 24, 417–428. https://doi.org/10.1093/plankt/24.5.417 (2002).
    Article  Google Scholar 

    12.
    Salmaso, N. & Padisák, J. Morpho-functional groups and phytoplankton development in two deep lakes (Lake Garda, Italy and Lake Stechlin, Germany). Hydrobiologia 578, 97–112. https://doi.org/10.1007/s10750-006-0437-0 (2007).
    Article  Google Scholar 

    13.
    Salmaso, N., Naselli-Flores, L. & Padisák, J. Functional classifications and their application in phytoplankton ecology. Freshw. Biol. 60, 603–619. https://doi.org/10.1111/fwb.12520 (2015).
    Article  Google Scholar 

    14.
    Borics, G., Tóthmérész, B., Lukács, B. A. & Várbíró, G. Functional groups of phytoplankton shaping diversity of shallow lake ecosystems. Hydrobiologia 698, 251–262. https://doi.org/10.1007/s10750-012-1129-6 (2012).
    Article  Google Scholar 

    15.
    Vellend, M. Conceptual synthesis in community ecology. Q. Rev. Biol. 85, 183–206. https://doi.org/10.1086/652373 (2010).
    Article  PubMed  PubMed Central  Google Scholar 

    16.
    Padisák, J., Vasas, G. & Borics, G. Phycogeography of freshwater phytoplankton: Traditional knowledge and new molecular tools. Hydrobiologia 764, 3–27. https://doi.org/10.1007/s10750-015-2259-4 (2016).
    CAS  Article  Google Scholar 

    17.
    Padisák, J. Cylindrospermopsis raciborskii (Woloszynska) Seenayya et Subba Raju, an expanding, highly adaptive cyanobacterium: worldwide distribution and review of its ecology. Arch. Hydrobiol. 107, 563–593 (1997).
    Google Scholar 

    18.
    Méndez, V., Assaf, M., Masó-Puigdellosas, A., Campos, D. & Horsthemke, W. Demographic stochasticity and extinction in populations with Allee effect. Phys. Rev. E. 99, 022101. https://doi.org/10.1103/PhysRevE.99.022101 (2019).
    ADS  Article  Google Scholar 

    19.
    Parvinen, K., Dieckmann, U., Gyllenberg, M. & Metz, J. A. Evolution of dispersal in metapopulations with local density dependence and demographic stochasticity. J. Evolut. Biol. 16, 143–153. https://doi.org/10.1046/j.1420-9101.2003.00478.x (2003).
    CAS  Article  Google Scholar 

    20.
    Borics, G., Abonyi, A., Salmaso, N. & Ptacnik, R. Freshwater phytoplankton diversity: Models, drivers and implications for ecosystem properties. Hydrobiologia https://doi.org/10.1007/s10750-020-04332-9 (2020).
    Article  PubMed  PubMed Central  Google Scholar 

    21.
    Naselli-Flores, L., Padisák, J., Dokulil, M. T. & Chorus, I. Equilibrium/steady-state concept in phytoplankton ecology. Hydrobiologia 502, 395–403. https://doi.org/10.1023/B:HYDR.0000004297.52645.59 (2003).
    Article  Google Scholar 

    22.
    Weiher, E. & Keddy, P.A. Assembly rules, null models, and trait dispersion: New questions from old patterns. Oikos 74, 159–164 (1995). https://www.jstor.org/stable/3545686.

    23.
    Coyle, J. R. et al. Using trait and phylogenetic diversity to evaluate the generality of the stress-dominance hypothesis in eastern North American tree communities. Ecography 37, 814–826. https://doi.org/10.1111/ecog.00473 (2014).
    Article  Google Scholar 

    24.
    Baastrup-Spohr, L., Sand-Jensen, K., Nicolajsen, S. V. & Bruun, H. H. From soaking wet to bone dry: Predicting plant community composition along a steep hydrological gradient. J. Veg. Sci. 26, 619–630. https://doi.org/10.1111/jvs.12280 (2015).
    Article  Google Scholar 

    25.
    Lhotsky, B. et al. Changes in assembly rules along a stress gradient from open dry grasslands to wetlands. J. Ecol. 104, 507–517. https://doi.org/10.1111/1365-2745.12532 (2016).
    Article  Google Scholar 

    26.
    Butterfield, B. J., Bradford, J. B., Munson, S. M. & Gremer, J. R. Aridity increases below-ground niche breadth in grass communities. Plant Ecol. 218, 385–394. https://doi.org/10.1007/s11258-016-0696-4 (2017).
    Article  Google Scholar 

    27.
    Gastauer, M., Saporetti-Junior, A. W., Valladares, F. & Meira-Neto, J. A. Phylogenetic community structure reveals differences in plant community assembly of an oligotrophic white-sand ecosystem from the Brazilian Atlantic Forest. Acta Bot. Bras. 31, 531–538. https://doi.org/10.1590/0102-33062016abb0442 (2017).
    Article  Google Scholar 

    28.
    Chapman, J. & McEwan, R. The role of environmental filtering in structuring appalachian tree communities: Topographic influences on functional diversity are mediated through soil characteristics. Forests 9, 19. https://doi.org/10.3390/f9010019 (2018).
    Article  Google Scholar 

    29.
    Lukács, B. A. et al. Carbon forms, nutrients and water velocity filter hydrophyte and riverbank species differently: A trait-based study. J. Veg. Sci. 30, 471–484 (2019).
    Article  Google Scholar 

    30.
    Kuczynski, L. & Grenouillet, G. Community disassembly under global change: Evidence in favor of the stress-dominance hypothesis. Glob. Chang. Biol. 24, 4417–4427. https://doi.org/10.1111/gcb.14320 (2018).
    ADS  Article  PubMed  Google Scholar 

    31.
    Patrick, L. E. & Stevens, R. D. Phylogenetic community structure of North American desert bats: Influence of environment at multiple spatial and taxonomic scales. J. Anim. Ecol. 85, 1118–1130. https://doi.org/10.1111/1365-2656.12529 (2016).
    Article  PubMed  Google Scholar 

    32.
    Lopez, B. et al. A new framework for inferring community assembly processes using phylogenetic information, relevant traits and environmental gradients. One Ecosyst. 1, e9501. https://doi.org/10.3897/oneeco.1.e9501 (2016).
    Article  Google Scholar 

    33.
    Ács, E. et al. Trait-based community assembly of epiphytic diatoms in saline astatic ponds: a test of the stress-dominance hypothesis. Sci. Rep. 9(1), 15749 (2019).
    ADS  Article  PubMed  PubMed Central  Google Scholar 

    34.
    Downing, J. A. & McCauley, E. The nitrogen:phosphorus relationship in lakes. Limnol. Oceanogr. 37, 936–945 (1992).
    ADS  CAS  Article  Google Scholar 

    35.
    Phillips, G. et al. A phytoplankton trophic index to assess the status of lakes for the Water Framework Directive. Hydrobiologia 704, 75–95. https://doi.org/10.1007/s10750-012-1390-8 (2013).
    Article  Google Scholar 

    36.
    Padisák, J. et al. Dominant species, functional assemblages and frequency of equilibrium phases in late summer phytoplankton assemblages in Hungarian small shallow lakes. Hydrobiologia 502, 157–168. https://doi.org/10.1023/B:HYDR.0000004278.10887.40 (2003).
    Article  Google Scholar 

    37.
    HilleRisLambers, J., Adler, P. B., Harpole, W. S., Levine, J. M. & Mayfield, M. M. Rethinking community assembly through the lens of coexistence theory. Annu. Rev. Ecol. Evol. Syst. 43, 227–248. https://doi.org/10.1146/annurev-ecolsys-110411-160411 (2012).
    Article  Google Scholar 

    38.
    Huisman, J., van Oostveen. P. & Weissing, F.J. Species dynamics in phytoplankton blooms: incomplete mixing and competition for light. Am. Nat. 154, 46–68, https://doi.org/10.1086/303220 (1999).

    39.
    Bird, D. F. & Kalff, J. Bacterial grazing by planktonic lake algae. Science 231, 493–495. https://doi.org/10.1126/science.231.4737.493 (1986).
    ADS  CAS  Article  Google Scholar 

    40.
    Stoecker, D. K. Mixotrophy among Dinoflagellates. J. Eukaryot. Microbiol. 46, 397–401. https://doi.org/10.1111/j.1550-7408.1999.tb04619.x (1999).
    Article  Google Scholar 

    41.
    Grime, J.P. Plant Strategies, Vegetation Processes, and Ecosystem Properties. (John Wiley and Sons, 2001). ISBN 0-471-49601-4.

    42.
    Navas, M. & Violle, C. Plant traits related to competition: How do they shape the functional diversity of communities?. Commun. Ecol. 10, 131–137. https://doi.org/10.1556/ComEc.10.2009.1.15 (2009).
    Article  Google Scholar 

    43.
    Reynolds, C. S. The Ecology of Phytoplankton (Cambridge University Press, Cambridge, 2006).
    Google Scholar 

    44.
    Borics, G., Grigorszky, I., Szabó, S. & Padisák, J. Phytoplankton associations in a small hypertrophic fishpond in East Hungary during a change from bottom-up to top-down control. Hydrobiologia 424, 79–90. https://doi.org/10.1023/A:1003948827254 (2000).
    Article  Google Scholar 

    45.
    Mason, N. W. H., de Bello, F., Doležal, J. & Lepš, J. Niche overlap reveals the effects of competition, disturbance and contrasting assembly processes in experimental grassland communities. J. Ecol. 99, 788–796. https://doi.org/10.1111/j.1365-2745.2011.01801.x (2011).
    Article  Google Scholar 

    46.
    Dobosi, Z. & Felméry, L. Climatology, ELTE TTK, Nemzeti Tankönyvkiadó, p. 500 (in Hungarian).

    47.
    Mihevc, A., Prelovšek, M. & Hajna, N.Z. Introduction to the Dinaric Karst. Inštitut za raziskovanje krasa ZRC SAZU. (2010).

    48.
    Borics, G. et al. Phytoplankton-based shallow lake types in the Carpathian basin: Steps towards a bottom-up typology. Fund. Appl. Limnol. 184, 23–34. https://doi.org/10.1127/1863-9135/2014/0518 (2014).
    Article  Google Scholar 

    49.
    Borics, G., Abonyi, A., Várbíró, G., Padisák, J. &T-Krasznai, E. Lake stratification in the Carpathian basin and its interesting biological consequences. Inland Waters 5, 173–186, https://doi.org/10.5268/IW-5.2.702 (2015).

    50.
    Utermöhl, H. Zur Vervollkommnung der quantitative Phytolankton-Methodik. Mitt. Int. Verein. Limnol. 9, 1–38. https://doi.org/10.1080/05384680.1958.11904091 (1958).
    Article  Google Scholar 

    51.
    Török, P. et al. Functional diversity supports the biomass–Diversity humped-back relationship in phytoplankton assemblages. Funct. Ecol. 30, 1593–1602. https://doi.org/10.1111/1365-2435.12631 (2016).
    Article  Google Scholar 

    52.
    Hillebrand, H., Dürselen, C. D., Kirschtel, D., Pollingher, U. & Zohary, T. Biovolume calculation for pelagic and benthic microalgae. J. Phycol. 35, 403–424. https://doi.org/10.1046/j.1529-8817.1999.3520403.x (1999).
    Article  Google Scholar 

    53.
    MSZ ISO 10260:1993. Water Quality. Measurement of Biochemical Parameters. Spectrometric Determination of the Chlorophyll-a Concentration.

    54.
    MSZ EN ISO 6878:2004. Water Quality. Determination of Phosphorus. Ammonium Molybdate Spectrometric Method.

    55.
    ISO 11905-1:1997. Water Quality. Determination of Nitrogen. Part 1: Method Using Oxidative Digestion with Peroxodisulfate.

    56.
    MSZ ISO 6060:1991.Water Quality. Determination of the Chemical Oxygen Demand.

    57.
    Botta-Dukát, Z. Cautionary note on calculating standardized effect size (SES) in randomization test. Commun. Ecol. 19, 77–83, https://doi.org/10.1556/168.2018.19.1.8 (2018).

    58.
    Botta-Dukát, Z. Rao’s quadratic entropy as a measure of functional diversity based on multiple traits. J. Veg. Sci. 16, 533–540. https://doi.org/10.1111/j.1654-1103.2005.tb02393.x (2005).
    Article  Google Scholar 

    59.
    Legendre, P. & Legendre, L.F. Numerical Ecology Vol. 24. (Elsevier, 2012).

    60.
    Götzenberger, L. et al. Which randomizations detect convergence and divergence in trait‐based community assembly? A test of commonly used null models. J. Veg. Sci. 27, https://doi.org/10.1111/jvs.12452, 1275–1287.

    61.
    Botta-Dukát, Z. & Czúcz, B. Testing the ability of functional diversity indices to detect trait convergence and divergence using individual-based simulation. Methods Ecol. Evol. 7, 114–126. https://doi.org/10.1111/2041-210X.12450 (2016).
    Article  Google Scholar 

    62.
    Garnier E. et al. Plant functional markers capture ecosystem properties during secondary succession. Ecology 85, e2630–2637, https://doi.org/10.1890/03-0799.

    63.
    Dray, S. & Dufour, A.B. The ade4 package: Implementing the duality diagram for ecologists. J. Stat. Softw. 22, 1–20, https://doi.org/10.18637/jss.v022.i04.

    64.
    Oksanen J. et al. Package ‘vegan’. Community Ecology Package, Version, 2(9) (2013).

    65.
    Team R. Core (2017) R: A Language and Environment for Statistical Computing (R Foundation for Statistical Computing, Vienna, 2017).
    Google Scholar  More

  • in

    DISPERSE, a trait database to assess the dispersal potential of European aquatic macroinvertebrates

    1.
    Bohonak, A. J. & Jenkins, D. G. Ecological and evolutionary significance of dispersal by freshwater invertebrates. Ecol. Lett. 6, 783–796 (2003).
    Article  Google Scholar 
    2.
    Clobert, J., Baguette, M., Benton, T. G. & Bullock, J. M. Dispersal Ecology and Evolution (Oxford Univ. Press, 2012).

    3.
    Heino, J. et al. Metacommunity organisation, spatial extent and dispersal in aquatic systems: Patterns, processes and prospects. Freshw. Biol. 60, 845–869 (2015).
    Article  Google Scholar 

    4.
    Barton, P. S. et al. Guidelines for using movement science to inform biodiversity policy. Environ. Manage. 56, 791–801 (2015).
    ADS  Article  Google Scholar 

    5.
    Heino, J. et al. Integrating dispersal proxies in ecological and environmental research in the freshwater realm. Environ. Rev. 25, 334–349 (2017).
    Article  Google Scholar 

    6.
    Rundle, S. D., Bilton, D. T. & Foggo, A. in Body Size: The Structure and Function of Aquatic Ecosystems (eds. Hildrew, A. G., Raffaelli, D. G. & Edmonds-Brown, R.) 186–209 (Cambridge Univ. Press, 2007).

    7.
    Macneale, K. H., Peckarsky, B. L. & Likens, G. E. Stable isotopes identify dispersal patterns of stonefly populations living along stream corridors. Freshw. Biol. 50, 1117–1130 (2005).
    Article  Google Scholar 

    8.
    Troast, D., Suhling, F., Jinguji, H., Sahlén, G. & Ware, J. A global population genetic study of Pantala flavescens. PLoS One 11, e0148949 (2016).
    Article  CAS  Google Scholar 

    9.
    French, S. K. & McCauley, S. J. The movement responses of three libellulid dragonfly species to open and closed landscape cover. Insect Conserv. Divers. 12, 437–447 (2019).
    Article  Google Scholar 

    10.
    Arribas, P. et al. Dispersal ability rather than ecological tolerance drives differences in range size between lentic and lotic water beetles (Coleoptera: Hydrophilidae). J. Biogeogr. 39, 984–994 (2012).
    Article  Google Scholar 

    11.
    Lancaster, J. & Downes, B. J. Dispersal traits may reflect dispersal distances, but dispersers may not connect populations demographically. Oecologia 184, 171–182 (2017).
    ADS  Article  Google Scholar 

    12.
    Lancaster, J. & Downes, B. J. Aquatic Entomology (Oxford Univ. Press, 2013).

    13.
    Stevens, V. M. et al. Dispersal syndromes and the use of life-histories to predict dispersal. Evol. Appl. 6, 630–642 (2013).
    Article  Google Scholar 

    14.
    Outomuro, D. & Johansson, F. Wing morphology and migration status, but not body size, habitat or Rapoport’s rule predict range size in North-American dragonflies (Odonata: Libellulidae). Ecography 42, 309–320 (2019).
    Article  Google Scholar 

    15.
    Tonkin, J. D. et al. The role of dispersal in river network metacommunities: patterns, processes, and pathways. Freshw. Biol. 63, 141–163 (2018).
    Article  Google Scholar 

    16.
    Brown, B. L. & Swan, C. M. Dendritic network structure constrains metacommunity properties in riverine ecosystems. J. Anim. Ecol. 79, 571–580 (2010).
    CAS  Article  Google Scholar 

    17.
    Wikelski, M. et al. Simple rules guide dragonfly migration. Biol. Lett. 2, 325–329 (2006).
    Article  Google Scholar 

    18.
    Schmidt-Kloiber, A. & Hering, D. An online tool that unifies, standardises and codifies more than 20,000 European freshwater organisms and their ecological preferences. Ecol. Indic. 53, 271–282, www.freshwaterecology.info (2015).
    Article  Google Scholar 

    19.
    Serra, S. R. Q., Cobo, F., Graça, M. A. S., Dolédec, S. & Feio, M. J. Synthesising the trait information of European Chironomidae (Insecta: Diptera): towards a new database. Ecol. Indic. 61, 282–292 (2016).
    Article  Google Scholar 

    20.
    Tachet, H., Richoux, P., Bournaud, M. & Usseglio-Polatera, P. Invertébrés d’Eau Douce: Systématique, Biologie, Écologie (CNRS Éditions, 2010).

    21.
    Vieira, N. K. M. et al. A Database of Lotic Invertebrate Traits for North America (U.S. Geological Survey Data Series 187, 2006).

    22.
    Chevenet, F., Dolédec, S. & Chessel, D. A fuzzy coding approach for the analysis of long-term ecological data. Freshw. Biol. 31, 295–309 (1994).
    Article  Google Scholar 

    23.
    Schmera, D., Podani, J., Heino, J., Erös, T. & Poff, N. L. R. A proposed unified terminology of species traits in stream ecology. Freshw. Sci. 34, 823–830 (2015).
    Article  Google Scholar 

    24.
    Lancaster, J., Downes, B. J. & Arnold, A. Lasting effects of maternal behaviour on the distribution of a dispersive stream insect. J. Anim. Ecol. 80, 1061–1069 (2011).
    Article  Google Scholar 

    25.
    Jenkins, D. G. et al. Does size matter for dispersal distance? Glob. Ecol. Biogeogr. 16, 415–425 (2007).
    Article  Google Scholar 

    26.
    Harrison, R. G. Dispersal polymorphisms in insects. Annu. Rev. Ecol. Syst. 11, 95–118 (1980).
    Article  Google Scholar 

    27.
    Graham, E. S., Storey, R. & Smith, B. Dispersal distances of aquatic insects: upstream crawling by benthic EPT larvae and flight of adult Trichoptera along valley floors. New Zeal. J. Mar. Freshw. Res. 51, 146–164 (2017).
    Article  Google Scholar 

    28.
    Hoffsten, P. O. Site-occupancy in relation to flight-morphology in caddisflies. Freshw. Biol. 49, 810–817 (2004).
    Article  Google Scholar 

    29.
    Bonada, N. & Dolédec, S. Does the Tachet trait database report voltinism variability of aquatic insects between Mediterranean and Scandinavian regions? Aquat. Sci. 80, 1–11 (2018).
    Article  Google Scholar 

    30.
    Sarremejane, R. et al. DISPERSE, a trait database to assess the dispersal potential of aquatic macroinvertebrates. figshare https://doi.org/10.6084/m9.figshare.c.5000633 (2020).

    31.
    Lévêque, C., Balian, E. V. & Martens, K. An assessment of animal species diversity in continental waters. Hydrobiologia 542, 39–67 (2005).
    Article  Google Scholar 

    32.
    Green, A. J. & Figuerola, J. Recent advances in the study of long-distance dispersal of aquatic invertebrates via birds. Divers. Distrib. 11, 149–156 (2005).
    Article  Google Scholar 

    33.
    Maasri, A. A global and unified trait database for aquatic macroinvertebrates: the missing piece in a global approach. Front. Environ. Sci. 7, 1–3 (2019).
    Article  Google Scholar 

    34.
    Cañedo-Argüelles, M. et al. Dispersal strength determines meta-community structure in a dendritic riverine network. J. Biogeogr. 42, 778–790 (2015).
    Article  Google Scholar 

    35.
    Datry, T. et al. Metacommunity patterns across three Neotropical catchments with varying environmental harshness. Freshw. Biol. 61, 277–292 (2016).
    Article  Google Scholar 

    36.
    Swan, C. M. & Brown, B. L. Metacommunity theory meets restoration: isolation may mediate how ecological communities respond to stream restoration. Ecol. Appl. 27, 2209–2219 (2017).
    Article  Google Scholar 

    37.
    Sarremejane, R., Mykrä, H., Bonada, N., Aroviita, J. & Muotka, T. Habitat connectivity and dispersal ability drive the assembly mechanisms of macroinvertebrate communities in river networks. Freshw. Biol. 62, 1073–1082 (2017).
    Article  Google Scholar 

    38.
    Jacobson, B. & Peres-Neto, P. R. Quantifying and disentangling dispersal in metacommunities: How close have we come? How far is there to go? Landsc. Ecol. 25, 495–507 (2010).
    Article  Google Scholar 

    39.
    Sarremejane, R. et al. Do metacommunities vary through time? Intermittent rivers as model systems. J. Biogeogr. 44, 2752–2763 (2017).
    Article  Google Scholar 

    40.
    Datry, T., Moya, N., Zubieta, J. & Oberdorff, T. Determinants of local and regional communities in intermittent and perennial headwaters of the Bolivian Amazon. Freshw. Biol. 61, 1335–1349 (2016).
    Article  Google Scholar 

    41.
    Cid, N. et al. A metacommunity approach to improve biological assessments in highly dynamic freshwater ecosystems. Bioscience 70, 427–438 (2020).
    Article  Google Scholar 

    42.
    Datry, T., Bonada, N. & Heino, J. Towards understanding the organisation of metacommunities in highly dynamic ecological systems. Oikos 125, 149–159 (2016).
    Article  Google Scholar 

    43.
    Hermoso, V., Cattarino, L., Kennard, M. J., Watts, M. & Linke, S. Catchment zoning for freshwater conservation: refining plans to enhance action on the ground. J. Appl. Ecol. 52, 940–949 (2015).
    Article  Google Scholar 

    44.
    Thuiller, W. et al. A road map for integrating eco-evolutionary processes into biodiversity models. Ecol. Lett. 16, 94–105 (2013).
    Article  Google Scholar 

    45.
    Mendes, P., Velazco, S. J. E., de Andrade, A. F. A. & De Marco, P. Dealing with overprediction in species distribution models: How adding distance constraints can improve model accuracy. Ecol. Model. 431, 109180 (2020).
    Article  Google Scholar 

    46.
    Willis, S. G. et al. Integrating climate change vulnerability assessments from species distribution models and trait-based approaches. Biol. Conserv. 190, 167–178 (2015).
    Article  Google Scholar 

    47.
    Cooper, J. C. & Soberón, J. Creating individual accessible area hypotheses improves stacked species distribution model performance. Glob. Ecol. Biogeogr. 27, 156–165 (2018).
    Article  Google Scholar 

    48.
    Markovic, D. et al. Europe’s freshwater biodiversity under climate change: distribution shifts and conservation needs. Divers. Distrib. 20, 1097–1107 (2014).
    Article  Google Scholar 

    49.
    Bush, A. & Hoskins, A. J. Does dispersal capacity matter for freshwater biodiversity under climate change? Freshw. Biol. 62, 382–396 (2017).
    Article  Google Scholar 

    50.
    Bohonak, A. J. Dispersal, gene flow, and population structure. Q. Rev. Biol. 74, 21–45 (1999).
    CAS  Article  Google Scholar 

    51.
    Dijkstra, K.-D. B., Monaghan, M. T. & Pauls, S. U. Freshwater biodiversity and aquatic insect diversification. Annu. Rev. Entomol. 59, 143–163 (2014).
    CAS  Article  Google Scholar 

    52.
    Múrria, C. et al. Local environment rather than past climate determines community composition of mountain stream macroinvertebrates across Europe. Mol. Ecol. 26, 6085–6099 (2017).
    Article  CAS  Google Scholar 

    53.
    Statzner, B. & Bêche, L. A. Can biological invertebrate traits resolve effects of multiple stressors on running water ecosystems? Freshw. Biol. 55, 80–119 (2010).
    Article  Google Scholar 

    54.
    Strayer, D. L. & Dudgeon, D. Freshwater biodiversity conservation: recent progress and future challenges. J. North Am. Benthol. Soc. 29, 344–358 (2010).
    Article  Google Scholar 

    55.
    Reid, A. J. et al. Emerging threats and persistent conservation challenges for freshwater biodiversity. Biol. Rev. 94, 849–873 (2019).
    Article  Google Scholar 

    56.
    R Core Team. R: A language and environment for statistical computing. https://www.r-project.org/ (2020).

    57.
    Dray, S. & Dufour, A.-B. The ade4 package: Implementing the duality diagram for ecologists. J. Stat. Softw. 1, 1–20 (2007).
    Google Scholar  More

  • in

    Detritivore conversion of litter into faeces accelerates organic matter turnover

    Detritivore and leaf litter collection
    We collected six phylogenetically diverse species of detritivores in various areas of the Scottish Lowlands in May and June 2018, including three millipede species (Diplopoda), two woodlouse species (Crustacea) and one snail species (Gastropoda). Millipede species include the common pill millipede (Glomeris marginata (Villers, 1789)) collected near Peebles, UK (55°38′45.8″N, 3°07′55.4″W), the striped millipede (Ommatoiulus sabulosus (Linnaeus, 1758)) collected near Dunfermline, UK (56°02′23.7″N, 3°19′49.2″W) and the white-legged millipede (Tachypodoiulus niger (Leach, 1815)) collected near Dundee, UK (56°32′08.5″N, 3°01′51.9″W). Woodlouse species include the common pill woodlouse (Armadillidium vulgare (Latreille, 1804)) collected near Dunfermline, UK (56°01′35.3″N 3°23′14.1″W) and the common rough woodlouse (Porcellio scaber (Latreille, 1804) collected in Stirling, UK (56°07′26.7″N, 3°55′51.2″W). The snail species was the brown-lipped snail (Cepaea nemoralis (Linnaeus, 1758)) collected in Stirling, UK (56°08′07.3″N, 3°55′16.3″W). These species are common in diverse ecosystems across Mediterranean and temperate ecosystems in Europe, where they feed on decomposing litter and produce large amounts of faeces16,38,39,40. Detritivores were kept in plastic boxes and fed with moist litter from various tree species from their respective collection sites before the start of the experiment.
    To obtain a gradient of leaf litter quality, we collected leaf litter from six deciduous broadleaf tree species in the Scottish Lowlands. These species include sycamore maple (Acer pseudoplatanus, L.), horse chestnut (Aesculus hippocastanum, L.), common hazel (Corylus avellana, L.), European beech (Fagus sylvatica, L.), English oak (Quercus robur, L.) from a woodland near Dundee, UK (56°32′08.5″N, 3°01′51.9″W) and lime (Tilia platyphyllos, L.) from a woodland in Stirling, UK (56°08′29.5″N, 3°55′14.2″W). Because detritivores are most active in spring and summer in these ecosystems, they feed on partially decomposed litter, which they prefer over freshly fallen litter (David and Gillon8). We thus collected leaf litter from the forest floor in May 2018, air-dried it and stored it in cardboard boxes until use.
    Faeces production
    To compare the quality and decomposability of leaf litter with faeces derived from the same litter and produced by diverse detritivore species, we set up two series of boxes for the production of the needed material. In the first of these series, we placed each detritivore species together with each litter species to produce the 36 different faeces types (Fig. 1; 6 litter species × 6 detritivore species = 36 faeces types). The second of these series contained the litter species only without any detritivores to produce intact litter from each tree species (6 litter species) under the same conditions for the same amount of time. In total, 42 different substrates were generated. To do so, we placed ca. 30 g of air-dry leaf litter from each species separately in plastic boxes (30 cm × 22 cm × 5.5 cm) to which we added ca. 50 individuals from each detritivore species separately or no detritivore for the intact litter treatment. We sprayed the litter with water to optimise litter moisture for detritivore consumption while avoiding water accumulation at the bottom of the boxes. We kept the boxes at room temperature (ca. 20 °C) for 4 weeks and collected the produced faeces/intact litter twice a week. For the faeces, we placed the content of each box in a large bucket and gently agitated to let detritivores and faeces fall to the bottom of the bucket. After collecting the faeces, we placed all the leaf litter and detritivores back into their boxes and sprayed the litter with water to keep moisture conditions constant. For the intact litter treatment, we followed the same procedure but collected just three random leaves out of the buckets. After each collection step, the combination-specific pools of leaf litter and faeces were dried at 30 °C. At the end of the faeces production period, we manually removed small leaf litter fragments from all combination-specific pools of faeces. Additionally, because detritivores feed on leaf lamina and leave leaf veins mostly uneaten6, we cut out the veins from the species-specific pools of intact leaf litter. This was done to ensure the comparability of quality and decomposability between faeces and intact litter.
    Litter and faeces quality
    To evaluate the effect of litter conversion into detritivore faeces on organic matter quality, we compared the quality of faeces to that of intact litter by measuring a series of physical and chemical quality parameters on all 42 substrates (6 litter species + 36 faeces types). Chemical characteristics included total carbon (C) and nitrogen (N) concentrations, DOC and TDN concentrations, total tannin concentrations, and 13C solid-state NMR spectra. Physical characteristics included WHC and specific area (surface area per unit of mass). Prior to these measurements, we drew three subsamples from each pool of substrate type. A part of each subsample was ground using a ball mill (TissueLyser II, Qiagen) to measure total C, N and tannin concentration and generate NMR spectra. The other part of each subsample was kept intact and used for all other measurements. All measurements were thus done on these three subsamples per substrate type, except for NMR spectra that were measured once per substrate type on a sample made by pooling all three ground subsamples. This pooling was necessary to obtain a sample large enough for the NMR analyses. Total C and N concentrations were measured with a flash CHN elemental analyser (Flash Smart, ThermoScientific). To measure DOC and TDN, we extracted leachates by placing ca. 30 mg of air-dried material with 25 ml of deionised water in 50 ml Falcon tubes and agitating the tubes horizontally on a reciprocal shaker for 1 h. Water extracts were then filtered through 0.45-μm cellulose nitrate filters to isolate the leachate fraction. Concentrations of DOC and TDN in leachates were measured with a TOC analyser (Shimadzu, Kyoto, Japan) equipped with a supplementary module for N. Tannin concentrations were measured with the protein-precipitable phenolics microplate assay, a microplate protocol adapted from Hagerman and Butler41. We obtained 13C-NMR spectra by applying 13C cross-polarisation magic angle spinning NMR spectroscopy using a 200 MHz spectrometer (Bruker, Billerica, USA). The samples were spun in 7 mm zirconium dioxide rotors at 6.8 kHz with an acquisition time of 0.01024 s. To avoid Hartmann–Hahn mismatches, a ramped 1H impulse was applied during a contact time of 1 ms. We applied a delay time of 2.0 s and the number of scans was set to 1500, yet some of the samples required longer measurements due to the low amount of sample material; in this case, we multiplied the number of scans to 3000, 6000 or 15000. As reference for the chemical shift, tetramethylsilane was used (0 ppm). We used the following chemical shift regions to integrate the spectra: −10–45 ppm alkyl C, 45–110 ppm O/N alkyl C, 110–160 ppm aromatic C, and 160–220 ppm carboxylic C. We measured the WHC by placing ca. 15 mg of air-dried intact material with 1.5 ml of deionised water in 2 ml Eppendorf tubes, agitating the tubes horizontally on a reciprocal shaker for 2 h, retrieving the material and placing it on a Whatman filter to remove excess water, weighing the wet material and reweighing it after drying at 65 °C for 48 h. We measured the specific area of leaf litter, faecal pellets and faeces particles from photographs using a stereomicroscope (ZEISS STEMI 508). For leaf litter and faecal pellets, we took photographs of ca. 20 mg of air-dried intact material. To visualise faeces particles, we weighed ca. 1 mg of air-dried faecal pellets and placed them in a beaker with 20 ml of deionised water for 2 h, allowing complete dissolution of the faecal pellets. We then filtered the faeces particles and photographed the filters under a stereomicroscope. Dimensions of each litter pieces and faecal pellets/faeces particles were measured using the image analysis software (ImageJ, version 1.46r). For all substrate types, we divided the calculated surface area by the dry mass of the sample to obtain the specific area.
    Faeces and litter decomposition parameters
    To evaluate the effect of litter conversion into detritivore faeces on C and N cycling, we compared the C and N loss of faeces to that of intact litter by incubating all 42 substrates in microcosms under controlled conditions for 6 months (180 days). Microcosms consisted of 250-ml plastic containers filled with 90 mg of air-dry soil collected from a temperate grassland (56°8′40.1″N, 3°54′50.9″W). We chose this soil to avoid any home-field advantage effect as this soil did not receive litter input from any of the studied tree species and none of the selected soil animals were present at this site. About 120 mg of each substrate were placed separately within a small polyvinyl chloride tube (30 mm diameter × 30 mm height) closed in the bottom with a 100-µm mesh and left open on the top. Each tube was then placed on top of the soil within the microcosm. Five replicates per substrate were prepared, resulting in a total of 210 microcosms (42 substrates × 5 replicates). Microcosms were watered by adding water directly over the tube containing faeces/litter so as to reach 70% of soil WHC and incubated at 22 °C and 70% relative humidity in a controlled environment chamber. To limit desiccation while ensuring gas exchange, we drilled four 3-mm holes in each microcosm cap. These microcosms were then weighed weekly and watered to their initial weight at 70% soil WHC. We placed replicates on separated shelves according to a randomised complete block design. Both block positions within the controlled environment chamber and microcosm positions within blocks were randomised weekly. After 180 days, remaining intact litter and faeces in microcosms were collected, dried at 30 °C for 48 h, weighed and ground with a ball mill (TissueLyser II, Qiagen). We measured C and N concentrations in all samples with a flash CHN Elemental Analyser (Flash Smart, ThermoScientific). The percentage of C and N lost after the incubation was calculated as:

    $$frac{{M_{rm{i}} times {rm{CN}}_{rm{i}} – M_{rm{f}} times {rm{CN}}_{rm{f}}}}{{M_{rm{i}} times {rm{CN}}_{rm{i}}}} times 100,$$

    where Mi and Mf are the initial and final 30 °C dry masses, respectively, and CNi and CNf are the initial and final C or N concentrations, respectively.
    Statistics and reproducibility
    To visualise how the 11 physicochemical characteristics were related and how their values differed between all substrates, we used a PCA, with all variables centred and standardised prior to ordination. Because NMR spectra were measured on a composite sample combining the three replicates of each substrate, a unique value was attributed to all replicates for each NMR region.
    To test our first hypothesis, we tested the overall effect of substrate form (faeces vs. intact litter) on quality (scores on PC1 and PC2) and decomposition (C and N losses) of all substrates using Student’s t tests. To identify the faeces types with significantly different quality (scores on PC1 and PC2) and decomposition (C and N losses) compared to that of the intact litter from which the faeces were derived, we tested the effect of substrate identity (all 42 substrates included as individual levels) on quality (scores on PC1 and PC2) and decomposition (C and N losses) using one-way ANOVAs. We then used Tukey’s honestly significant difference tests to determine significant differences between each faeces type and the corresponding intact litter.
    To test our second hypothesis, we expressed the changes in quality and decomposition following litter conversion into detritivore faeces as net differences in quality (scores on PC1 and PC2) and decomposition (C and N losses) between faeces and the litter from which faeces were derived. We then compared the hypothesised role of intact litter quality/decomposition (PC1 and PC2 scores, C and N losses) and the role of detritivore species on changes in quality/decomposition (net differences in PC1 and PC2 scores, C and N losses) by performing ANCOVAs with intact litter quality/decomposition as the continuous variable and detritivore species as categorical variable (all six detritivore species as individual levels). For all ANVOCAs, the variance associated with each term (intact litter quality/decomposition; detritivore species; interaction) was computed by dividing the sum of squares by the total sum of squares.
    To evaluate the relation between quality parameters (PC1 and PC2 scores) and C and N losses from intact litter and faeces separately, we determined the relations between intact litter and faeces C and N losses and their scores on PC1 and PC2 with simple linear regressions and visualised these relations by fitting these variables as supplementary variables on the PCA.
    For all statistical tests on C and N losses, block was included in the model as a random variable. All data were checked for normal distribution and homoscedasticity of residuals. All analyses were performed using the R software (version 3.5.3).
    Reporting summary
    Further information on research design is available in the Nature Research Reporting Summary linked to this article. More

  • in

    An annotated dataset of bioacoustic sensing and features of mosquitoes

    We conducted a laboratory study in the facilities provided by the Natural History Museum of Funchal (Mosquito Lab). Three species of mosquitoes were recorded to determine their dominant frequencies and spectral behaviors. The species used for this collection and study were A. Aegypti, C. Quinquefasciatus & Pipiens and Culiseta, which came from a lab colony established from captures collected in Funchal city in 2019.
    The mosquitoes were kept in an environmental room simulating natural conditions, with 60 ± 10% relative humidity and temperature of 20–25 °C. Mosquitoes were housed individuals in boxes (25 × 25 × 25 cm) covered with a mesh cap. They were fed with 20% sucrose solution supplemented with 1 g aquarium fish food mixed daily from the brand “Sera Guppy Gran”. The duration of the study was approximately 48 days. All mosquitoes used in these experiments were 7–25 days old. For the recording process, sensors were incorporated into the boxes and the tests conducted on 12–18 specimens for Aedes Aegypti, 7–12 specimens for Culex and 4 specimens for Culiseta. The duration of the extracted sequences ranged from 0 to 300 ms. To generate samples closer to real-world acquisition conditions we added environmental noise in some mosquito samples.
    Uncompressed audio of real sound waves was converted to digital format without any further processing. This means that recordings are exact copies of the source audio, recorded in WAV files.
    The acoustic sensor uses a low-noise omnidirectional microphone capsule2. The microphone converts sound into electrical signals with a specific signal to noise ratio (80 dB), self-noise, and residual noise. All these parameters influence the quality of the acquired sound.
    Noise can be a significant problem when acquiring physical signals as voltages. Signal smoothing attempts to capture the essential information in the signal while leaving out the noise. This is done by interpolating the raw signal to estimate the original one17.
    To collect samples, we used three devices: one of them was our prototype comprising a Teensy 3.2 audio board, microphone and environmental sensor for 44.1 kHz sampling rate. The other two were general-purpose smartphones (Huawei P20 Lite and IPhone 4) used to record samples with a 8 and 48 kHz sampling rate, respectively.
    To start a colony for our experience, we installed traps and buckets of water to catch eggs and adult mosquitoes. The female Aedes mosquitoes require a blood meal before each egg-laying18. The eggs are deposited individually on the inner walls of any container capable of storing water. This work was conducted jointly with the Natural History Museum of Madeira and IASaude (the regional health authority of Madeira islands) as part of a plan to control the spread of mosquitoes in the city of Funchal (Fig. 1).
    Fig. 1

    Location and number of traps in the city of Funchal, Madeira, Portugal.

    Full size image

    A. Aegypti mosquitoes, lay the most eggs in the velcro tape, while Culex and Culiseta prefer to lay directly in rafts on still water or in other substances19. Traps with a ventilation system were also used to capture adult mosquitoes, especially Culex and Culiseta.
    Figure 2 shows the procedure from egg collection to mosquito germination, and also the boxes that are used for further acquisition of sound samples. It is noteworthy that after 25–30 days the mosquitoes die due to the conditions imposed in the study.
    Fig. 2

    Procedure for collecting audio samples for different species.

    Full size image

    Step A comprises the gathering of eggs and mosquitoes. The figures show a bucket inside which mosquitoes lay eggs on a velcro tape, and also a trap. These traditional methods allow a fine assessment of the distribution of mosquito populations over time and space (periodically summarized in epidemiological bulletins). In step B the collected eggs are germinated to create a colony. Then, (step C) mosquitoes are placed in boxes and fed with a sugar solution and fish food20. Finally, in step D, audio samples are collected by the devices: mobile phones and low-cost IoT. This procedure is repeated when the colony dies after 25 days, starting from step B.
    Audio was recorded inside boxes (25 × 25 × 25 cm) where the mosquitoes were located at a maximum distance of 27 cm from the microphone placed in the center of the box. The signal amplitude fluctuates significantly over time as the mosquitoes in free flight approach the microphone or move away.
    Continuous recordings were then split into 300 millisecond (ms) snippets. Since mosquitoes have a very short flight, it was necessary to apply a slight stimulus on the wall of the boxes (covered by a net) to force them to fly.
    To analyze each mosquito recording, 34 features were extracted taking into account several parameters of the signal belonging to three different domains: time (1–3), frequency (4–8, 22–34) and cepstrum (9–21), analyzed below in the Technical Validation section21,22.
    These features are often used for speech signal classification, but are useful when handling non-speech signals as well. They enable a comprehensive analysis of the mosquito sounds in terms of amplitude, energy, zero crossing rate, power, frequency variation in the audio file, tonality, loudness, etc. The features are included in the dataset23 and their computation is demonstrated in the Code Availability section. More