BASE
====== Crescimento denso-independente com estocasticidade demográfica ======
{{:ecovirt:roteiro:den_ind:dodos.jpg?300 |}}
Os modelos determinísticos de dinâmica populacional não consideram a variação de aptidão entre os indivíduos. Por exemplo, quando usamos o modelo de crescimento discreto
$$N_{t+1} = 1,5 \times N_t$$
supomos que a cada intervalo de tempo o saldo médio entre nascimentos e mortes é de três para dois, causando um crescimento de 50% na população. Isso poderia acontecer se metade dos indivíduos morrer sem deixar descendentes e a outra metade sobreviver e tiver dois filhotes, cada. Também é possível se todos morrerem, mas apenas um ter $1,5 \times N_t$ filhotes antes.
O raciocínio é o mesmo para outros modelos determinísticos. No modelo de crescimento exponencial
$$N(t) = N_0e^{rt}$$
por exemplo, a população cresce por um fator de $e^{rt}$, devido à taxa de crescimento instantânea $r$, que nada mais é que o saldo das taxas de nascimentos e mortes.
Em resumo, as taxas populacionais são médias que resultam de uma infinidade de arranjos de mortes e nascimentos na população, a maioria com variação de aptidão. A simples natureza fracionária das taxas já implica em variação. Por exemplo, uma taxa de natalidade de 0,5 indivíduo/indivíduo.ano indica que alguns se reproduzem e outros não, pois filhotes não vêm em metades!
A //estocasticidade demográfica// é o efeito da variação de aptidão individual sobre a dinâmica populacional. O objetivo desse roteiro é entender a lógica de construção de modelos de dinâmica populacional com esses efeitos, e conhecer suas principais propriedades.
===== Apenas mortes =====
Vamos começar com uma população de $N_0$ indivíduos em que não há nascimentos nem migrações. As mortes ocorrem a uma [[ecovirt:roteiro:den_ind:di_rcmdr#taxa_instantanea_de_crescimento|taxa instantânea]] //per capita// de $\mu = 0,693 \ \text{ano}^{-1}$. O modelo mais simples para calcular os tamanhos dessa população no tempo é o [[ecovirt:roteiro:den_ind:di_rcmdr|exponencial]]:
$$N(t)=N_0e^{(\text{nascimentos}-\text{mortes})t} \ = \ N_0 e^{-0,693t}$$
que para essa taxa de mortalidade prevê que a população é reduzida à metade a cada ano((ou seja, que a meia-vida da população é de um ano. A meia-vida é calculada da mesma forma que o [[ecovirt:roteiro:math:exponencial#tempo_de_duplicacao|tempo de duplicação]].)).
Para isso acontecer, metade do indivíduos deve morrer e a outra sobreviver. Isso quer dizer que a taxa de mortalidade não é a mesma para todos? Para manter essa premissa de homogeneidade (e nosso modelo simples) podemos dizer que a //probabilidade// de morrer é igual para todos. No caso, todos os indivíduos têm 50% de chance de sobreviver ao próximo ano. Se começamos com $N_0=100$, após um ano em média teremos $50$, após dois $25$ e assim por diante, exatamente como prevê o modelo exponencial.
Mas algo importante mudou em nosso modelo: agora o acaso faz o tamanho populacional variar em torno da média, por exemplo:
Vamos supor que temos apenas dois indivíduos. Cada um deles tem 50% de chance sobreviver ao próximo ano. Assumindo que as probabilidades são independentes, há três resultados possíveis((para entender os cálculos a seguir lembre-se das duas regras básicas de probabilidades: 1-a probabilidade de eventos independentes é o produto das probabilidades de cada evento; 2- a probabilidade de eventos alternativos é a soma de suas probabilidades. Para uma explicação completa veja [[https://www.khanacademy.org/math/probability/independent-dependent-probability|aqui]].)):
* Os dois indivíduos morrem, com probabilidade $0,5 \times 0,5 = 0,25$
* Um indivíduos morre e o outro sobrevive, com probabilidade $2 \times 0,5 \times 0,5 = 0,5$((dobramos o produto das probabilidades porque este resultado pode acontecer de duas maneiras: o indivíduo A morre e B sobrevive, **ou** A vive e B morre, cada maneira com probabilidade de $0,25$.))
* Os dois indivíduos sobrevivem, com probabilidade $0,5 \times 0,5 = 0,25$
Isso mostra que em nosso modelo estocástico((o mesmo que probabilístico aqui)) há mais de um valor possível para a população no futuro. Portanto, agora temos incerteza nas projeções, que pode ser muito grande. Em nosso exemplo de uma população com dois indivíduos, a chance de termos o valor esperado é de apenas 50%!
Mas [[http://www.lage.ib.usp.br/rserve/panic.jpg|nada de pânico]]. Nosso exemplo também mostra que a incerteza é mensurável: é possível calcular a probabilidade de cada tamanho populacional ocorrer. No caso da nossa dinâmica apenas com mortes, a probabilidade de um indivíduo sobreviver até o tempo $t$ é:
//**__Probabilidade de Sobreviver__**//
$$p(t)=e^{-\mu t}$$
Assim, esperamos ter $p(t)N_0$ no tempo $t$, ou seja, o tamanho esperado da população((o mesmo que a média das projeções, que chamaremos de $E[N(t)]$)) continua o mesmo do modelo sem estocasticidade:
$$E[N(t)]\ = \ p(t)N_0 \ = \ N_0e^{-\mu t}$$
O que mostra que **em média** o modelo com estocasticidade resulta nas mesmas projeções do modelo determinístico. Mas quanta variação há em torno dessa média? Ou seja, qual a chance de outros valores ocorrerem?
==== Distribuição de probabilidades dos tamanhos populacionais ====
Como calcular a chance de cada tamanho populacional ocorrer? Isso nos leva ao conceito de //distribuição de probabilidades//. Vamos começar com um cálculo simples: a probabilidade de todos sobreviverem até o tempo $t$, em nosso modelo estocástico apenas com mortes. Chamamos essa probabilidade de $P(N(t)\!=\!N_0)$. Como supomos que as probabilidades de morte são independentes entre indivíduos, seu valor é:
//**__Probabilidade de todos sobreviverem__**//
$$P(N(t)\!=\!N_0) \ = \ p(t)^{N_0}$$
Para tamanhos populacionais pequenos essa probabilidade pode ser alta como em nosso exemplo com $N_0=2$ e $p(t\!=\!1)=0,5$:
$$P(N(t)=2)\ = \ 0,5^2 \ = \ 0,25$$
Já quando a população é grande, as chances de que todos sobrevivam são muito pequenas. O mesmo acontece para a probabilidade de todos morrerem, que é
$$P(N(t)=0) \ = \ (1-p(t))^{N_0}$$
É o raciocínio análogo a se perguntar qual a chance de ter apenas caras ou apenas coroas em um certo número de lançamentos de uma moeda. Todos os outros valores entre estes extremos são possíveis, e a cada um deles corresponde uma probabilidade, dada por:
$$P(N(t)\!=\!n) \ = \ \binom{N_0}{n} \ p(t)^n(1-p(t))^{(N_0-n)}$$
Esta é a [[http://en.wikipedia.org/wiki/Binomial_distribution|distribuição binomial de probabilidades]]. Dado um certo número inicial de indivíduos $N_0$ com probabilidades iguais e independentes de morrer após um tempo $t$, essa distribuição dá a probabilidade de $n$ indivíduos sobreviverem. De uma maneira mais geral, a binomial dá a probabilidade de $n$ sucessos em $N_0$ tentativas, dada uma probabilidade constante de sucesso por tentativas.
======A distribuição binomial======
Você verá uma janela com o gráfico com o número de sucessos (no caso sobreviventes) no eixo **//x//**, que vai de zero a $N_0$, e suas probabilidades no eixo **//y//**, de acordo com a distribuição binomial.
Avalie o efeito de mudar os dois parâmetros da binomial: número de tentativas e a probabilidade de sucesso. Experimente alguns valores e proponha regras gerais sobre seus efeitos. Sugestões:
* Mantenha o número de tentativas em 10 e faça a probabilidade de sucesso ir de $0$ a $1$ a intervalos de $0,2$
* Mantenha a probabilidade de sucesso em $0,5$ e aumente o número de tentativas em $2$, $5$, $10$, $100$ , $1000$.
==== Pergunta ====
Para uma população sob dinâmica estocástica apenas de mortes com taxa de mortalidade $\mu=0,693$ e tamanho inicial $N_0=10$:
- Calcule as probabilidades de sobrevivência para $t=1$, $t=2$ e $t=3$
- Faça os gráficos das distribuições de probabilidades dos tamanhos populacionais nestes 3 tempos
====== Simulação no computador ======
Até agora vimos algumas propriedades teóricas da dinâmica populacionais com estocasticidade demográfica:
- Há mais de um tamanho populacional possível a cada tempo;
- Quando há apenas mortes, a probabilidade dos tamanhos populacionais a cada tempo segue uma distribuição binomial;
- A média dos tamanhos populacionais a cada tempo corresponde ao valor previsto pelo modelo sem estocasticidade (determinístico).
Vamos agora testar na prática essas propriedades, e descobrir mais algumas, simulando populações com a dinâmica estocástica de mortes.
====== Parâmetros ======
Os parâmetros controlam as simulações de populações sob estocasticidade demográfica em tempo contínuo:
^ Opção ^ parâmetro ^ O que faz ^
^''Enter name for last simulation data set'' | objeto no R |nome para salvar os resultados da simulação em um objeto no R |
^''Maximum time''| tmax | tempo máximo da simulação na escala de tempo das taxas|
^''Number of simulations''| nsim | número de populações a simular|
^''Initial size''| N0 |tamanho inicial das populações |
^''birth rate''| b |taxa instantânea de nascimentos|
^''death rate''| d |taxa instantânea de mortes|
Vamos simular dez populações de nosso primeiro exemplo, até o tempo 5. Para isso mude as opções de simulação para:
tmax = 5
nsim = 10
N0 = 2
b = 0
d = 0.693
Você deve ver um gráfico parecido com este:
{{:ecovirt:roteiro:den_ind:simpleD-plot.png|}}
As linhas coloridas são as trajetórias das dez populações, e a linha preta a trajetória esperada. No canto superior direito do gráfico está a média e o desvio-padrão do tempo que as populações levaram para diminuir pela metade (meia-vida da população).
==== Meia-vida da população ====
A meia-vida esperada em nossa simulação é de um ano, mas note como algumas populações demoraram bem mais do que outras para cair de dois para um indivíduo, ou para se extinguir. Vamos verificar se ainda assim o tempo médio para a população cair para metade corresponde ao valor teórico.
Testamos isso simulando muitas populações de tamanho inicial $N0=20$ e calculando a média dos tempos que levaram para passar para $N=10$. Ajuste as opções de simulação para
tmax = 3
nsim = 1000
N0 = 20
b = 0
d = 0;693
O gráfico vai ficar bem cheio de linhas, mas o que nos interessa é o valor em //Halving time//. A média está próxima do valor teórico? Agora aumente o tamanho inicial da população para $80$, mantendo os outros parâmetros.
=== Perguntas ===
- Qual foi o efeito do tamanho populacional inicial sobre a média e a variação da meia-vida das populações?
- Qual a explicação para o resultado que você encontrou?
======Distribuição dos tamanhos populacionais======
Veja que os tamanhos das populações ao final da simulações ($t=2$) variam. Sabemos que os valores possíveis vão de $N_0$ (no caso $20$) até zero. A distribuição probabilidade esperada para esses valores é a de uma binomial com $N_0=20$ tentativas e probabilidade de sucesso de $p(t)=e^{-0,693 \times 2}=0,25$.
Faça o gráfico da distribuição com esses parâmetros usando a função da distribuição binomial, como explicamos na seção [[ecovirt:roteiro:den_ind:di_edrcmdrdistribuicao_de_probabilidades_dos_tamanhos_populacionais|sobre distribuições de probabilidades]], acima.
Agora vamos comparar esse gráfico da distribuição teórica com os resultados das simulações. Faça um gráfico da proporção das simulações que terminou com cada tamanho executando os seguintes comandos:
## tamanhos finais de cada populacao
sim1.Nt <- sapply(sim1, function(x) x[sum(x[,1]<2, na.rm=TRUE),2])
## Tabela de frequencia dos tamanhos populacionais
sim1.tab <- table(factor(sim1.Nt, levels=0:20))
## abre uma nova janela grafica
x11()
## Grafico com as proporcoes esperadas de cada tamanho populacional
plot(sim1.tab/1000, xlab="N(t=2)", ylab="Proporcao das populacoes", lwd=5)
Compare os dois gráficos. Há uma boa correspondência? Se quiser sobrepor os valores da binomial no gráfico das simulações copie e execute estes comandos:
## probabilidades esperadas pela binomial
(sim1.esp <- dbinom(0:20, size=20 ,prob=.25))
lines(0:20,sim1.esp, col="blue", type="b")
===== Tamanho médio da população =====
Já guardamos o tamanho das mil populações simuladas em um objeto no R. Para obter a médias desses tamanhos copie e execute o comando
mean(sim1.Nt)
====Pergunta====
Esta média é compatível com o valor esperado?
====== Nascimentos e mortes ======
O que esperar de uma população com nascimentos e mortes estocásticos? O modelo resultante é uma extensão do anterior, com o acréscimo de que a probabilidade de aumento na população não é nula. Vamos usar o **Ecovirtual** para investigar o que muda.
===== Simulação em computador =====
Simule 200 populações com tamanho inicial 1 e com taxa de nascimento duas vezes maiores que a de morte. Execute as simulações com os valores:
# salve o resultado da função no objeto sim2
tmax = 20
nsim = 200
N0 = 1
b = 0.2
d = 0.1
Os tamanhos populacionais agora oscilam em uma [[ecovirt:roteiro:math:bebadorcmdr|caminhada aleatória]], devido à sucessão de nascimentos e mortes. Como a taxa de nascimentos é o dobro da de mortes, um nascimento é duas vezes mais provável do que uma morte, a qualquer instante. O senso comum diria então que as populações não correm risco de extinção. Será?
===== Tamanho populacional médio =====
De forma similar ao modelo anterior, o tamanho populacional esperado é
$$E[N(t)]=N_0e^{rt}$$
onde $r$ é a taxa instantânea de crescimento, que é a diferença entre a taxa de nascimentos e de mortes.
Novamente guardamos os resultados em um objeto do R, do qual podemos calcular os tamanhos populacionais médios. Para isso copie e execute os comandos abaixo
sim2.Nt <- sapply(sim2, function(x)x[nrow(x),2])
mean(sim2.Nt)
=== Pergunta ===
Verifique se a média dos tamanhos populacionais das simulações é compatível com a média teórica.
===== Distribuição dos tamanhos populacionais =====
Mas já sabemos que a média não conta toda a estória. Como em qualquer modelo com estocasticidade, agora não temos mais um único valor possível para o tamanho populacional a cada tempo, e sim um conjunto de valores possíveis e suas respectivas probabilidades. Faça um histograma dos tamanhos populacionais no tempo final da simulação com o código:
sim2.tab <- table(factor(sim2.Nt,
levels=0:max(sim2.Nt)))/length(sim2.Nt)
plot(sim2.tab, xlab="N(tmax)", ylab="Proporcao das populacoes", lwd=5)
No processo de mortes estocásticas vimos que as probabilidades de cada tamanho populacional a cada tempo seguem a distribuição binomial. Para um processo de nascimentos e mortes, as probabilidades seguem uma outra distribuição, chamada [[http://en.wikipedia.org/wiki/Negative_binomial_distribution|binomial negativa]].
O que mais nos interessa aqui é que há uma probabilidade de um tamanho populacional zero, ou seja, de extinção das populações. Em nossas simulações, este $N(t)=0$ foi o valor mais frequente, o que é fácil de entender: como o tamanho inicial é de um indivíduo, há uma chance de um em dois de uma morte antes do primeiro nascimento, o que extingue a população.
===Perguntas===
* 1. Para uma dinâmica estocástica com taxa de nascimentos maior que de mortes, qual é o efeito sobre a probabilidade de extinção de:
* Tempo de simulação?
* Tamanho inicial da população?
* Razão entre as taxas?
**__Utilize o EcoVirtual__**
- anote a cada simulação a proporção de populações extintas;
- faça um gráfico da probabilidade de extinção e as variáveis (tempo, tamanho e razão entre taxas);
*2. Compare suas conclusões com as obtidas se nascimentos e mortes se equivalem, como no roteiro de [[ecovirt:roteiro:math:bebadorcmdr|caminhada aleatória]].
====== Para saber mais ======
* Renshaw, E. (1991). Modelling biological populations in space and time Cambridge University Press. //Este roteiro segue a lógica do segundo capítulo desse livro, que é uma ótima introdução a modelos estocásticos de nascimentos e mortes//.
* Akçakaya H.R., Burgman M.A & Ginzburg, L.V. (1999). [[http://www.ramas.com/apppopn.htm|Applied population ecology - Principles and computer exercises using RAMAS EcoLab]]. //Outro livro muito didático, com exercícios em computador com o programa proprietário [[http://www.ramas.com/ecolab.htm|RAMAS ecolab]]. O capítulo 2 é uma excelente introdução às fontes de estocasticidade em dinâmicas populacionais.//
* [[http://cmq.esalq.usp.br/BIE5781/doku.php?id=01-discretas:01-discretas|Tutoriais em R]] sobre distribuições de probablilidades discretas. Da disciplina [[http://cmq.esalq.usp.br/BIE5781|Modelagem Estatística para Ecologia e Recursos Naturais]], dos Programas de Pós-Graduação em Ecologia e em Recursos Florestais da USP. Incluem tutoriais sobre a distribuição binomial e binomial negativa.