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:ce227-2014-01:historico [2014/04/09 22:00] paulojus |
disciplinas:ce227-2014-01:historico [2016/05/18 09:51] (atual) paulojus |
||
---|---|---|---|
Linha 14: | Linha 14: | ||
===== Conteúdos das Aulas ===== | ===== Conteúdos das Aulas ===== | ||
- | ^ ^^ B & M ^^ Online ^ | ||
^ Data ^ Conteúdo ^ Leitura ^ Exercícios ^ Tópico ^ | ^ Data ^ Conteúdo ^ Leitura ^ Exercícios ^ Tópico ^ | ||
| 12/02 Qua |Teorema de Bayes: revisão, interpretações e generalização. Expressão probabilística de informação subjetiva, estimação baseada nos dados, estimação combinando informação prévia (subjetiva) e dados. Exemplo Binomial-Beta | [[#12/02|Ver abaixo]] | [[#12/02|Ver abaixo]] | | | | 12/02 Qua |Teorema de Bayes: revisão, interpretações e generalização. Expressão probabilística de informação subjetiva, estimação baseada nos dados, estimação combinando informação prévia (subjetiva) e dados. Exemplo Binomial-Beta | [[#12/02|Ver abaixo]] | [[#12/02|Ver abaixo]] | | | ||
Linha 33: | Linha 32: | ||
| 04/04 Sex |Extensões e fundamentos do algorítmo de Gibbs. Gibbs metropolizado. Predição. Introdução do [[http:mcmc-jags.sourceforge.net|JAGS]] e rjags | | |[[#04/04|Ver abaixo]] | | | 04/04 Sex |Extensões e fundamentos do algorítmo de Gibbs. Gibbs metropolizado. Predição. Introdução do [[http:mcmc-jags.sourceforge.net|JAGS]] e rjags | | |[[#04/04|Ver abaixo]] | | ||
| 09/04 Qua |Inferência Bayesiana em modelos lineares generalizados (ex distribuição de Poisson). Introdução aos modelos de efeitos aleatórios - comparações com modelos de efeitos fixos e verossimilhança | | |[[#09/04|Ver abaixo]] | | | 09/04 Qua |Inferência Bayesiana em modelos lineares generalizados (ex distribuição de Poisson). Introdução aos modelos de efeitos aleatórios - comparações com modelos de efeitos fixos e verossimilhança | | |[[#09/04|Ver abaixo]] | | ||
+ | | 11/04 Sex |Modelos de efeitos aleatórios/hierárquicos e inferência Bayesiana. Algorítmos de inferência e INLA | | |[[#11/04|Ver abaixo]] | | ||
+ | | 16/04 Qua |Avaliação semanal. Discussão da avaliação sobre modelos declarados em JAGS. Distribuições marginais <latex>[Y] = \int [Y|\theta][\theta] d\theta</latex> e interpretação. Exemplo: caso Poisson | |[[#16/04|Ver abaixo]] | | | ||
+ | | 23/04 Qua |sem aula expositiva. Preparar e discutir [[#16/04|atividades propostas na última aula]] para apresentação na próxima aula. | | | | | ||
+ | | 25/04 Sex |Avaliação: apresentação das análises de dados correspondentes aos {{:disciplinas:ce227:ce227-av04.pdf|modelos da quarta avaliação semanal}} | | | | | ||
+ | | 30/04 Sex |continuação das apresentações e discussões | | | | | ||
+ | | 07/05 Qua |Modelos para efeitos aleatórios dependentes. Estrutura do modelo, comparação com outras estratégias/modelos. Exemplo: série temporal para dados binários. | | |[[#07/05|Ver abaixo]] | | ||
+ | | 09/05 Sex |Fundamentos do INLA | | |{{:disciplinas:ce227:apresentacao-inla.pdf|Apresentação}} | | ||
+ | | 14/05 Sex |Fundamentos do INLA - II - Comentários sobre a (excelente!) apresentação de Gianluca Baio | | |[[http://www.statistica.it/gianluca/Talks/INLA.pdf|A apresentação]] | | ||
+ | |||
Linha 195: | Linha 203: | ||
* Incluir bandas de predição usuais e bayesianas no gráfico | * Incluir bandas de predição usuais e bayesianas no gráfico | ||
* Generalizar código para família de priori Normal-GammaInversa | * Generalizar código para família de priori Normal-GammaInversa | ||
- | * Verificar o efeito/sensitividade à espeficicação das prioris | + | * Verificar o efeito/sensitividade à especificação das prioris |
**Exemplos JAGS/rjags** | **Exemplos JAGS/rjags** | ||
Linha 309: | Linha 317: | ||
=== 09/04 === | === 09/04 === | ||
- | - Fazer uma análise Bayesiana de uma regressão de Poisson e comparar com resultados de um MLG (modelo linear generalizado) usual | + | - Fazer uma análise Bayesiana de uma regressão de Poisson e comparar com resultados de um MLG (modelo linear generalizado) usual. Utilize um conjunto de dados qualquer da literatura (sugestão: conjunto ''gala'' do pacote "faraway") |
- Use os dados do artigo //Confiabilidade e Precisão na Estimação de Médias// de Singer et al na [[http://www.rbes.ibge.gov.br/c/document_library/get_file?uuid=007b1342-7715-4539-9a26-e53cd9dc220b&groupId=37690208|Revista Brasileira de Estatística, v73, n. 236, jan./jun 2012]] e: | - Use os dados do artigo //Confiabilidade e Precisão na Estimação de Médias// de Singer et al na [[http://www.rbes.ibge.gov.br/c/document_library/get_file?uuid=007b1342-7715-4539-9a26-e53cd9dc220b&groupId=37690208|Revista Brasileira de Estatística, v73, n. 236, jan./jun 2012]] e: | ||
- ajuste um modelo "ingênuo" de média e resíduo | - ajuste um modelo "ingênuo" de média e resíduo | ||
Linha 316: | Linha 324: | ||
- ajuste um modelo Bayesiano. | - ajuste um modelo Bayesiano. | ||
Compare e discuta os resultados. | Compare e discuta os resultados. | ||
+ | |||
+ | === 11/04 === | ||
+ | <code R> | ||
+ | ## | ||
+ | ## Dados simulados do modelo: | ||
+ | ## Y_{ij} \sim N(\mu_{i}, \sigma^2_y) | ||
+ | ## mu_{i} = theta + b_{i} | ||
+ | ## b_{i} \sim N(0, \sigma^2_b) | ||
+ | ## que, por ser normal (com ligação identidade) | ||
+ | ## pode ser escrito por: | ||
+ | ## Y_{ij} = \beta_0 + b_{i} + \epsilon_{ij} | ||
+ | ## | ||
+ | Ngr <- 25 | ||
+ | Nobs <- 10 | ||
+ | set.seed(12) | ||
+ | sim <- data.frame(id = Ngr*Nobs, | ||
+ | gr = rep(1:Ngr, each=Nobs), | ||
+ | bs = rep(rnorm(Ngr, m=0, sd=10), each=Nobs), | ||
+ | eps = rnorm(Ngr*Nobs, m=0, sd=4) | ||
+ | ) | ||
+ | sim <- transform(sim, y = 100 + bs + eps) | ||
+ | sim | ||
+ | |||
+ | ## estimativas "naive" | ||
+ | resumo <- function(x) c(media=mean(x), var=var(x), sd=sd(x), CV=100*sd(x)/mean(x)) | ||
+ | (sim.res <- aggregate(y~gr, FUN=resumo, data=sim)) | ||
+ | var(sim.res$y[,1]) | ||
+ | mean(sim.res$y[,2]) | ||
+ | mean(sim$y) | ||
+ | |||
+ | ## Modelo de efeitos aleatórios | ||
+ | require(lme4) | ||
+ | fit.lme <- lmer(y ~ 1|gr, data=sim) | ||
+ | summary(fit.lme) | ||
+ | ranef(fit.lme) | ||
+ | coef(fit.lme)$gr - fixef(fit.lme) | ||
+ | print(VarCorr(fit.lme), comp="Variance") | ||
+ | |||
+ | ## JAGS | ||
+ | require(rjags) | ||
+ | |||
+ | sim.lst <- as.list(sim[c("gr","y")]) | ||
+ | sim.lst$N <- nrow(sim) | ||
+ | sim.lst$Ngr <- length(unique(sim$gr)) | ||
+ | mean(sim.lst$y) | ||
+ | |||
+ | cat("model{ | ||
+ | for(j in 1:N){ | ||
+ | y[j] ~ dnorm(mu[gr[j]], tau.e) | ||
+ | } | ||
+ | for(i in 1:Ngr){ | ||
+ | mu[i] ~ dnorm(theta, tau.b) | ||
+ | } | ||
+ | theta ~ dnorm(0, 1.0E-6) | ||
+ | tau.b ~ dgamma(0.001, 0.001) | ||
+ | sigma2.b <- 1/tau.b | ||
+ | tau.e ~ dgamma(0.001, 0.001) | ||
+ | sigma2.e <- 1/tau.e | ||
+ | cci <- sigma2.e/(sigma2.e+sigma2.b) | ||
+ | }", file="sim.jags") | ||
+ | |||
+ | sim.jags <- jags.model(file="sim.jags", data=sim.lst, n.chains=3, n.adapt=1000) | ||
+ | ## inits = ... | ||
+ | |||
+ | fit.jags <- coda.samples(sim.jags, c("theta", "sigma2.b", "sigma2.e", "cci"), 10000, thin=10) | ||
+ | |||
+ | summary(fit.jags) | ||
+ | plot(fit.jags) | ||
+ | |||
+ | ## | ||
+ | require(INLA) | ||
+ | |||
+ | fit.inla <- inla(y ~ f(gr) , family="gaussian", data=sim) | ||
+ | summary(fit.inla) | ||
+ | sqrt(1/fit.inla$summary.hyperpar[,1]) | ||
+ | </code> | ||
+ | |||
+ | === 16/04 === | ||
+ | - (Individual ou em duplas) Encontre um conjunto de dados para cada um dos casos na avaliação semanal de 16/04. Proceda análises não Bayesianas e Bayesianas e discuta os resultados. Trazer //scripts// para discussão na aula de 23/04 e apresentação em 25/04. | ||
+ | - Considere o modelo de verossimilhança <latex>[Y|\mu, \sigma^2] \sim N(\theta, \sigma^2)</latex> e a priori <latex>\tau = 1/\sigma^2 \sim Ga(a, b)</latex>. Mostre como obter a densidade: \\ <latex>[Y|\theta, a, b] = \frac{\Gamma((n/2)+a)}{\pi^{n/2} \Gamma(a) (\sum_i (x_i - \theta)^2 + 2b)^{(n/2)+a}}</latex>. \\ Como este resultado pode ser interpretado? | ||
+ | |||
+ | |||
+ | === 07/05 === | ||
+ | <code R> | ||
+ | require(INLA) | ||
+ | ## | ||
+ | ## Visualizado dados | ||
+ | ## | ||
+ | data(Tokyo) | ||
+ | head(Tokyo) | ||
+ | plot(y ~ time, data=Tokyo) | ||
+ | ## colocando na forma de proporção de dias com chuva | ||
+ | plot(y/2 ~ time, data=Tokyo) | ||
+ | ## | ||
+ | ## 1. Modelo "Nulo": só intercepto | ||
+ | ## estimando a probabilidade de chuva como uma constante: | ||
+ | fit.glm <- glm(cbind(y, n-y) ~ 1, family=binomial, data=Tokyo) | ||
+ | abline(h=exp(coef(fit.glm))/(1+exp(coef(fit.glm))), col=2, lty=3, lwd=3) | ||
+ | ## como, como este modelo, todos os valores preditos são iguais bastaria fazer: | ||
+ | abline(h=fitted(fit.glm)[1], col=2, lty=3, lwd=3) | ||
+ | ## | ||
+ | ## 2. Modelo com probabilidades variando no tempo (variável/processo latente) | ||
+ | ## modelando o logito(probabilidade) como um efeito aleatório correlacionado no tempo | ||
+ | ## segundo um "random walk" cíclico de ordm 2 | ||
+ | modelo = y ~ 0 + f(time, model="rw2", cyclic=T, param=c(1, 0.0001)) | ||
+ | fit <- inla(modelo, data=Tokyo, family="binomial", Ntrials=n, control.predictor=list(compute=TRUE)) | ||
+ | ## | ||
+ | names(fit) | ||
+ | head(fit$summary.fitted.values) | ||
+ | ## sobrepondo ao gráfico dos dados (moda, mediana e médi são praticamente indistinguíveis | ||
+ | with(fit, matlines(summary.fitted.values[,c(1,3:6)], lty=c(1,2,2,2,3), col=1)) | ||
+ | ## | ||
+ | ## 3. Agora o mesmo modelo nulo anterior porém jusado pelo INLA | ||
+ | ## | ||
+ | modelo0 = y ~ 1 | ||
+ | fit0 <- inla(modelo0, data=Tokyo, family="binomial", Ntrials=n, control.predictor=list(compute=TRUE)) | ||
+ | summary(fit0) | ||
+ | fit0$summary.fitted.values | ||
+ | with(fit0, matlines(summary.fitted.values[,c(1,3:6)], lty=c(1,2,2,2,3), col=2)) | ||
+ | ## | ||
+ | ## 4. Modelando usando GAM (generalised additive model) | ||
+ | ## | ||
+ | require(mgcv) | ||
+ | fit.gam <- gam(cbind(y, n-y) ~ s(time), family=binomial, data=Tokyo) | ||
+ | names(fit.gam) | ||
+ | fitted(fit.gam, se=T) | ||
+ | pred.gam <- predict(fit.gam, type="response", se=T, newdata=Tokyo["time"]) | ||
+ | names(pred.gam) | ||
+ | with(pred.gam, matlines(cbind(fit, fit+2*se.fit, fit-2*se.fit), lty=c(1,2,2), col=4)) | ||
+ | </code> | ||
+ |