https://github.com/cran/GAS
Tip revision: 5b6062863bb2ac5de9bd14604256e6be00bdf36c authored by Leopoldo Catania on 19 August 2024, 09:00:52 UTC
version 0.3.4.1
version 0.3.4.1
Tip revision: 5b60628
gamma.cpp
#include <RcppArmadillo.h>
using namespace Rcpp;
using namespace arma;
double dGAMMA(double dY, double dAlpha, double dBeta, bool bLog=false) {
double dLPDF = Rf_dgamma(dY, dAlpha, 1.0/dBeta, 1);
if(!bLog) dLPDF=exp(dLPDF);
return dLPDF;
}
//
double pGAMMA(double dY, double dAlpha, double dBeta) {
double dP = Rf_pgamma(dY, dAlpha, 1.0/dBeta, 1, 0);
return dP;
}
double qGAMMA(double dP, double dAlpha, double dBeta){
double dQ = Rf_qgamma(dP, dAlpha, 1.0/dBeta, 1 ,0 );
return dQ;
}
double rGAMMA(double dAlpha, double dBeta){
double dY = Rf_rgamma(dAlpha, 1.0/dBeta);
return dY;
}
arma::vec mGAMMA(double dAlpha, double dBeta){
arma::vec vMoments(4);
vMoments(0) = dAlpha/dBeta;
vMoments(1) = dAlpha/pow(dBeta,2.0);
vMoments(2) = 2.0*pow(dAlpha,0.5);
vMoments(3) = 6.0/dAlpha + 3.0;
return vMoments;
}
arma::vec gamma_Score(double dY, arma::vec vTheta){
double dAlpha = vTheta(0);
double dBeta = vTheta(1);
arma::vec vScore(2);
double dAlpha_s = log(dBeta) + log(dY) - Rf_digamma(dAlpha);
double dBeta_s = dAlpha/dBeta - dY;
vScore(0) = dAlpha_s;
vScore(1) = dBeta_s;
return vScore;
}
arma::mat gamma_IM(arma::vec vTheta){
double dAlpha = vTheta(0);
double dBeta = vTheta(1);
arma::mat mIM(2,2);
mIM(0,0) = - Rf_trigamma(dAlpha);
mIM(0,1) = 1.0/dBeta;
mIM(1,0) = mIM(0,1);
mIM(1,1) = -dAlpha/pow(dBeta,2.0);
return mIM;
}