########################################################################################## ## Solução da lista II - Conglomerados em 2 estágios ##################################### ########################################################################################## y <- c(4,6,2,5,9,5,8,9,8,4,5,1,0,6,4,4,1,4,7,8) cong <- rep(1:10,each=2) dados <- data.frame(y,cong) B.alpha <- rep(50,10) ## Primeiro estimador y.2c1 <- function(dados, B.alpha){ y.alpha <- tapply(dados$y,dados$cong, mean) a <- length(unique(dados$cong)) B.barra <- mean(B.alpha) y2c1 <- (1/a)*sum( (B.alpha/B.barra)*y.alpha) return(y2c1) } var.y2c1 <- function(dados,B.alpha){ B.barra <- mean(B.alpha) y.alpha <- tapply(dados$y, dados$cong, mean) a <- length(unique(dados$cong)) s2 <- (1/(a-1))*sum( ( (B.alpha/B.barra)*y.alpha - y.2c1(dados=dados,B.alpha=B.alpha))^2) s2y2c1 <- s2/a return(s2y2c1) } ########################################################################################## ## Resultados com y.2c1 ################################################################## ########################################################################################## pontual.1 <- y.2c1(dados,B.alpha) var.1 <- var.y2c1(dados,B.alpha) IC.1 <- pontual.1 + c(-1,1)*qnorm(0.975)*sqrt(var.1) IC.1 ########################################################################################## ## Segundo estimador ##################################################################### ########################################################################################## y.2c2 <- function(dados){ b.alpha <- tapply(dados$y, dados$cong, length) b.barra <- mean(b.alpha) y.alpha <- tapply(dados$y,dados$cong, mean) y2c2 <- sum(b.alpha*y.alpha)/sum(b.alpha) return(y2c2) } ########################################################################################## ## Variância para o segundo estimador #################################################### ########################################################################################## var.y2c2 <- function(dados){ b.alpha <- tapply(dados$y, dados$cong, length) b.barra <- mean(b.alpha) y.alpha <- tapply(dados$y,dados$cong, mean) a <- length(unique(dados$cong)) s2 <- (1/(a-1))*sum( (b.alpha/b.barra)^2*(y.alpha - y.2c2(dados=dados))^2) s2y2c2 <- s2/a return(s2y2c2) } pontual.2 <- y.2c2(dados) var.2 <- var.y2c2(dados) IC.2 <- pontual.2 + c(-1,1)*qnorm(0.975)*sqrt(var.2) IC.2 ########################################################################################## ## Estimador 3 - ######################################################################### ########################################################################################## y.2c3 <- function(dados){ y.alpha <- tapply(dados$y,dados$cong, mean) y2c3 <- mean(y.alpha) return(y2c3) } ########################################################################################## ## Variancia estimador 3 ################################################################# ########################################################################################## var.y2c3 <- function(dados){ y.alpha <- tapply(dados$y,dados$cong, mean) a <- length(unique(dados$cong)) s2 <- (1/(a-1))*sum( (y.alpha - y.2c3(dados))^2) s2y2c3 <- s2/a return(s2y2c3) } pontual.3 <- y.2c3(dados) var.3 <- var.y2c3(dados) IC.3 <- pontual.3 + c(-1,1)*qnorm(0.975)*sqrt(var.3) IC.3