Vamos investigar nesse exercício modelos de sucessão e incorporar neles distúrbio. No exercício de coexistência simulamos distúrbio através da taxa de extinção. Vamos partir desse modelo e avaliar o que acontece com a coexistência conforme aumento da intensidade do distúrbio. Vamos relembrar o modelo de coexistência e o exercício que fizemos. Tínhamos duas espécie, sendo a primeira melhor competidora pois podia colonizar manchas já ocupadas pela outra, enquanto a segunda só colonizava manchas vazias. A variação na ocupação de manchas era dada por:

$$ (df_1)/dt = i_1f_1(1-f_1)-p_e f_1 $$

$$ (df_2)/dt = i_2f_2(1-f_1-f_2)-i_1f_1f_2-p_e f_2 $$

onde:

• f = fração de machas ocupadas
• pe = probabilidade de extinção por mancha
• i = taxa de incremento da probabilidade de colonização com o aumento de f
• A taxa de colonização, portanto, é o produto i.f, e varia com a fração de manchas ocupadas (quanto mais ocupação, mais propágulos).

Para que a espécie dois persistisse no sistema era necessário que:

$$ e/i_(1) > i_(1)/i_(2) $$

Abaixo o código da função utilizada na aula de coexistência, copie e cole no R para podermos utilizá-la.

meta.comp<-function(tf,l,c,fi1,fi2,i1,i2,pe,plot.eq=FALSE,D=0){
  F1 <- 1-(pe/i1)
  F2 <- pe/i1-i1/i2
  if(F1<=0) F2 <- 1-(pe/i2)
  Nt <- l*c
  N <- floor(Nt*(1-D))
  resultado=matrix(nrow=tf,ncol=3)
  n1 <- floor(fi1*N)
  n2 <- floor(fi2*N)
  antes <- sample(rep(c(2,1,0),c(n2,n1,N-(n1+n2))))
  resultado[,1] <- 1:tf
  resultado[1,2:3] <- c(sum(antes==1),sum(antes==2))/N
  for(t in 2:tf){
    depois <- rep(0,N) 
    pi1=i1*sum(antes==1)/Nt
    pi2=i2*sum(antes==2)/Nt
    depois[antes==1]<-sample(c(0,1),sum(antes==1),replace=T,prob=c(pe,1-pe))
    depois[antes==2]<-sample(c(0,2),sum(antes==2),replace=T,prob=c(pe,1-pe))
    depois[antes==0] <- sample(c(0,2),sum(antes==0),replace=T,prob=c(1-pi2,pi2))
    d1<-sample(c(0,1),sum(antes!=1),replace=T,prob=c(1-pi1,pi1))
    depois[antes!=1][d1==1] <- 1
    resultado[t,2:3]=c(sum(depois==1),sum(depois==2))/Nt
    antes <- depois
  }
  plot(1:tf,resultado[,2],type="l",xlab="Tempo",ylab="Fração de manchas ocupadas",
       ylim=c(0,1),main=paste("Competição com Colonização Interna","\n c=",c," l=",l," fi1=",fi1," fi2=",
                     fi2,"\n i1=",i1," i2=",i2," pe=",pe," D=",D),font.lab=2,lwd=2, col="red")
  lines(1:tf,resultado[,3],col="blue", lwd=2)
  if(plot.eq==T){
    abline(h=F1,col="red",lwd=1.5,lty=2)
    abline(h=F2,col="blue",lwd=1.5,lty=2)
  }
  if(D>0)abline(h=1-D,lty=2)
  legend("topright", c("Melhor competidora", "Pior competidora"),col=c("red","blue"),lty=1)
  invisible(resultado)
}

Os argumentos da função são:

• o número de linhas e colunas da matriz (l, c),
• o número de intervalos de tempo para reiterar a simulação (tf),
• a fração inicial de manchas ocupadas por cada espécie (fi1, fi2) e
• os parâmetros do modelo (i1, i2, pe ) ¹⁾.

Vamos simular um aumento crescente da probabilidade de extinção:

par(mfrow=c(2,3)) ## dividindo a tela gráfica em seis partes 
meta.comp(tf=100,c=100,l=100,fi1=0.1,fi2=0.1,i1=0.4,i2=0.8,pe=0.10, plot.eq=FALSE)
meta.comp(tf=100,c=100,l=100,fi1=0.1,fi2=0.1,i1=0.4,i2=0.8,pe=0.20, plot.eq=FALSE)
meta.comp(tf=100,c=100,l=100,fi1=0.1,fi2=0.1,i1=0.4,i2=0.8,pe=0.25, plot.eq=FALSE)
meta.comp(tf=100,c=100,l=100,fi1=0.1,fi2=0.1,i1=0.4,i2=0.8,pe=0.30, plot.eq=FALSE)
meta.comp(tf=100,c=100,l=100,fi1=0.1,fi2=0.1,i1=0.4,i2=0.8,pe=0.40, plot.eq=FALSE)
meta.comp(tf=100,c=100,l=100,fi1=0.1,fi2=0.1,i1=0.4,i2=0.8,pe=0.50, plot.eq=FALSE)

