4 Resultados e Discussão

Neste capítulo apresentaremos resultados numéricos para a metodologia estudada. Dividimos o capítulo em duas partes: na primeira parte, ilustraremos a metodologia em duas aplicações com dados reais de captura-recaptura; e na segunda parte realizaremos um estudo de simulação Monte Carlo para verificar o comportamento dos estimadores de máxima verossimilhança – com foco no estimador do parâmetro N – dos dois modelos apresentados, Mt e Mtb.

4.1 Aplicação com dados reais

Nesta seção vamos aplicar os dois modelos de captura-recaptura estudados, Mt e Mtb, em dois conjuntos de dados reais da literatura. Para fins comparativos, iremos escolher como melhor ajuste aos dados o modelo que apresentar os menores valores dos critérios de informação AIC (Critério de Informação de Akaike) e BIC (Critério de Informação Bayesiano), definidos respectivamente por


AIC=2logL(ˆθ|x)2dBIC=2logL(ˆθ|x)dlog(n)


onde ˆθ é a estimativa de máxima verossimilhança do vetor paramétrico θ de dimensão d, e x é o conjunto de dados contendo n observações . Em nossos modelos, n é o número total de distintos capturados, denotado neste projeto por r.

4.1.1 Primeira aplicação

Como parte de um estudo de ciências populacionais, Coulombe, Ridgway, e Evans (1965) conduziram um estudo com armadilhas para ratos domésticos selvagens Mus musculus. Um quadrado 10×10 foi usado com armadilhas espaçadas de 3 metros entre elas e as armadilhas foram verificadas pelos pesquisadores duas vezes ao dia, de manha e a noite, durante 5 dias, totalizando k=10 épocas de captura. Um total de r=173 ratos diferentes foram capturados. Os dados estão disponíveis em Otis et al. (1978) e são exibidos na Tabela 4.1.


Tabela: 4.1 Dados reais de captura-recaptura de ratos Otis et al. (1978) .
j 1 2 3 4 5 6 7 8 9 10
Animais marcados antes da ocasião Mj 0 68 102 128 140 156 159 171 171 173
Animais capturados (não marcados) uj 68 34 26 12 16 3 12 0 2 0
Animais recapturados (marcados) mj 0 27 36 40 58 38 64 36 74 39
Animais selecionados nj 68 61 62 52 74 41 76 36 76 39

Observe a baixa frequência de animais não marcados (uj) nas últimas ocasiões de captura, indicando que poucos animais da população ainda não foram vistos no estudo. Os tamanhos amostrais (nj) diferentes também indicam heterogeneidade ao longo do tempo.

A Tabela 4.2 a seguir exibe as estimativas dos dois modelos, Mt e Mtb.


Tabela: 4.2 Estimativas de máxima verossimilhança dos parâmetros dos modelos Mt e Mtb na Aplicação 1.
Estimativas Mt Mtb
ˆN 175 176
ˆc 0.0858
ˆp1 0.3886 0.3864
ˆp2 0.3486 0.3391
ˆp3 0.3543 0.3410
ˆp4 0.2971 0.2826
ˆp5 0.4229 0.4039
ˆp6 0.2343 0.2196
ˆp7 0.4343 0.4129
ˆp8 0.2057 0.1913
ˆp9 0.4343 0.4115
ˆp10 0.2229 0.2074


Primeiramente, destacamos que as estimativas para o tamanho populacional dos dois modelos foram bem próximas, 175 animais para Mt e 176 para Mtb, e não muito superiores ao número total de animais distintos capturados, r=173 animais. Isto é, estima-se que apenas dois ou três animais não foram vistos durante o estudo de acordo com os modelos Mt e Mtb, respectivamente. A heterogeneidade entre as probabilidades de captura é observada nas estimativas dos parâmetros pj, que variaram entre 0,2057 e 0,4343 para o modelo Mt e de 0,1913 a 0,4129 no modelo Mtb. O parâmetro de efeito comportamental à marcação foi estimado em ˆc=0,0858, bem próximo de zero que indica não efeito da marcação na chance de recapturar o animal. Essa evidência foi confirmada pelos critérios AIC e BIC exibidos na Tabela 4.3 onde escolhemos o modelo Mt como sendo de melhor ajuste por apresentar menores valores dos critérios.


