Atividade: Refazer o algoritmo apresentado em Rcpp. A função deve retornar um dataframe com as variâncias para os casos de N 10, 100, 200, 1000. Tanto para o tradicional e com usando a técnica de redução de variância (Cap 20. Seção de variáveis antitéticas do Jones.)
#include <RcppArmadillo.h>
using namespace Rcpp;
// [[Rcpp::depends(RcppArmadillo)]]
// [[Rcpp::export]]
arma::mat g(arma::mat m) {
arma::mat X = m.transform( [](double val) { return (1 - pow(val,2)); } );
return(X);
}
// [[Rcpp::export]]
arma::vec var(int N=5000,int n =50){
arma::mat u1 = arma::randu<arma::mat>(2*n,N);
arma::mat theta1 = arma::mean(g(u1),0);
arma::vec aux = arma::var(theta1,0,1);
double var1 = aux(0);
//
arma::mat ua = arma::randu<arma::mat>(n,N);
arma::mat thetaa = 0.5*( (arma::mean(g(ua),0)) + (arma::mean(g(1-ua),0)) );
arma::vec aux2 = arma::var(thetaa,0,1);
double vara = aux2(0);
//
double reduction = (100*(var1 - vara)/var1);
arma::vec variances = {var1,vara,reduction};
return (variances);
}
// [[Rcpp::export]]
Rcpp::DataFrame variancias(){
arma::Col<int> Nsize = {10,100,200,1000};
arma::vec size1 = var(Nsize(0),50);
arma::vec size2 = var(Nsize(1),50);
arma::vec size3 = var(Nsize(2),50);
arma::vec size4 = var(Nsize(3),50);
arma::vec v1 = {size1(0),size2(0),size3(0),size4(0)};
arma::vec va = {size1(1),size2(1),size3(1),size4(1)};
arma::Col<int> redu = {(int)round(size1(2)),(int)round(size2(2)),(int)round(size3(2)),(int)round(size4(2))};
Rcpp::DataFrame result = Rcpp::DataFrame::create(Rcpp::Named("N") = Nsize,
Rcpp::Named("Var.1") = v1,
Rcpp::Named("Var.A") = va,
Rcpp::Named("Reduction") = redu);
return (result);
}
O código escrito acima é a versão em Rcpp da mesma função descrita no Jones da técnica de redução de variancia usando Variáveis Antitéticas. O código original em R se encontra no livro do Jones. Como saída, geramos um data frame como foi solicitado. A função ‘g’ transforma os dados. A função ‘var’ estima as variancias e a redução ao utilizar a técnica. Já a função ‘variancias’ gera o dataframe com as informações calculadas. Ao chamar a função variancias obtemos o dataframe que pode ser visto a seguir.
# variancias()
# N Var.1 Var.A Reduction
# 1 10 0.0019028962 4.744384e-05 98
# 2 100 0.0009137750 1.024020e-04 89
# 3 200 0.0008403715 9.802160e-05 88
# 4 1000 0.0008166650 1.112221e-04 86