writing a function for columnwise computation of sample skewness in Rcpp I'm in trouble using sqrt()-function. I know that sqrt(x) works for NumericVector types (tested it in a separate file), however in my code (where I'm trying to pass doubles) it does not work.
Here's my code:
#include <Rcpp.h>
#include <cmath>
using namespace Rcpp;
// [[Rcpp::export]]
NumericVector colSkew(NumericMatrix x) {
int nc = x.ncol();
int nr = x.nrow();
NumericVector colS(nc);
for(int i = 0; i < nc; i++){
double cMean = mean(x( _ ,i));
double xSq = 0;
double cSt = 0;
for(int j = 0; j < nr; j++){
xSq += std::pow(x(j,i), 2.0);
cSt += std::pow(x(j,i) - cMean, 3.0);
}
double colMsq = nr * std::pow(cMean, 2.0);
double cTT = sqrt((xSq - colMsq)) / (nr - 1);
double colVar = cTT / (nr - 1);
double colNew = nr * std::pow(colVar, 3);
colS[i] = cSt / colNew;
}
return(colS);
}
I tried std::sqrt()
instead, giving me the "call to sqrt is ambiguous" error. For that one, neither Getting: error C2668: 'sqrt' : ambiguous call to overloaded function, nor http://dirk.eddelbuettel.com/code/rcpp.examples.html's
inline static double sqrt_double( double x ){ return ::sqrt( x ); }
helped. (The latter one has itself either the "no matching function" or the "ambiguous" problem). The code itself compiles, but does not work properly (I can call the function but results are not fine).
Reproducible example:
# First load the cpp function, however you want to;
Create example data:
A = matrix(rchisq(1000, 5), nrow = 100)
library(timeDate)
skew.1 = apply(A, 2, skewness)
skew.2 = colSkew(A)
# A custom R-function which should mimic the cpp function
colSkew.r = function(x){
nc = ncol(x)
nr = nrow(x)
colS = numeric(nc)
cMean = colMeans(x)
xSq = colSums(x^2)
cSt = 0
for(i in 1:nc){
cSt = sum((x[,i]-cMean[i])^3)
colMsq = nr * cMean[i]^2
cTT = sqrt((xSq[i] - colMsq) / (nr - 1))
colNew = nr * cTT^3
colS[i] = cSt / colNew
}
return(colS)
}
skew.3 = colSkew.r(A)