Tabela: 4.3 Critérios AIC e BIC para os modelos Mt e Mtb na Aplicação 1.
Critério Mt Mtb
AIC 112.74 114.55
BIC 147.42 152.39


Para fins ilustrativos, reajustamos os dados considerando somente as k=5 primeiras épocas de captura, o que totaliza r=156 animais capturados. A Tabela 4.4 exibe as novas estimativas dos dois modelos, Mt e Mtb, de acordo com esses dados restritos. Observe que a estimativa de N para o modelo Mt se manteve consistente com ˆN=173, enquanto que o modeo Mtb elevou a estimativa para ˆN=184. Os critérios AIC e BIC neste caso foram de 57.20 e 75.50 para o modelo Mt e 58.55 e 79.90 para o modelo Mtb, ainda favorecendo o modelo Mt para ajuste.


Tabela: 4.4 Estimativas de máxima verossimilhança dos parâmetros dos modelos Mt e Mtb na Aplicação 1 (k=5).
Estimativas Mt Mtb
N 173 184
c 0.3298
p1 0.3886 0.3864
p2 0.3526 0.3042
p3 0.3584 0.2965
p4 0.3006 0.2376
p5 0.4277 0.3431


4.1.2 Segunda aplicação

Nesta seção vamos utilizar dados reais de captura-recaptura de ratos veados que foram coletados por S. Hoffman em meados de julho de 1974, inicialmente utilizado por Otis et al. (1978) e discutidos no contexto bayesiano por Wang et al. (2015). Os dados foram coletados da seguinte forma: a captura foi em 5 manhãs consecutivas e 110 ratos distintos foram capturados, onde presumiu-se que exista resposta comportamental à marcação. Os resultados são apresentados na Tabela 4.5.


Tabela: 4.5 Dados reais de captura-recaptura de ratos veados Wang et al. (2015).
j 1 2 3 4 5
Animais marcados antes da ocasião Mj 0 37 68 77 98
Animais capturados (não marcados) uj 37 31 9 21 12
Animais recapturados (marcados) mj 0 23 49 44 57
Animais selecionados nj 37 54 58 65 69


A 4.6 exibe as estimativas de máxima verossimilhança dos parâmetros, onde nota-se a diferença significativa entre as estimativas para ˆN=113 para o modelo Mt e de ˆN=215 para o modelo Mtb. De fato, o efeito comportamental à marcação nesta aplicação foi de ˆc=2,5195, considerada alta, o que pode estar influenciando as estimativas dos demais parâmetros entre os dois ajustes.


Tabela: 4.6 Estimativas de máxima verossimilhança dos parâmetros dos modelos Mt e Mtb na Aplicação 2.
Estimativas Mt Mtb
N 113 215
c 2.5196
p1 0.3274 0.1721
p2 0.4779 0.1579
p3 0.5133 0.1121
p4 0.5752 0.1201
p5 0.6106 0.1012


De fato, temos que a probabilidade de capturar um animal não marcado (sendo a mesma para o animal marcado no modelo Mt) na j-ésima ocasião é: pj=eηj1+eηj,   j=1,2,,k e que a probabilidade de capturar um animal marcado pelo modelo Mtb na j-ésima ocasião é: cj=ec+ηj1+ec+ηj,   j=2,3,,k.

Logo, pela propriedade da invariância dos estimadores de máxima verossimilhança, temos pelo modelo Mtb que as probabilidades estimadas de recapturas são ˆc2= 0.6996, ˆc3= 0.6106, ˆc4= 0.6290 e ˆc5= 0.5831, muito superiores que os ˆpj exibidos na tabela anterior. De fato, temos que a razão de chances está de ˆw=exp(ˆc)=exp(2.5196)=12,4236 de recaptura para a probabilidade de captura, evidenciando um forte efeito (feliz’’) no comportamento do animal devido à marcação. Isso é consistente com os dados, uma vez que há mais recapturas do que as capturas iniciais na maioria das ocasiões de amostragem.

O efeito significativo no parâmetro comportamental no ajuste final é descatado também pelos critérios AIC e BIC apresentados na Tabela 4.7, onde observados menores valores de ambos critérios para o modelo Mtb e, assim, sendo escolhido como melhor ajuste aos dados da Tabela 4.5.


