Wishart {mixAK} | R Documentation

## Wishart distribution

### Description

Wishart distribution

Wishart(nu, S),

where nu are degrees of freedom of the Wishart distribution and S is its scale matrix. The same parametrization as in Gelman (2004) is assumed, that is, if W~Wishart(nu,S) then

E(W) = nu*S.

### Usage

    dWishart(W, df, S, log=FALSE)

    rWishart(n, df, S)

### Arguments

`W` | Either a matrix with the same number of rows and columns as `S` (1 point sampled from the Wishart distribution) or a matrix with `ncol` equal to `ncol`*`(ncol+1)`/2 and `n` rows (`n` points sampled from the Wishart distribution for which only lower triangles are given in rows of the matrix `W`).
`n` | number of observations to be sampled.
`df` | degrees of freedom of the Wishart distribution.
`S` | scale matrix of the Wishart distribution.
`log` | logical; if `TRUE`, log-density is computed

### Details

The density of the Wishart distribution is the following

f(W) = (2^{nu*p/2} * pi^{p*(p-1)/4} * prod[i=1]^p Gamma((nu + 1 - i)/2))^{-1} * |S|^{-nu/2} * |W|^{(nu - p - 1)/2} * exp(-0.5*tr(S^{-1}*W)),

where p is number of rows and columns of the matrix W.

In the univariate case, Wishart(nu,S) is the same as Gamma(nu/2, 1/(2*S)).

Generation of random numbers is performed by the algorithm described in Ripley (1987, pp. 99).

### Value

Some objects.

### Value for dWishart

A numeric vector with evaluated (log-)density.

### Value for rWishart

If `n` equals 1 then a sampled symmetric matrix W is returned.

If `n` > 1 then a matrix with sampled points (lower triangles of W) in rows is returned.

### Author(s)

Arnošt Komárek arnost.komarek[AT]mff.cuni.cz

### References

Gelman, A., Carlin, J. B., Stern, H. S., and Rubin, D. B. (2004). Bayesian Data Analysis, Second edition. Boca Raton: Chapman and Hall/CRC.

Ripley, B. D. (1987). Stochastic Simulation. New York: John Wiley and Sons.

### Examples

    set.seed(1977)
    ### The same as gamma(shape=df/2, rate=1/(2*S))
    df <- 1
    S  <- 3

    w <- rWishart(n=1000, df=df, S=S)
    mean(w)    ## should be close to df*S
    var(w)     ## should be close to 2*df*S^2

    dWishart(w[1], df=df, S=S)
    dWishart(w[1], df=df, S=S, log=TRUE)

    dens.w <- dWishart(w, df=df, S=S)
    dens.wG <- dgamma(w, shape=df/2, rate=1/(2*S))
    rbind(dens.w[1:10], dens.wG[1:10])

    ldens.w <- dWishart(w, df=df, S=S, log=TRUE)
    ldens.wG <- dgamma(w, shape=df/2, rate=1/(2*S), log=TRUE)
    rbind(ldens.w[1:10], ldens.wG[1:10])

    ### Bivariate Wishart
    df <- 2
    S <- matrix(c(1,3,3,13), nrow=2)

    print(w2a <- rWishart(n=1, df=df, S=S))
    dWishart(w2a, df=df, S=S)

    w2 <- rWishart(n=1000, df=df, S=S)
    print(w2[1:10,])
    apply(w2, 2, mean)                ## should be close to df*S
    (df*S)[lower.tri(S, diag=TRUE)]

    dens.w2 <- dWishart(w2, df=df, S=S)
    ldens.w2 <- dWishart(w2, df=df, S=S, log=TRUE)
    cbind(w2[1:10,], data.frame(Density=dens.w2[1:10], Log.Density=ldens.w2[1:10]))

    ### Trivariate Wishart
    df <- 3.5
    S <- matrix(c(1,2,3,2,20,26,3,26,70), nrow=3)

    print(w3a <- rWishart(n=1, df=df, S=S))
    dWishart(w3a, df=df, S=S)

    w3 <- rWishart(n=1000, df=df, S=S)
    print(w3[1:10,])
    apply(w3, 2, mean)                ## should be close to df*S
    (df*S)[lower.tri(S, diag=TRUE)]

    dens.w3 <- dWishart(w3, df=df, S=S)
    ldens.w3 <- dWishart(w3, df=df, S=S, log=TRUE)
    cbind(w3[1:10,], data.frame(Density=dens.w3[1:10], Log.Density=ldens.w3[1:10]))

[Package mixAK version 0.6 Index]