Revision 1618 trunk/contrib/statslib/dnbinom.cpp

dnbinom.cpp (revision 1618)
1
#include "statsLib.h"
2 1
#include <df1b2fun.h>
3 2
#include <adrndeff.h> 
4

  
5
/**
6
* 
7
* \file dnbinom.cpp
8
* \brief Negative binomial density functions.
9
* \ingroup STATLIB
10
* \author Steven Martell and Mollie Brooks
11
* \date 2/05/2011
12
* 
13
* 
14
* This file contains the negative loglikelihood 
15
* functions for the negative binomial distribution. The function
16
* dnbinom is overloaded to accomodate variable and data objects. 
17
* 
18
* The negative log function is implemented as:
19
* \f[
20
*  \ln \Gamma (n) + \ln x! - \ln \Gamma(x+n) - n \ln(p) - x \ln(1-p)
21
* \f]
22
* where \f$p\f$ is the negative binomial probability and \f$n\f$
23
* number of trials, and \f$x\f$ is the number of successes.
24
* 
25
* 
26
* 
27
*/
28
/** \brief Negative binomial with size and mean
29
	
30
	
31
	
32
	\author Steven James Dean Martell UBC Fisheries Centre
33
	\date 2011-06-23
34
	\param  x This represents the number of failures which 
35
	* occur in a sequence of Bernoulli trials before a target
36
	*  number of successes is reached. 
37
	\param  size target for number of successful trials, or
38
	*  dispersion parameter (the shape parameter of the gamma
39
	*  mixing distribution). Must be strictly positive, need not be integer.
40
	\param mu is the expected mean number of successful trials.
41
	\return the negative log likelihood for the negative binomial distribution.
42
	\sa
43
    \ingroup STATLIB
3
/** negative log likelihood of negative binomial with mean=mu and variance = mu + mu^2 /k 
4
\brief Negative binomial with mean and quadratic variance
5
\author Steven Martell and Mollie Brooks
6
\param x observed count
7
\param mu is the predicted mean
8
\param k is the overdispersion parameter, i.e. shape parameter of underlying gamma heterogeneity (different from tau). should be >0
9
\return negative log likelihood \f$ -( \ln(\Gamma(x+k))-\ln(\Gamma(k))-\ln(x!)+k\ln(k)+x\ln(\mu)-(k+x)\ln(k+\mu) )\f$
10
\ingroup STATLIB
44 11
**/
45

  
46
dvariable dnbinom(const double& x,const prevariable& mu, const prevariable& size)
12
dvariable dnbinom(const double& x, const prevariable& mu, const prevariable& k)
47 13
{
48
	if (value(size)<0.0)
14
	//x is the observed count
15
	//mu is the predicted mean
16
	//k is the overdispersion parameter
17
	if (value(k)<0.0)
49 18
	{
50
		cerr<<"size is <=0.0 in"
51
			"dnbinom(const double& x,const prevariable& mu, const prevariable& size)\n";
52
		return(0);
19
		cerr<<"k is <=0.0 in dnbinom()";
20
		return(0.0);
53 21
	}
54 22
	RETURN_ARRAYS_INCREMENT();
55
	//dvariable p=size/(size+mu);
56
	dvariable tmp;
57
	//tmp = gammln(x+size) - gammln(size) - factln(x)
58
	//	+size*log(p)+x*log(1.-p);
59
	dvariable tau = (1.0 + size*mu);
60
	tmp = log_negbinomial_density(x,mu,tau);
61
		
23
	dvariable loglike;
24

  
25
	loglike = gammln(k+x)-gammln(k)-gammln(x+1)+k*log(k)-k*log(mu+k)+x*log(mu)-x*log(mu+k);
26

  
62 27
	RETURN_ARRAYS_DECREMENT();
63
	return(-tmp);
28
	return(-loglike);
64 29
}
65 30

  
66
/** negative log likelihood of negative binomial with mean and size 
67
\brief Negative binomial with size and mean
31
/** negative log likelihood of negative binomial with mean=mu and variance = mu + mu^2 /k 
32
\brief Negative binomial with mean and quadratic variance
68 33
\author Mollie Brooks
69 34
\param x observed count
70 35
\param mu is the predicted mean
71
\param k is the overdispersion parameter, i.e. shape parameter of underlying heterogeneity (different from tau). should be >0
36
\param k is the overdispersion parameter, i.e. shape parameter of underlying gamma heterogeneity (different from tau). should be >0
72 37
\return negative log likelihood \f$ -( \ln(\Gamma(x+k))-\ln(\Gamma(k))-\ln(x!)+k\ln(k)+x\ln(\mu)-(k+x)\ln(k+\mu) )\f$
73 38
\ingroup STATLIB
... This diff was truncated because it exceeds the maximum size that can be displayed.

Also available in: Unified diff