Tabela: 4.7 Estimativas de máxima verossimilhança dos parâmetros dos modelos Mt e Mtb na Aplicação 2.
Critério Mt Mtb
AIC 89.10 67.37
BIC 105.31 86.28


4.2 Estudo de simulação

Para estudar o comportamento dos estimadores de máxima verossimilhança, utilizaremos as medidas de desempenho: Viés, Viés Relativo e Erro Quadrático Médio das estimativas obtidas dos parâmetros para um número M de amostras de captura-recaptura simuladas.

Genericamente, suponha que temos um parâmetro θ e para cada uma das M amostras geradas, temos uma estimativa de máxima verossimilhança ˆθm para θ, m=1,2,,M. Suponha que o verdadeiro valor de θ usado na geração das amostras seja θverd. Assim, definimos o Viés, Viés Relativo e Erro Quadrático Médio, respectivamente, por

Viés=1MMm=1ˆθmθverd.ViésR=1MMm=1ˆθmθverd.θverd.EQM=1MMm=1(ˆθmθverd.)2

4.2.1 Modelo Mt

Para o estudo do comportamento do estimador de máxima verossimilhança do parâmetro N no modelo Mt, fixamos o verdadeiro tamanho populacional nas simulações em N=100,300 e 500 e consideramos números diferentes de épocas de capturas, sendo k=5,6,,15. Para os parâmetros pj, consideramos quatro cenários diferentes na geração destes valores, ilustrados na Tabela 4.8. A ideia é avaliar as estimativas variando cenários onde a probabilidade de captura é baixa (P1) até uma probabilidade de captura relativamente alta (P4).


Tabela: 4.8 Cenários para geração das probabilidades de captura pj.
Cenário Distribuição
P1 p1,...,pkUniforme(0,05;0,15)
P2 p2,...,pkUniforme(0,10;0,20)
P3 p3,...,pkUniforme(0,15;0,25)
P4 p3,...,pkUniforme(0,20;0,30)


Observe que, com a combinação dos diferentes tamanhos populacionais N, épocas de capturas k e probabilidades de captura pj’s, temos um total de 3×11×4=132 diferentes cenários para avaliarmos no estudo de simulação. Para cada um dos 132 cenários, foram geradas 10.000 amostras de captura-recaptura.

A função implementada em código R para geração dos dados é exibida a seguir. Observe que a função é genérica para ambos os modelos (Mt e Mtb) sendo que, quando estivermos gerando dados do modelo Mt, basta fixar c=0 durante a geração dos dados.

# --- Funcao para gerar as amostras : 

gera_dados <- function(N,K,pj,c,seed){
  set.seed(seed)
  u = m = n = M = c()
  etaj = log(pj/(1-pj))
  cj   = exp(etaj+c)/(1+exp(etaj+c))  # probabilidade de recaptura
  u[1] = rbinom(1,N,pj[1])            # numero de animais capturados 
  m[1] = 0                            # numero de animais marcados
  n[1] = u[1]+m[1]                    # numero de animais selecionados
  M[1] = 0                            # numero de animais marcados na populacao
  # segunda ocasiao em diante:
  for(i in 2:K){
    M[i] = M[i-1] + u[i-1]
    u[i] = rbinom(1,N-M[i],pj[i])
    m[i] = rbinom(1,M[i],cj[i])
    n[i] = u[i] + m[i]
  }
  r = M[K] + u[K]
  return(list(u,m,M,n,r))
}

Na Tabela 4.9 temos a relação em porcentagem do número médio de animais distintos observados dividido pelo tamanho da população em cada um dos cenários. Por exemplo, no cenário com k=5, N=100 e Cenário P1, temos que, em média, o número de animais distintos observados nas 10.000 amostras está em torno de 41,4% do total da população, isto é, 41,4 animais em N=100. Já para o caso de k=15, N=500 e Cenário P4, captura-se em torno de 98,7% dos animais (média de 493,5 animais), isto é, quase todos eles são vistos durante o estudo. Em geral, observamos que com o crescimento de k ou cenários onde a probabilidade de captura é maior, tem-se em média mais animais capturados, como esperado.


