Normalizing constant for the hyperdirichlet distribution

Numerical techniques for calculating the normalizing constant for the hyperdirichlet distribution

B(H, disallowed=NULL, give=FALSE, ...)
probability(H, disallowed=NULL, ...)
mgf(H, powers, ...) 
dhyper2(ip,H,...)
dhyper2_e(e,H,include.Jacobian=TRUE)
mean_hyper2(H, normalize=TRUE, ...)
Jacobian(e)
e_to_p(e)
p_to_e(p)

Arguments

H

Object of class hyper2

powers

Vector of length dim(x) whose elements are the powers of the expectation; see details section

disallowed

Function specifying a subset of the simplex over which to integrate; default NULL means to integrate over the whole simplex. The integration proceeds over p with disallowed(p) evaluating to FALSE

e,p

A vector; see details

ip

A vector of probabilities corresponding to indep(p) where p is vector with unit sum

include.Jacobian

Boolean, with default TRUE meaning to include the Jacobian transformation in the evaluation, and FALSE meaning to ignore it; use FALSE for likelihood work and TRUE for probability densities

give

Boolean, with default FALSE meaning to return the value of the integral and TRUE meaning to return the full output of adaptIntegrate()

normalize

Boolean, indicates whether return value of mean_hyper2() is normalized to have unit sum

...

Further arguments passed to adaptIntegrate()

Details

  • Function B() returns the normalizing constant of a hyperdirichlet likelihood function. Internally, \(p\) is converted to e (by e_to_p()) and the integral proceeds over a hypercube. This function can be very slow, especially if disallowed is used.

  • Function dhyper2(ip,H) is a probability density function on the independent components of a unit-sum vector, that is, ip=indep(p). This function calls B() each time so might be a performance bottleneck.

  • Function probability() gives the probability of an observation from a hyperdirichlet distribution satisfying !disallowed(p).

  • Function mgf() is the moment generating function, taking an argument that specifies the powers of p needed: the expectation of \(\prod_{i=1}^n {p_i}^{{\rm powers}[i]}\) is returned.

  • Function mean_hyper2() returns the mean value of the hyperdirichlet distribution. This is computationally slow (consider maxp() for a measure of central tendency). The function takes a normalize argument, not passed to adaptIntegrate(): this is Boolean with FALSE meaning to return the value found by integration directly, and default TRUE meaning to normalize so the sum is exactly 1

Value

  • Function B() returns a scalar: the normalization constant

  • Function dhyper2() is a probability density function over indep(p)

  • Function mean() returns a \(k\)-tuple with unit sum

  • Function mgf() returns a scalar equal to the expectation of p^power

  • Functions is.proper() and validated() return a Boolean

  • Function probability() returns a scalar, a (Bayesian) probability

Author

Robin K. S. Hankin

Note

The adapt package is no longer available on CRAN; from 1.4-3, the package uses adaptIntegrate of the cubature package.

See also

loglik

Examples


# Two different measures of central tendency:
# mean_hyper2(chess,tol=0.1)   # takes ~10s to run
maxp(chess)                    # faster
#>   Topalov     Anand    Karpov 
#> 0.4036108 0.3405168 0.2558723 

# Using the 'disallowed' argument typically results in slow run times;
# use high tol for speed:

# probability(chess,disallowed=function(p){p[1]>p[2]},tol=0.5)
# probability(chess,disallowed=function(p){p[1]<p[2]},tol=0.5)

# Above should sum to 1 [they are exclusive and exhaustive events]