integralForBayesian.cpp

#include <Rcpp.h>

using namespace Rcpp;

double getProb(double p, int success, int fail){
	double e = success * log(p) + fail * log(1-p);
	return e;
}


//[[Rcpp::export]]
NumericVector probOfProb(NumericVector p, int success, int fail){
	int size = p.size();
	IntegerVector a = seq_len(success), b = seq_len(fail), c = seq_len(fail+success);
	NumericVector result(size), b1 = log(a), b2 = log(b), b3 = log(c);
	double c1 = sum(b1), c2 = sum(b2), c3 = sum(b3);
	double d = c1 + c2 - c3;
	for (int i = 0 ; i < size; i++){
		result[i] = exp(getProb(p[i],success,fail) - d);
	}
	return result;
}