Normalizing constant for the hyperdirichlet distribution

Description

Uses numerical techniques for calculating the normalizing constant for the hyperdirichlet distribution

Usage

1
2
3
4
5
6
7
8
B(x, ...)
NC(x)
calculate_B(x, disallowed=NULL, give=FALSE, ...)
probability(x, disallowed, ...)
mgf(x, powers, ...)
mean(x, ...)
is.proper(x,irregardless)
validated(x)

Arguments

x

Object of class “hyperdirichlet” (or coerced thereto)

powers

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

irregardless

Boolean; 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

give

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

...

Further arguments passed to adaptIntegrate()

Details

  • Function B() is the user-friendly version. It accesses the NC slot. If not NA, the value is returned; if NA, the normalizing constant is calculated using adaptIntegrate() of the cubature() package, via calculate_B().

  • Function NC() is not intended for the user. It is used internally as an accessor method for the NC slot, and this value is returned indiscriminately.

  • Function calculate_B() is the engine which actually does the work. Observe how p is converted to e (by e_to_p()) and the integral proceeds over a hypercube. Function dirichlet() and gd() do not use this as the normalizing constant has an analytical expression and this is used instead.

  • 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 p^powers is returned.

  • Function mean() returns the mean value of the hyperdirichlet distribution. This is computationally slow (consider maximum_likelihood() 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

  • Function is.proper() checks a hyperdirichlet distribution for being normalizable: a “proper” hyperdirichlet object has a finite integral and therefore can be normalized. This function is quite time-consuming for hyperdirichlet distributions of large dimension.

    The irregardless argument to function is.proper() is Boolean, with TRUE meaning to carry out the checks whatever the value of slot @validated [that is, validated(x)]. Default FALSE means that function is.proper() returns TRUE if @validated is TRUE and to carry out the check otherwise. Use this argument to force is.proper() to carry out a check even if not strictly necessary.

  • Function validated() is an accessor method for the @validated slot of hyperdirichlet object. It returns a Boolean variable with TRUE meaning that the object is known to be “proper” (ie is.proper(x) returns TRUE), so it is normalizable, even if the normalization constant is not known. This flag is present because many hyperdirichlet objects of interest are known a priori to be proper, so executing is.proper() would be unnecessary.

Value

Functions B(), NC(), calculate_NC() notionally return a scalar: the normalization constant

Functions mean() and mgf() return a k-tuple

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

Function probability() returns a scalar, a probability.

Note

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

Author(s)

Robin K. S. Hankin

See Also

hyperdirichlet

Examples

 1
 2
 3
 4
 5
 6
 7
 8
 9
10
a <- hyperdirichlet(c(4,3,6,5,4,3,2,1))
## Not run: 
B(a)                                    # Not recommended
a <- as.hyperdirichlet(a,TRUE)          # Recommended

is.proper(a)

mgf(a,powers=1:3)    # expectation of p1^1 * p2^2 * p3^3

## End(Not run)

Questions? Problems? Suggestions? or email at ian@mutexlabs.com.

Please suggest features or report bugs with the GitHub issue tracker.

All documentation is copyright its authors; we didn't write any of that.