1. O que está acontecendo com o sistema conforme aumentamos a intensidade do distúrbio (pe)?
2. Em que cenário o sistema é mais diverso?
3. Interprete o comportamento da espécie pior competidora (azul) no início das quatro primeiras simulações. O que está acontecendo?

Para criar um modelo onde exista demanda conflitante (tradeoff) entre colonização e competição em um sistema multiespécies, precisamos primeiro criar uma cenário onde há um variação na habilidade de colonização/competição das espécies. Para tanto vamos criar um sistema onde a habilidade competitiva é inversamente relacionada a habilidade de colonização, uma clássica demanda conflitante. No equilíbrio ou em intensidades baixas de distúrbios, nossa comunidade virtual teria a espécie melhor competidora com maior proporção de manchas ocupadas, assim como no modelo de coexistência de duas espécies. Vamos estabelecer que essa proporção seja 20% de manchas ocupadas e em seguida estabelecer uma sequência ordenada hierárquica onde a próxima melhor competidora ocupará 20% das manchas restantes e a seguinte 20% do que restou depois, até nossa última espécie (pior competidora) ocupar uma fração de 20% daquelas últimas manchas não ocupadas por nenhuma das outras espécies. Nesse caso, teríamos o seguinte cenário na distribuição de abundâncias das nossas espécies:

$$p_(i)=d(1-d)^(i-1) $$

onde: $$p_(i)$$ = proporção de manchas ocupadas pela espécie cuja ordem na classificação de abundância é i d= proporção de manchas ocupadas pela espécie mais abundante

Vamos fazer um gráfico desse cenário de abundâncias para 10 espécies sendo a proporção de manchas ocupadas pela melhor competidora 20%:

x11() # abre uma nova janela gráfica
s=10 # número de espécies no sistema no início da simulaçao
i=1:s # posição da espécie no classificação de competição
i
d=0.2 # potencial de manchas ocupadas pela melhor competidora
p=d*(1-d)^(i-1) # proporção de manchas ocupadas por cada espécie no equilíbrio
p
plot(i, p, type="b", ylab="Proporção de manchas ocupadas", xlab="Ordem das espécies")

O Pesquisador David Tilman demonstrou que esse cenário de distribuição de abundâncias é possível quando todas as espécies experienciam a mesma taxa de extinção/mortalidade e quando a taxa de colonização de cada espécie é dada por:

$$ c_(i)= e/(1-d)^(2i-1)$$

Vamos adicionar os valores de taxas de colonização no nosso gráfico anterior, calculados para uma taxa de extinção igual a 0.04 para todas as espécies:

pe=0.04
c= pe/(1-d)^(2*i-1)
c
par(mfrow=c(1,2))
plot(i, p, type="b", ylab="Abundância Proporcional", xlab="Ordem de classificação da espécie")
plot(i, c, type="b", col="red", ylab="Taxa de Colonização", xlab="Ordem de classificação da espécie" ) 

Essa é nossa situação de equilíbrio na comunidade, uma distribuição geométrica de abundâncias com a mais abundante chegando a 20% e com um taxa de extinção de 4% para todas as espécies. Além disso, o mesmo pesquisador generalizou a equação de variação da proporção de manchas ocupadas ao longo do tempo para n espécies, como sendo:

$$ (dp_(i))/dt = c_(i)p_(i)(1-sum_(j=1)^i(p_i)) - m_ip_i - sum_(j=1)^(i-1)(c_jp_jp_i)$$

Vamos montar essa comunidade simulando sua dinâmica ao longo do tempo. A nossa função tem os seguintes passos:

1. cria a ordenação das espécies
2. calcula a coeficiente de colonização como acima
3. se o fi (proporção de manchas ocupadas no início) é um valor único, preenche essa proporção de manchas com (s) espécies ao acaso
4. se fi é um vetor de s valores, esses são utilizados para o preenchimento inicial das manchas sendo as proporções no vetor fi relacionadas à ordem das espécies, começando com a espécie melhor competidora para a pior
5. inicia ciclos até tf
6. a cada ciclo calcula a extinção para cada espécie e sua colonização respeitando a hierarquia de competição para decidir qual espécie é deslocada da mancha e qual permanece
7. guarda a proporção de manchas ocupadas por cada espécie a cada ciclo
8. produz o gráfico e retorna o resultado de número de manchas ocupadas por cada espécie a cada tempo

