Bootstrap

Objetivo

Construir simulação bootstrap usando o Rcpp e comparar velocidade com a versao em R

Funções

# Função do R
ic_bootR <- function(R,nx,nboot,Rboot,make.graph = TRUE){
Io = Ib = numeric(R)
for(i in 1:R){
# -------- Amostra original -------- #
xo <- rnorm(nx)# xo eh o x original
icxo <- c(mean(xo) -1.96/sqrt(nx),mean(xo) + 1.96/sqrt(nx))
Io[i] <- ifelse(0>=icxo[1]&0<=icxo[2],1,0)
# ------------ BootStrap ----------- #
xb <- numeric(Rboot)
for(j in 1:Rboot){
xb[j] <- mean(sample(xo,size = nboot,replace = TRUE))
}
icxb <- quantile(xb,probs = c(0.025,0.975))
Ib[i] <- ifelse(0>=icxb[1]&0<=icxb[2],1,0)
}
if(make.graph==TRUE){
plot(c(1:R),cumsum(Ib)/c(1:R),
type = "l",
xlab = "Nº de réplicas",
ylab = "Nível de confiança",col = "blue",
ylim = c(0,1),
main = paste("Simulação para uma amostra de tamanho",nx,
"com ",R,"réplicas de monte carlo e",
"\n","utilizando um bootstrap de tamanho",nboot,
"e ",Rboot,"réplicas."))
points(c(1:R),cumsum(Io)/c(1:R),type = "l")
abline(h = 0.95,lty = 2)
legend("bottomright",legend = c("Bootstrap",
"Amostra Aleátoria",
"Nível de conf. = 95%"),
col = c("blue","black","black"),lty = c(1,1,2))
}
return(data.frame(original = Io,bootstrap = Ib))
}


#// Funcao em C++
#include <RcppArmadillo.h>
#using namespace Rcpp;
#// [[Rcpp::depends(RcppArmadillo)]]
#// Calculo dos quantis
#// [[Rcpp::export]]
#double quantileC(NumericVector x,double q){
#NumericVector y = clone(x);
#std::sort(y.begin(),y.end());
#return y[x.size()*(q - 0.000000001)];
#}
#// [[Rcpp::export]]
#DataFrame IcBootstrap(int R,int nx,int nboot,int Rboot){
#NumericVector Io(R);
#NumericVector Ib(Rboot);
#NumericVector xb(Rboot);
#double icxo1, icxo2;
#double icxb1, icxb2;
#for(int i=0;i<R;i++){
#// -------- Amostra original --------- //
#NumericVector xo = rnorm(nx); // xo eh o x original
#icxo1 = mean(xo) - (1.96/std::sqrt(nx));
#icxo2 = mean(xo) + (1.96/std::sqrt(nx));
#if((0>=icxo1)&(0<=icxo2)){
#Io(i) = 1;
#}
#else{
#Io(i) = 0;
#}
#// -------- Bootstrap -------- //
#for(int j = 0;j<Rboot;j++){
#xb(j) = mean(sample(xo,nboot,TRUE));
#icxb1 = quantileC(xb,0.025);
#icxb2 = quantileC(xb,0.975);
#if((0>=icxb1)&(0<=icxb2)){
#Ib(i) = 1;
#}
#else{
#Ib(i) = 0;
#}
#}
#}
#return DataFrame::create(Named("Io") = Io,
#Named("Ib") = Ib);
#}
# Funcao em C passada para os mesmos moldes da feita em R
#ic_bootC <- function(R,nx,nboot,Rboot,make.graph = TRUE){
#boot <- IcBootstrap(R,nx,nboot,Rboot)
#Io <- boot$Io
#Ib <- boot$Ib

#if(make.graph==TRUE){
#plot(c(1:R),cumsum(Ib)/c(1:R),
#type = "l",
#xlab = "Nº de réplicas",
#ylab = "Nível de confiança",col = "blue",
#ylim = c(0,1),
#main = paste("Simulação para uma amostra de tamanho",nx,
#"com ",R,"réplicas de monte carlo e",
#"\n","utilizando um bootstrap de tamanho",nboot,
#"e ",Rboot,"réplicas."))
#points(c(1:R),cumsum(Io)/c(1:R),type = "l")
#abline(h = 0.95,lty = 2)
#legend("bottomright",legend = c("Bootstrap",
#"Amostra Aleátoria",
#"Nível de conf. = 95%"),
#col = c("blue","black","black"),lty = c(1,1,2))
#}
#return(data.frame(original = Io,bootstrap = Ib))
#}

Teste das funções

Simulação para uma amostra de tamanho 30 com 300 réplicas de monte carlo e utilizando um bootstrap de tamanho 30 e 300 réplicas.

com_r <- ic_bootR(R = 300,nx = 30,nboot = 30,Rboot = 300)

#com_c <- ic_bootC(R = 300,nx = 30,nboot = 30,Rboot = 300)

Teste das velocidades

#library(microbenchmark)
#tempo <- microbenchmark(ic_bootR(R = 300,nx = 30,nboot = 30,Rboot = 300,make.graph = F),
#ic_bootC(R = 300,nx = 30,nboot = 30,Rboot = 300,make.graph = F),
#times = 100)
#tempo

Conclusão

Observa-se que a função do R é mais rápida comparada com a do C++, porém a variabilidade também é maior.