How to implement non-random extinction scenarios

Generation of probabilistic numerical simulations that test how species loss may affect ecosystem functioning can be achieved using the R statistical and programming environment available from . This process requires information on relevant species traits (abundance, A., and biomass, B_i ), an indication of each species contribution to ecosystem functioning, and reference values pertaining to extinction risk for any extinction drivers of interest (Table 10.1). Following the study of

Solan et al. (2004) for marine benthic invertebrates, the latter can be achieved by characterizing the per capita contributions of each species using an appropriate trait-based index e.g. bioturbation potential, BP., that differentially scores each species with respect to their contribution to the relevant ecosystem process (in this case, bioturbation). These scores can then be collated at the population (BP_p = A. x BP) and community level (BP. = 2 BP_p) and related to measured values of ecosystem functioning e.g. mixed depth of sediments, MD; chlorophyll a, Chl a; or organic carbon, Org-C. This relationship can then be used to predict the level of ecosystem functioning (i.e. MD, Chl a or Org-C) for novel communities following a scenario-based manipulation of the community i.e. recalculated levels of BP. . Scenarios of extinction may be random = 1/n, i.e. extinction risk is a reciprocal function of the number of species in the assemblage and, therefore, is not linked to species traits, as in Table 10.1, or will match likely drivers of change, including extinctions ordered by body size (B.), rarity (A.), trophic position, or particular pressures, such as an anthropogenic activity or an element of environmental change. Extinction probabilities for the latter can be estimated, for example, by quantifying the proportional change in species abundance along a gradient of increasing pressure. Once the scenario is assembled, additional biotic conditions can be added. For example, as the functional consequences of extinction are known to depend on the response of surviving species, the simulation can include compensatory responses expressed as changes in population size in the surviving population following an extinction event. This can represent a range of assumed community responses, from no compensation i.e. loss of a species does not affect the population size of any surviving species following, for example, predation or competitive release to full numerical, or biomass, compensation i.e. the best-case outcome that assumes compensation is additive and substitutions of abundance, or biomass, maintain total community density. The following shows how to implement these basic functions in R; for convenience, we supply the implementation code as a Tinn-R file (extinc- tion_code.r, available at http://ukcatalogue.oup. com/product/9780199642267.do).

To call your data (tab delimited text file required), copy and paste the following code into R:

You can check that your data has imported correctly by entering:

Table 10.1 Example dataset for the generation of extinction scenarios detailing species identity (here, benthic invertebrates), species traits (abundance, A,, and biomass, fi), an indication of each species contribution to ecosystem functioning (bioturbation potential of the individual, BP,, and at population level, BP), and reference values pertaining to extinction risk (Extinction Probability, here assumed to be random, 1/n).

A summary of the structure of the data confirming that you have the correct number of observations (x rows = x species) with their relevant traits in columns (y columns = y traits plus an additional column providing any extinction scenario information), followed by a list of your column names will be returned.

The sample dataset we have provided to accompany this chapter (sampledata.csv, available at http:// ukcatalogue.oup.com/product/9780199642267.do) is for a subset of benthic invertebrates and has six columns: Species—the name of the species or taxon; A.—the mean abundance of individuals within

i

each species population; B.—the mean biomass or other appropriate size metric of individuals within each species population; BP.—the bioturbation potential of an individual; BPp—the bioturbation potential of the population; and ExtinctProb—the probability of extinction. In the example dataset, the extinction probabilities are equal indicating that in this scenario the extinctions are random. Alternatively, extinction probabilities in this column could be set to match specific traits, such as A. B. or other information—e.g. vulnerability to a specific extinction driver—such that the most numerous, largest or most susceptible species are more likely to become extinct. To change the risk of extinction to relate to B., for example, use the following code:

Before the simulation can be implemented, it is necessary to generate two extra columns in the data set to record the abundances and extinction probabilities during each simulation. Their values are set to NA to show that no data has yet been generated:

For convenience, a couple of variables are set that store the number of species in the community and the original total biomass of the community—this is needed in case biomass compensation is used:

Two parameters then need to be set to determine how many simulations will be run (for brevity, we use 100 here, but > 1000 simulations is more appropriate for full datasets) and whether or not biomass compensation should be included (TRUE or FALSE):

Finally, an output data frame is set up to record the results of the simulation. This has the following columns to record the respective information:

Simulation—the simulation number,

Nsp—the number of species remaining in the community,

ExtinctSpecies—the species that will go extinct at the next iteration,

Measure—the community function that is being calculated; e.g. Org-C or Chl a,

Value—the value of that community function,

Note that the ExtinctSpecies column contains the species that will go extinct at the next iteration, i.e. the first row records the function of the full community and describes which species will be the next to go extinct. This is to save staggering the ExtinctSpecies column when it comes to the final species going extinct.

The dataset is now established to run the simulations. The simulations are run inside a large for loop which counts down the number of simulations:

At the start of each simulation the abundances and extinction probabilities need to be reset—otherwise they will hold the results left over from any previous simulations:

We now set up another for loop—still inside the main simulation for loop—which counts down the number of species in the community:

Each time around the species loop, we record the value of the community functions, select a species for extinction, apply biomass compensation, if chosen, and then remove the species that has become extinct from the community.

Following Solan et al. (2004), we first calculate BPc (see section 10.3.2 below) using the current abundances and calculated bioturbation potentials of each species:

We then relate measured levels of ecosystem functioning to the value of BPc using an appropriate relationship, such as a linear regression. Here, we use calculated values for Org-C and Chl a—an equivalent equation will relate to your own data:

We then store these values in the output data frame:

Now we select the species to remove based on the extinction probabilities. This is achieved by picking a random number from a uniform distribution, with a minimum of 0 and maximum of 1, and comparing it to the current extinction probabilities. The value of Extinct is the row number of the species that will go extinct:

I f the extinction scenario included biomass compensation, then the biomass that will be lost during each species extinction is calculated. To implement biomass compensation, the abundances of the surviving species are increased by the same proportion as the lost biomass to the total biomass. This ensures that the total biomass is kept constant:

The abundance and the extinction probability of the expired species is set to 0, thereby removing it from the community:

The remaining extinction probabilities need to be normalized so that they sum to 1:

Finally, we record the identification number—the row number in the original data—of the expired species in the output data frame:

We are now at the end of the species loop so it returns to the beginning, removes another species, records the corresponding community function values, and so on, until all species have become extinct. When the final species is reached, the simulation ends and we return to the start of the simulation loop, ready for another simulation with a full community. After all the simulations have been run, the results in the output data frame can be explored.