Copie e cole o código abaixo no R para rodar a função com.compete

com.compete = function(tf,l,c, s, fi, d, pe,fr=0,int=0)
{
rank=1:s
vetor.dist=rep("n", tf)
  if(fr>0 & int>0)
  {
  n.dist=floor(tf/fr)
  max.dist=n.dist*fr
  vetor.dist[seq(0,max.dist, by=fr)[-1]]="d"
  }
vetor.dist=vetor.dist[-1]
ci= pe/(1-d)^(2*rank-1)
N <- l*c
resulta=matrix(nrow=s,ncol=tf)
attributes(resulta)=c(attributes(resulta),list(tempo=tf, riqueza=s, n.manchas=N, inicial= fi, comp=d, coloniza=ci, freq.dist=fr, int=int))
n.ocup= sum(fi)*N
t=1
	if(length(fi)==s)
	{ 
	resulta[,1]=fi
	antes=sample(c(rep(0,N-n.ocup),rep(rank,each=fi*N)))
	}
	if(length(fi)==1)
	{
	antes <- sample(c(1:s, sample(1:s, (n.ocup-s),replace=TRUE), rep(0, N-n.ocup)))
	t.antes=table(antes)
	n.sp=rep(0,s)
	names.antes=match(as.numeric(names(t.antes)), rank)
	rank.antes=match(rank,as.numeric(names(t.antes)))
	n.sp[na.omit(names.antes)]=t.antes[na.omit(rank.antes)]
	resulta[,1]<-n.sp/N
	}
	for (f in vetor.dist)
		{
		t=t+1
		depois <- rep(0,N)
		ncum=cumsum(resulta[,1])-resulta[1,1]
		pi=ci*resulta[,t-1]
		pi[pi>1]=0.999
			for(r in s:1)
			{
			depois[antes==r]<-sample(c(0,r),sum(antes==r),replace=T,prob=c(pe,1-pe))
			d1<-sample(c(0,r),sum(antes>r | antes==0),replace=T,prob=c(1-pi[r],pi[r]))
			depois[antes>r | antes==0][d1==r] <- r
			}
		n.sp=rep(0,s)
		if(f=="d")
		{
		depois[sample(1:N,N*int)]=0
		}
		t.depois=table(depois)
		names.match=match(as.numeric(names(t.depois)), rank)
		rank.match=match(rank,as.numeric(names(t.depois)))
		n.sp[na.omit(names.match)]=t.depois[na.omit(rank.match)]
		resulta[,t]<-n.sp/N
		antes<-depois
		}
matplot(1:tf,t(resulta),type="l",lty=1:s, col=rainbow(s),ylim=c(0,max(resulta)+0.2), xlab="Tempo", ylab="Fração de manchas ocupadas",main=paste("Competição x Colonização \n d=",d," pe=",pe, "fr= ",fr, "int= ", int))
legend("topright", legend=paste("sp", 1:s, sep="_"), lty=1:s, col=rainbow(s), bty="n", cex=0.5)
invisible(resulta)
}

Nossa função tem os seguintes argumentos:

1. tf= tempo final da simulação,
2. l= número de linhas no espaço de manchas,
3. c= número de colunas no espaço de manchas,
4. s= número de espécies no tempo inicial,
5. fi= proporção de manchas ocupadas no tempo inicial,
6. d= proporção de manchas ocupadas pela espécie melhor competidora no equilibrio,
7. pe= taxa de extinção das manchas,
8. fr= frequencia de distúrbio (usada somento mais a frente no exercício),
9. int= intensidade do distúrbio (idem)

Vamos testar o nosso modelo com 10 espécies (s=10), 10000 manchas (l = 100, c= 100 ), com a espécie competidora tendo o potencial de ocupar 20% das manchas, probabilidade de extinção pe=0.04. Vamos rodar primeiro com 200 ciclos de tempo e depois com 1000. Em ambas as simulações vamos começar com uma proporção de manchas ocupadas no início baixa (fi=0.2), simulando uma situação onde o sistema está fora do equilíbrio.

teste1=com.compete(tf=200,l=100,c=100, s=10, fi=0.2, d=0.2, pe=0.04)
x11()
teste2=com.compete(tf=1000,l=100,c=100, s=10, fi=0.2, d=0.2, pe=0.04)

