dgamma.cpp

Go to the documentation of this file.

#include "statsLib.h"

dvariable dgamma( const dvariable &x, const double& a, const double& b )
{
  //E(x)=a/b;
  //V(x)=a/b^2
  if ( a<= 0 || b<=0 )
  {
    cerr<<"a or b parameters are less than or equal to 0 in"
      "dgamma( const dvariable &x, const double a, const double b )\n";
    return 0;
  }
  if ( x<=0 )
  {
    cerr<<"The x parameter is less than or equal to 0 in"
      "dgamma( const dvariable &x, const double a, const double b )\n";
    return 0;
  }
  
  return -a*log(b)+gammln(a)-(a-1.)*log(x)+b*x;
}

dvariable dgamma( const dvariable &x, const dvariable& a, const dvariable &b )
{
  //E(x)=a/b;
  //V(x)=a/b^2
  if ( a<= 0 || b<=0 )
  {
    cerr<<"a or b parameters are less than or equal to 0 in"
      "dgamma( const dvariable &x, const double a, const double b )\n";
    return 0;
  }
  if ( x<=0 )
  {
    cerr<<"The x parameter is less than or equal to 0 in"
      "dgamma( const dvariable &x, const double a, const double b )\n";
    return 0;
  }
  
  return -a*log(b)+gammln(a)-(a-1.)*log(x)+b*x;

}

dvariable dgamma( const dvector &x, const prevariable& a, const prevariable& b )
{
  //E(x)=a/b;
  //V(x)=a/b^2
  if ( a<= 0 || b<=0 )
  {
    cerr<<"a or b parameters are less than or equal to 0 in"
      "dgamma( const dvariable &x, const double a, const double b )\n";
    return 0;
  }
  for(int i=x.indexmin(); i<=x.indexmax(); i++)
  {
    if ( x(i)<=0 )
    {
      cerr<<"The x("<<i<<") parameter is less than or equal to 0 in"
        "dgamma( const dvariable &x, const double a, const double b )\n";
      return 0;
    }
  }
  return sum(-a*log(b)+gammln(a)-(a-1.)*log(x)+b*x);

}
dvariable dgamma( const dvector &x, const prevariable& a, const dvar_vector& b )
{
  //E(x)=a/b;
  //V(x)=a/b^2
  for(int i=x.indexmin(); i<=x.indexmax(); i++)
  {
    if ( a<= 0 || b(i)<=0 )
    {
      cerr<<"a or b ("<<i<<") parameters are less than or equal to 0 in"
        "dgamma( const dvariable &x, const double a, const double b )\n";
      return 0;
    }
  } 
  for(int i=x.indexmin(); i<=x.indexmax(); i++)
  {
    if ( x(i)<=0 )
    {
      cerr<<"The x ("<<i<<") parameter is less than or equal to 0 in"
        "dgamma( const dvariable &x, const double a, const double b )\n";
      return 0;
    }
  }
  return sum(-a*log(b)+gammln(a)-(a-1.)*log(x)+elem_prod(b, x));
}
dvariable dgamma( const dvector &x, const dvar_vector& a, const dvar_vector& b )
{
  //E(x)=a/b;
  //V(x)=a/b^2
  for(int i=x.indexmin(); i<=x.indexmax(); i++)
  {
    if ( a(i)<= 0 || b(i)<=0 )
    {
      cerr<<"a or b ("<<i<<") parameters are less than or equal to 0 in"
        "dgamma( const dvariable &x, const double a, const double b )\n";
      return 0;
    }
  } 
  for(int i=x.indexmin(); i<=x.indexmax(); i++)
  {
    if ( x(i)<=0 )
    {
      cerr<<"The x ("<<i<<") parameter is less than or equal to 0 in"
        "dgamma( const dvariable &x, const double a, const double b )\n";
      return 0;
    }
  }
  return sum(-elem_prod(a, log(b))+gammln(a)-elem_prod((a-1.), log(x))+elem_prod(b, x));
}
dvariable dgamma( const dvector &x, const dvar_vector& a, const prevariable& b )
{
  //E(x)=a/b;
  //V(x)=a/b^2
  for(int i=x.indexmin(); i<=x.indexmax(); i++)
  {
    if ( a(i)<= 0 || b<=0 )
    {
      cerr<<"a or b ("<<i<<") parameters are less than or equal to 0 in"
        "dgamma( const dvariable &x, const double a, const double b )\n";
      return 0;
    }
  } 
  for(int i=x.indexmin(); i<=x.indexmax(); i++)
  {
    if ( x(i)<=0 )
    {
      cerr<<"The x ("<<i<<") parameter is less than or equal to 0 in"
        "dgamma( const dvariable &x, const double a, const double b )\n";
      return 0;
    }
  }
  return sum(-a*log(b)+gammln(a)-elem_prod((a-1.), log(x))+b*x);
}