Exploring Spatial Patterns and Coexistance

June 11, 2016

(This article was first published on R – biologyforfun, and kindly contributed to R-bloggers)

Today is a rainy day and I had to drop my plans for going out hiking, instead I continued reading “Self-Organization in Complex Ecosystems” from Richard Solé and Jordi Bascompte.

As I will be busy in the coming weeks with spatial models at the iDiv summer school I was closely reading chapter three on spatial self-organization and decided to try and implement in R one of the Coupled Map Lattice Models that is described in the book.

The model is based on a discrete time version of the Lotka-Volterra competition model, for example for two competing species we model the population normalized abundances (the population abundance divided by the carrying capacity) N1 and N2 using the following equations:

N_{1,t+1} = \phi_{N_{1}}(N_{1,t},N_{2,t}) = N_{1,t} * exp[r_{N_{1}} * (1 - N_{1} - \alpha_{12} * N_{2})]

N_{2,t+1} = \phi_{N_{2}}(N_{1,t},N_{2,t}) = N_{2,t} * exp[r_{N_{2}} * (1 - N_{2} - \alpha_{21} * N_{1})]

Where r is the growth rate and the \alpha are the interspecific competition coefficients. Note that since we use normalized abundance the intraspecific competition coefficient is equal to one.

Space can be introduced in this model by considering dispersion between neighboring cells. Say we have a 50 x 50 lattice map and that individuals disperse between neighboring cells with a probability D, as a result the normalized population abundances in cell k at time t+1 will be:

N_{1,t+1,k} = (1 - D) * \phi_{N_{1}}(N_{1,t,k},N_{2,t,k}) + \frac{D}{8} * \sum_{j=1}^{8} \phi_{N_{1}}(N_{1,t,j},N_{2,t,j})

Where j is an index of all 8 neighboring cells around cell k. We have a similar equation for the competing species.

Now let’s implement this in R (you may a nicer version of the code here).

I start by defining two functions: one that will compute the population dynamics based on the first two equations, and the other one that will return for each cell in lattice the row numbers of the neighboring cells:

#Lotka-Voltera discrete time competition function
#@Ns are the population abundance
#@rs are the growth rate
#@as is the competition matrix
  #return population abundance at time t+1

#an helper function returning line numbers of all neighbouring cells
  lnb<-which(dat[,"x"]%in%(vec["x"]+c(-1,-1,-1,0,1,1,1,0)) &
+ dat[,"y"]%in%(vec["y"]+c(-1,0,1,1,1,0,-1,-1)))
  #remove own line

Then I create a lattice, fill it with the two species and some random noise.


The next step is to initialize all the parameters, in the example given in the book the two species have equal growth rate, equal dispersal and equal interspecific competition.

#competition parameters
#growth rate
#dispersal rate
#infos on neighbouring cells
neigh<-apply(dat,1,function(x) find_neighbour(x,dat))

Now we are ready to start the model (the code is rather ugly, sorry for that, saturday coding makes me rather lazy …)

#model the dynamic assuming absorbing boundary (ie ind that go out of the grid are lost to the system)
for(i in generation){  
    #population dynamics within one grid cell
    #grab the neighbouring cells
    #the number of immigrant coming into the focal grid cell
    imm<-(ds/8)*rowSums(apply(neigh_cell,1,function(y) {
    + +imm[2],lnb=x["lnb"])

First let’s reproduce figure 3.8 from the book:

#look at coexistence between the two species within one cell
cell525<-ldply(list_dat,function(x) c(x[525,"N1"],x[525,"N2"]))
+ measure.vars=1:2,variable.name="Species",value.name="Pop")

#look at coexistence across the system
all<-ldply(list_dat,function(x) c(sum(x[,"N1"]),sum(x[,"N2"])))
+ variable.name="Species",value.name="Pop")



And now figure 3.9 plus an extra figure looking at the correlation in population abundance between the two species at the beginning and the end of the simulation:

#compare species-species abundance correlation at the beginning and at the end of the simulation
exmpl$Gen<-rep(c("Gen : 1","Gen : 60"),each=2500)


#look at variation in spatial patterns at the beginning and the end of the simulation
+ measure.vars=c("N1","N2"),variable.name="Species",value.name="Pop")



Neat stuff, I am looking forward to implement other models in this book!


Filed under: Biological Stuff, R and Stat Tagged: ecology, R

To leave a comment for the author, please follow the link and comment on their blog: R – biologyforfun.

R-bloggers.com offers daily e-mail updates about R news and tutorials on topics such as: Data science, Big Data, R jobs, visualization (ggplot2, Boxplots, maps, animation), programming (RStudio, Sweave, LaTeX, SQL, Eclipse, git, hadoop, Web Scraping) statistics (regression, PCA, time series, trading) and more...

If you got this far, why not subscribe for updates from the site? Choose your flavor: e-mail, twitter, RSS, or facebook...

Comments are closed.


Never miss an update!
Subscribe to R-bloggers to receive
e-mails with the latest R posts.
(You will not see this message again.)

Click here to close (This popup will not appear again)