Não foi possível enviar o arquivo. Será algum problema com as permissões?
Diferenças
Aqui você vê as diferenças entre duas revisões dessa página.
Ambos lados da revisão anterior Revisão anterior Próxima revisão | Revisão anterior | ||
disciplinas:ce089-2014-02 [2014/11/10 21:22] walmes |
disciplinas:ce089-2014-02 [2014/12/15 14:16] (atual) walmes |
||
---|---|---|---|
Linha 22: | Linha 22: | ||
---- | ---- | ||
+ | /* | ||
==== Histórico das Aulas do Curso ==== | ==== Histórico das Aulas do Curso ==== | ||
Linha 67: | Linha 67: | ||
---- | ---- | ||
+ | |||
+ | */ | ||
<code R> | <code R> | ||
Linha 137: | Linha 139: | ||
plot(profile(m0)) | plot(profile(m0)) | ||
+ | |||
+ | |||
+ | ##-------------------------------------------- | ||
+ | ## Implementação conforme sugestão do Wikipedia. | ||
+ | |||
+ | require(rootSolve) | ||
+ | |||
+ | n <- rbinom(1, size=200, p=0.8) | ||
+ | y <- c(rpois(n, lambda=exp(2)), rep(0, 200-n)) | ||
+ | length(y) | ||
+ | |||
+ | barx <- mean(y) | ||
+ | p <- sum(y==0)/length(y) | ||
+ | |||
+ | f <- function(lambda, L){ | ||
+ | L$barx*(1-exp(-lambda))-lambda*(1-L$p) | ||
+ | } | ||
+ | |||
+ | L <- list(barx=barx, p=p) | ||
+ | gradient(f, x=2, L=L) | ||
+ | |||
+ | curve(f(x, L=L), 0, 15); abline(h=0) | ||
+ | |||
+ | ## Newton-Raphson. | ||
+ | maxiter <- 50; i <- 1 ## Número máximo de iterações e contador. | ||
+ | tol <- 1e-5; error <- 100*tol ## Tolerância e erro inicial. | ||
+ | theta <- matrix(NA, nrow=1, ncol=1) | ||
+ | theta[1,] <- barx | ||
+ | while(i <= maxiter & error>tol){ | ||
+ | theta <- rbind(theta, rep(NA,1)) | ||
+ | G <- f(theta[i,], L) | ||
+ | H <- gradient(f=f, x=theta[i,], L=L) | ||
+ | theta[i+1,] <- theta[i,]-solve(H)%*%G | ||
+ | error <- sum(abs((theta[i+1,]-theta[i,])/theta[i+1,])) | ||
+ | i <- i+1 | ||
+ | ## print(c(theta[i,], G)) | ||
+ | print(cbind(H, G, theta[i,])) | ||
+ | } | ||
+ | |||
+ | lam <- theta[i,] ## Estimativa do lambda. | ||
+ | pii <- 1-barx/theta[i,] ## Estimativa do pi. | ||
+ | |||
+ | ##----------------------------------------------------------------------------- | ||
+ | ## Vendo os contornos da verossimilhança. | ||
+ | |||
+ | llmax <- ll(th=c(log(pii/(1-pii)), log(lam)), y=y) | ||
+ | |||
+ | th1 <- seq(-7, 5, l=50) | ||
+ | th2 <- seq(-1, 4, l=50) | ||
+ | lla <- outer(th1, th2, llv, y=y) | ||
+ | |||
+ | contour(th1, th2, lla, | ||
+ | ## levels=seq(from=llmax-30, to=llmax, by=10), | ||
+ | nlevels=20, | ||
+ | xlab="theta1: logit(p)", | ||
+ | ylab="theta2: log(lambda)") | ||
+ | abline(v=log(pii/(1-pii)), h=log(lam), col=2) | ||
</code> | </code> |