• Interprete o resultado do modelo em termos de:
  1. demanda conflitante de colonização x competição,
  2. substituição de espécies e
  3. sucessão ecológica.

## incluir aqui mais alguns exemplos variando parâmetros

Vamos simular agora um sistema com distúrbios. Nos modelos de metapopulação vimos que a probabilidade de extinção ($$p_e$$) podia ser interpretada como distúrbio, já que atingia uma parte das manchas tornando-as vagas para posterior colonização. Aqui, vamos interpretar o $$p_e$$ como a mortalidade basal das populações, não mais como distúrbio. Essa probabilidade de morte é constante e igual para todas as espécies no nosso modelo, somente para simplificar a simulação. Vamos criar distúrbios na nossa comunidade que acontecem com uma certa frequência de intervalos regulares (fr) e uma certa intensidade (int) relacionada à proporção de manchas afetadas.

Antes de incluir o distúrbio no sistema, vamos fazer um gráfico do sistema de referência sem distúrbio. No exemplo 2 o resultado foi guardado no objeto teste2. Vamos usá-lo para não ter que rodar novamente a função ²⁾:

par(mfrow=c(2,2))
matplot(1:attributes(teste2)$tempo,t(teste2),type="l",ylim=c(0,max(teste2)+0.1),lty=1:dim(teste2)[1], col=rainbow(dim(teste2)[1], start=.7, end=.4), xlab="Tempo", ylab="Fração de manchas ocupadas",main=paste("Competição x Colonização" ))
legend("topright", legend=paste("sp", 1:10, sep="_"), lty=1:10, col=rainbow(10, start=.7, end=.4), bty="n", cex=0.5)

## testando distúrbios com frequência baixa (intervalo grande entre eventos fr=100) e intensidade baixa (0.2)
teste3=com.compete(tf=1000,l=100,c=100, s=10, fi=0.2, d=0.2, pe=0.04, fr=100, int=0.2)

## distúrbio com frequência baixa (intervalo grande fr=100) e intensidade média (0.5)
teste4=com.compete(tf=1000,l=100,c=100, s=10, fi=0.2, d=0.2, pe=0.04, fr=100, int=0.5)

## distúrbio com frequência alta (intervalo pequeno entre eventos, 10) e intensidade média (0.2)
teste5=com.compete(tf=1000,l=100,c=100, s=10, fi=0.2, d=0.2, pe=0.04, fr=10, int=0.2)

• Responda as seguintes questões:
  1. Em que situações as melhores competidoras se dão melhor?
  2. E as piores competidoras?
  3. Consegue detectar algum padrão na coexistência de espécies com essas poucas simulações?

É creditado ao pesquisador Joseph Connell (1978) a teoria do distúrbio intermediário, apesar dela ter sido proposta por J. Phillip Grime alguns anos antes em 1973, tratando de exclusão competitiva em plantas. O trabalho de Connel (veja referência) pode ter tido preferência na citação da teoria por ter contrastando dois ambientes muito distintos e reconhecidamente entre os mais diversos do planeta (recifes de corais e florestas tropicais). Nele o pesquisador advoga que o principal fator relacionado à manutenção da alta diversidade nesse ecossistemas é a presença de distúrbios em frequências e intensidade intermediárias. Faça simulações para testar essa teoria.

<box red 70% | Dicas > Salve cada simulação em um objeto de nome diferente. Nos exemplos acima fizemos isso denominando cada simulação por um nome (teste1, teste2…) à direta do sinal de igualdade ou com o símbolo (←), veja abaixo:

teste5 <- com.compete(tf=1000,l=100,c=100, s=10, fi=0.2, d=0.2, pe=0.04, fr=10, int=0.2)

Desta forma o resultado da simulação é salvo em objetos diferentes e podemos depois usar esses resultados. Por exemplo, podemos estar interessados em ver quais as espécies mais abundantes no final da simulação, por exemplo nos 10 últimos ciclos:

tf5=dim(teste5)[2]
teste5[,(tf5-10):tf5]

Podemos contar o número de espécies ao longo do tempo, da seguinte forma:

teste2a=teste2
teste2a[teste2>0]=1
apply(teste2, 2, sum)

Caso queire fazer alguma outra operação com os dados resultantes das simulações, contate um monitor ou professor. </box>

• Connell, J. H. 1978 Diversity in tropical rain forests and coral reefs. Science 199(4335):1302 - 1310
• Hastings, A. (1980) Disturbance, coexistence, history and competition for space. Theoretical Population Biology, 18:363–373.
• Stevens, M.H.H. (2009) A primer in ecology with R. New York, Springer. capítulo 9