Tabela: 4.9 Relação em porcentagem do número médio de animais distintos observados dividido pelo tamanho da população no primeiro estudo de.
Cenário P1 Cenário P2 Cenário P3 Cenário P4
k N=100 300 500 100 300 500 100 300 500 100 300 500
5 41.4 40.5 40.9 55.8 55.6 55.6 67.3 67.3 67.3 76.3 76.3 76.3
6 47.2 46.5 46.8 62.4 62.3 62.3 73.8 73.8 73.8 82.1 82.2 82.2
7 52.5 51.8 52.1 68.0 67.9 67.9 79.0 79.0 79.1 86.6 86.7 86.7
8 57.1 56.7 56.9 72.8 72.7 72.7 83.2 83.2 83.2 89.9 90.0 90.0
9 61.4 61.0 61.2 76.9 76.8 76.8 86.5 86.6 86.6 92.4 92.5 92.5
10 65.2 64.9 65.0 80.3 80.3 80.3 89.2 89.3 89.3 94.3 94.4 94.4
11 68.7 68.4 68.5 83.2 83.3 83.3 91.4 91.4 91.4 95.7 95.8 95.8
12 71.8 71.5 71.7 85.7 85.8 85.8 93.1 93.1 93.1 96.8 96.8 96.8
13 74.6 74.3 74.5 87.9 87.9 87.9 94.4 94.5 94.5 97.6 97.6 97.6
14 77.1 76.9 77.0 89.7 89.7 89.7 95.5 95.6 95.6 98.2 98.2 98.2
15 79.4 79.2 79.4 91.2 91.3 91.3 96.4 96.5 96.5 98.6 98.7 98.7


Primeiramente, vale ressaltar que nem todas as amostras resultaram em estimativas válidas para os parâmetros. Especificamente, quando N=100 e Cenário P1, obteve-se um total de 42 amostras sem estimativas quando k=5 e cinco amostras sem estimativas quando K=6, número relativamente baixo se comparado ao total de amostras simuladas (10.000). É conhecido em modelos de captura-recaptura que quando não se observa animais marcados em nenhuma época, a estimativa de máxima verossimilhança de N é infinito, ou extremamente elevada quando o número de recapturados é muito baixo que, em geral, leva a não convergência dos processos de maximização. Geralmente, baixos números de animais recapturados são observados em cenários com N e k pequenos e baixa probabilidade de captura, em concordância com os cenários observados sem estimativas no nosso estudo de simulação. De fato, as medidas de desempenho foram calculadas utilizando somente as estimativas válidas para os parâmetros.

Discutemos em sequência o comportamento do estimador de máxima verossimilhança de N.
A Figura 4.1 abaixo apresenta os Viéses Relativo e o EQM para as estimativas de N do modelo Mt. Visualmente, observamos que conforme aumentam o número de épocas de captura, a estimação do tamanho da população torna-se mais precisa no sentido que diminui o viés relativo e o erro quadrático médio. Como esperado, o Cenário P1 acarreta em viéses e EQM’s maiores se comparado aos demais cenários pois tem probabilidades de captura pj mais baixas, variando de 0,05 a 0,15. Em geral, podemos assumir uma boa precisão das estimativas para k10 para todos os cenários, até mesmo para o cenário P1, o que corresponde aproximadamente pelo menos 60% dos animais capturados da população (veja Tabela 4.9).


Figura: 4.1 Medidas de desempenho das estimativas do tamanho da população, N para o modelo Mt.

grafico 01


4.2.2 Modelo Mtb

Com o propósito de analisar o comportamento das estimativas dos parâmetros N e c para o modelo Mtb, foi realizado um estudo de simulação com 10000 amostras simuladas para cada valor verdadeiro fixado do parâmetro de efeito comportamental, c=(1.5,1.0,0.5,0,0.5,1.0,1.5), número de épocas de captura, k=(5,6,7,8,9,10,11,12,13,14,15), e cada tamanho populacional, N=(100,300,500). Fixamos a geração dos pjUniforme(0,10; 0,20) para todos os cenários acima. Considerando a combinação das quantidades acima, temos 3×7×11=231 cenários diferentes.

