**Código em construção, ainda não é a versão final!** # Dados para teste x <- c(0,0,0,1,0,0,NA,NA,1,0,0,1,NA,1,NA,NA,0,NA,NA,0,0,0,0,NA,0,0,0,1,0,NA,1,NA,0,1,0,1,NA,0,NA,0,0,0,1,0,0,0,NA,1,0,NA,0,0,NA,1,0,1,1,0,0,0,1,1,1,0,NA,NA,NA,1,NA,NA,NA,NA,NA,NA,0,0,0,0,NA,0,1,0,0,1,0,0,0,1,1,1,1,1,NA,NA,0,0,0,NA,0,NA,0,1,0,1,0,0,1,0,0,NA,NA,0,1,NA,NA,0,0,1,0,1,1,1,1,1,NA,0,1,NA,0,1,1,0,1,0) y <- c(0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,1,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,1,0,0,0,1,0,0,0,0,0,0,0,0,0,0,0,1,0,0,0,1,1,0,1,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,1,0,0,0,0,0,0,1,0,1,1,0,0,0,0,0,0,0,0,0,0) prev <- function(x,y, datatype=c("rrt", "qd", "RRT.qd"), probs=c(0.25, 0.5833), nsim=1000) # função prev com argumentos: x, y, probs (probabilidades de "sim" forçado e "diga a verdade", ver Hox & Lensvelt-Mulders 2005), teste (se usuário quer realizar teste de permutação), nsim (número de simulações) { if(datatype=="rrt") { rrt <- function(x) ((sum(na.omit(x))/length(na.omit(x)))-probs[1])/probs[2] # fórmula de Hox & Lensvelt-Mulders (2005) para estimar prevalência de comportamento com base em dados RRT a <- rrt(x) # criando objeto a com estimativa de prevalência para os dados de respostas randomizadas rrt.boot <- rep(NA, nsim) # criando objeto rrt.boot para receber valores de simulação for(i in 1:nsim) # criando contador { sample.x <- sample(x, replace=TRUE) # reamostrando vetor de dados de respostas randomizadas com substituição rrt.boot[i] <- rrt(sample.x) # preenchendo objeto rrt.boot } b <- quantile(rrt.boot, probs=c(0.025, 0.975), na.rm=T) # criando objeto b com quantis referentes ao intervalo de confiança de 95% n.NAx <- (length(x))-(length(na.omit(x))) # calcula numero de NAs no vetor de dados x cat("\n\t", n.NAx," valores NA omitidos de x\n") # insere texto informando numero de valores de NA omitidos de x ret.rrt <- list("Prevalencia RRT" = a, "intervalo de confiança de 95% RRT" = b) # cria objeto com lista de itens a serem exibidos no return return(ret.rrt) } if(datatype=="qd") { qd <- function(y) sum(na.omit(y))/length(na.omit(y)) # fórmula para estimar prevalência de comportamento com base em dados de questionamento direto a2 <- qd(y) # criando objeto a2 com estimativa de prevalência para os dados de questionamento direto qd.boot <- rep(NA, nsim) # criando objeto qd.boot para receber valores de simulação for(i in 1:nsim) # criando contador { sample.y <- sample(y, replace=TRUE) # reamostrando vetor de dados de questionamento direto com substituição qd.boot[i] <- qd(sample.y) # preenchendo objeto qd.boot } b2 <- quantile(qd.boot, probs=c(0.025, 0.975), na.rm=T) # criando objeto b2com quantis referentes ao intervalo de confiança de 95% n.NAy <- (length(y))-(length(na.omit(y))) # calcula numero de NAs no vetor de dados y cat("\t", n.NAy," valores NA omitidos de y\n\n") ret.qd <- list("Prevalencia QD" = a2, "intervalo de confiança de 95% QD" = b2) # cria objeto com lista de itens a serem exibidos no return return(ret.qd) } if(datatype=="rrt.qd") { rrt <- function(x) ((sum(na.omit(x))/length(na.omit(x)))-probs[1])/probs[2] # fórmula de Hox & Lensvelt-Mulders (2005) para estimar prevalência de comportamento com base em dados RRT a <- rrt(x) # criando objeto a com estimativa de prevalência para os dados de respostas randomizadas rrt.boot <- rep(NA, nsim) # criando objeto rrt.boot para receber valores de simulação for(i in 1:nsim) # criando contador { sample.x <- sample(x, replace=TRUE) # reamostrando vetor de dados de respostas randomizadas com substituição rrt.boot[i] <- rrt(sample.x) # preenchendo objeto rrt.boot } b <- quantile(rrt.boot, probs=c(0.025, 0.975), na.rm=T) # criando objeto b com quantis referentes ao intervalo de confiança de 95% qd <- function(y) sum(na.omit(y))/length(na.omit(y)) # fórmula para estimar prevalência de comportamento com base em dados de questionamento direto a2 <- qd(y) # criando objeto a2 com estimativa de prevalência para os dados de questionamento direto qd.boot <- rep(NA, nsim) # criando objeto qd.boot para receber valores de simulação for(i in 1:nsim) # criando contador { sample.y <- sample(y, replace=TRUE) # reamostrando vetor de dados de questionamento direto com substituição qd.boot[i] <- qd(sample.y) # preenchendo objeto qd.boot } b2 <- quantile(qd.boot, probs=c(0.025, 0.975), na.rm=T) # criando objeto b2com quantis referentes ao intervalo de confiança de 95% boxplot(rrt.boot, qd.boot, xlab="Metodo de amostragem", ylab="Prevalência", xaxt="n") # boxplot com estimativas de prevalência obtidas pelos dois métodos axis(side=1, at=(1:2), labels=c("RRT", "QD")) # adicionando labels ao boxplot dif <- rrt.boot-qd.boot # calculando a diferença entre as estimativas de prevalência obtidas pelos dois métodos n <- sum(abs(dif) >= a-a2) # calculando quantas vezes o módulo do valor observado foi obtido ao acaso prob <- n/length(dif) # calculando a probabilidade de o valor observado se dever ao acaso x11() # abrindo novo graphic device hist(dif, ylab="Frequência", xlab="RRT-QD", main=NULL) # plota histograma com diferenças simuladas entre as estimativas de prevalência obtidas pelos dois métodos abline(v=a-a2, col="red") # indica a posição do valor observado no histograma n.NAx <- (length(x))-(length(na.omit(x))) # calcula numero de NAs no vetor de dados x cat("\n\t", n.NAx," valores NA omitidos de x\n") # insere texto informando numero de valores de NA omitidos de x n.NAy <- (length(y))-(length(na.omit(y))) # calcula numero de NAs no vetor de dados y cat("\t", n.NAy," valores NA omitidos de y\n\n") ret.rrt.qd <- list("Prevalencia RRT" = a, "intervalo de confiança de 95% RRT" = b, "Prevalencia QD" = a2, "intervalo de confiança de 95% QD" = b2, "Probabilidade de diferença entre estimativas RRT e QD ter sido obtida ao acaso"=prob) # cria objeto com lista de itens a serem exibidos no return return(ret.rrt.qd) } } # testando: prev(x, datatype="rrt") prev(,y, datatype="qd") prev(x, y, datatype="rrt.qd")