Aqui você vê as diferenças entre duas revisões dessa página.
biometria:verossim:03b-model [2009/06/10 00:47] prado |
biometria:verossim:03b-model [2015/08/10 20:48] |
||
---|---|---|---|
Linha 1: | Linha 1: | ||
- | ^[[start|BIE 5187 - Modelos Estatísticos na Pesquisa em Ecologia e Recursos Naturais]]^ | ||
- | |||
- | \\ | ||
- | < | ||
- | <font face=" | ||
- | 4. Modelos com Parâmetros Constantes | ||
- | </ | ||
- | </ | ||
- | |||
- | \\ | ||
- | |||
- | |||
- | ====== Tutoriais ====== | ||
- | |||
- | ===== Poisson e Binomial Negativa ===== | ||
- | |||
- | Vamos retomar [[biometria: | ||
- | |||
- | < | ||
- | set.seed(42) | ||
- | cx <- runif(20, | ||
- | cy <- runif(20, | ||
- | px <- rnorm(2000) | ||
- | py <- rnorm(2000) | ||
- | x1 <- cx+px | ||
- | y1 <- cy+py | ||
- | x2 <- x1[x1> | ||
- | y2 <- y1[x1> | ||
- | x2.parc <- cut(x2, | ||
- | y2.parc <- cut(y2, | ||
- | </ | ||
- | |||
- | O numero de plantas por parcela é obtido com: | ||
- | |||
- | < | ||
- | cont2 <- data.frame(table(x2.parc, | ||
- | </ | ||
- | |||
- | A avaliação visual feita no tutorial anterior indica que a variável binomial negativa é uma descrição mais adequada destes dados. Vamos ajustar estes dois modelos e fazer a comparação. | ||
- | |||
- | Primeiro definimos funções de log-verossimilhança negativa para cada modelo: | ||
- | |||
- | < | ||
- | LL.pois <- function(lam){ | ||
- | -sum(dpois(cont2, | ||
- | } | ||
- | |||
- | LL.nbin <- function(media, | ||
- | -sum(dnbinom(cont2, | ||
- | } | ||
- | </ | ||
- | |||
- | Em seguida minimizamos estas funções. Para isto carregue o pacote //bbmle// e use a função [[http:// | ||
- | < | ||
- | library(bbmle) # basta uma vez por seção | ||
- | mod1 <- mle2(LL.pois, | ||
- | mod2 <- mle2(LL.nbin, | ||
- | </ | ||
- | |||
- | Como esperado, a binomial negativa é um modelo muito mais plausível: | ||
- | |||
- | < | ||
- | > logLik(mod1) | ||
- | 'log Lik.' -1786.502 (df=1) | ||
- | > logLik(mod2) | ||
- | 'log Lik.' -1014.969 (df=2) | ||
- | </ | ||
- | |||
- | Você pode obter um resumo dos modelo com o comando '' | ||
- | < | ||
- | summary(mod1) | ||
- | summary(mod2) | ||
- | </ | ||
- | |||
- | E podemos fazer um gráfico dos valores previsto pelos dois modelos com: | ||
- | |||
- | < | ||
- | ##MLEs de cada modelo | ||
- | cf1 <- coef(mod1) | ||
- | cf2 <- coef(mod2) | ||
- | |||
- | ##grafico | ||
- | cont2.f <- factor(cont2, | ||
- | plot(table(cont2.f)/ | ||
- | points(x=0: | ||
- | points(x=0: | ||
- | |||
- | </ | ||
- | |||
- | Por fim, avalie o perfil de verossimilhança dos dois parâmetros da binomial negativa ((para isto você vai precisar da função {{: | ||
- | |||
- | < | ||
- | p.mod2 <- profile(mod2) | ||
- | par(mfrow=c(1, | ||
- | plot.profmle(p.mod2) | ||
- | par(mfrow=c(1, | ||
- | </ | ||
- | |||
- | ====== Recursos para Estudo ====== | ||
- | ===== Na Internet ===== | ||
- | |||
- | * Um roteiro de ajuste de modelos na página do [[http:// | ||
- | * Excelentes exercícios de simulação e ajustes de distribuições no [[http:// | ||
- | |||
- | |||
- | ====== Pesquisa ====== | ||
- | |||
- | Indique os tutoriais e partes da leitura básica que merecem mais atenção na discussão [[http:// | ||
- | |||
- | |||
- | |||
- | |||
- | |||