A Tabela 4.10 exibe a porcentagem média de animais distintos observados em cada cenário. Note que são bem parecidos ao Cenário P2 do estudo de simulação anterior, pois utilizamos a mesma variação para os pj. Outro fato que observados é que não há diferença no número de animais distintos observados de acordo com o parâmetro c e por isso exibimos na tabela apenas quatro cenários diferentes para c. De fato, esse parâmetro modifica a probabilidade de capturar os animais marcados, que não influencia a probabilidade de capturar novos animais.


Tabela: 4.10 Relação em porcentagem do número médio de animais distintos observados dividido pelo tamanho da população no primeiro estudo de.
c=1,5 c=0,5 c=0,5 c=1,5
k N=100 300 500 100 300 500 100 300 500 100 300 500
5 55.8 55.6 55.7 55.8 55.6 55.7 55.8 55.6 55.7 55.8 55.6 55.6
6 62.4 62.3 62.3 62.4 62.3 62.3 62.4 62.3 62.3 62.4 62.3 62.3
7 68.0 67.9 67.9 68.0 67.9 67.9 68.0 67.9 68.0 68.0 67.9 67.9
8 72.8 72.7 72.8 72.8 72.7 72.8 72.8 72.7 72.8 72.8 72.7 72.7
9 76.9 76.8 76.8 76.9 76.8 76.8 76.9 76.8 76.8 76.9 76.8 76.8
10 80.3 80.3 80.3 80.3 80.3 80.3 80.3 80.3 80.3 80.3 80.3 80.3
11 83.2 83.2 83.2 83.2 83.2 83.2 83.2 83.2 83.3 83.2 83.3 83.3
12 85.7 85.7 85.8 85.7 85.7 85.8 85.7 85.7 85.8 85.7 85.8 85.8
13 87.9 87.9 87.9 87.9 87.9 87.9 87.9 87.9 87.9 87.9 87.9 87.9
14 89.7 89.7 89.7 89.7 89.7 89.7 89.7 89.7 89.7 89.7 89.7 89.7
15 91.2 91.2 91.3 91.2 91.2 91.3 91.2 91.2 91.3 91.2 91.3 91.3


A Figura 4.2 exibe a quantidade de estimativas que convergiram para cada cenário. Diferente do que se observou no primeiro estudo de simulação, o modelo Mtb apresentou maior proporção de amostras que não tiveram estimativas válidas, principalmente para valores de k pequenos. Não se observou diferença aparente para os diferentes valores do parâmetro N. Contudo, se observou no geral que pelo menos 90% das amostras convergiram e que este número aumenta significativamente com o aumento de k.


Figura: 4.2 Quantidade das estimativas válidas de cada cenário para o modelo Mtb.

grafico 02


Analisando a Figura 4.3, podemos verificar um comportamento interessante do estimador do tamanho populacional: quando N=100, tem-se subestimativas para N quando k é pequeno, passando para estimativas viesadas (viés positivo) quando k está em torno de 8 e quando k está acima de 10 tem-se que praticamente os viéses relativos convergem para zero. Quando N=300 ou N=500, as subestimativas são geradas praticamente só quando c verdadeiro é negativo. No geral, observa-se a forte dependência do estimador de N com o valor verdadeiro do parâmetro c, e que o viés relativo tende a zero e EQM tende a diminuir quando k aumenta, principalmente para k10.


Figura: 4.3 Medidas de desempenho das estimativas do tamanho da população, N, para o modelo Mtb.

grafico 03


A Figura 4.4 mostra um comportamento de subestimativas para o parâmetro c, principalmente para valores de c negativos. Contudo, observa-se para crescimento de N e k, o estimador torna-se mais preciso no sentido que viés e EQM tendem a zero.


Figura: 4.4 Medidas de desempenho das estimativas do parâmetro de efeito comportamental à marcação, c, para o modelo Mtb.

grafico 04


Referências

Coulombe, Harry N, Sam H Ridgway, e William E Evans. 1965. “Respiratory water exchange in two species of porpoise”. Science 149 (3679): 86–88.
Otis, David L, Kenneth P Burnham, Gary C White, e David R Anderson. 1978. “Statistical inference from capture data on closed animal populations”. Wildlife monographs, nº 62: 3–135.
Wang, Xiaoyin, Zhuoqiong He, Dongchu Sun, et al. 2015. “Bayesian Estimation of Population Size via Capture-Recapture Model with Time Variation and Behavioral Response”. Open Journal of Ecology 5 (01): 1.