Estimate the scale parameter and shape parameters of the Mathai and Moschopoulos (1992) multivariate gamma distribution by maximum likelihood estimation.

gammaff.mm(lscale = "loglink", lshape = "loglink",
           iscale = NULL, ishape = NULL, imethod = 1,
           eq.shapes = FALSE, sh.byrow = TRUE, zero = "shape")

Arguments

lscale, lshape

Link functions applied to the (positive) parameters \(b\), and \(s_1\), ..., \(s_Q\) respectively. See Links for more choices. In the future, lshapes might be used instead; of course, this link applies to all the shape parameters.

iscale, ishape, sh.byrow

Optional initial values. The default is to compute them internally. Argument sh.byrow is fed into byrow in matrix and concerns the ordering of the initial shape parameters; a matrix of dimension \(n\) by \(Q\) is ultimately constructed. See also CommonVGAMffArguments.

eq.shapes

Logical. Constrain the shape parameters to be equal? See also CommonVGAMffArguments.

imethod, zero

See CommonVGAMffArguments.

Details

This distribution has the bivariate gamma distribution bigamma.mckay as a special case. Let \(Q > 1\) be the number of columns of the response matrix y. Then the joint probability density function is given by $$f(y_1,\ldots,y_Q; b, s_1, \ldots, s_Q) = y_1^{s_1} (y_2 - y_1)^{s_2} \cdots (y_Q - y_{Q-1})^{s_Q} \exp(-y_Q / b) / [b^{s_Q^*} \Gamma(s_1) \cdots \Gamma(s_Q)]$$ for \(b > 0\), \(s_1 > 0\), ..., \(s_Q > 0\) and \(0<y_1< y_2<\cdots<y_Q<\infty\). Also, \(s_Q^* = s_1+\cdots+s_Q\). Here, \(\Gamma\) is the gamma function, By default, the linear/additive predictors are \(\eta_1=\log(b)\), \(\eta_2=\log(s_1)\), ..., \(\eta_M=\log(s_Q)\). Hence \(Q = M - 1\). The marginal distributions are gamma, with shape parameters \(s_1\) up to \(s_Q\), but they have a common scale parameter \(b\).

The fitted value returned is a matrix with columns equalling their respective means; for column \(j\) it is sum(shape[1:j]) * scale.

The correlations are always positive; for columns \(j\) and \(k\) with \(j < k\), the correlation is sqrt(sum(shape[1:j]) /sum(shape[1:k])). Hence the variance of column \(j\) is sum(shape[1:j]) * scale^2.

Value

An object of class "vglmff" (see vglmff-class). The object is used by modelling functions such as vglm and vgam.

References

Mathai, A. M. and Moschopoulos, P. G. (1992). A form of multivariate gamma distribution. Ann. Inst. Statist. Math., 44, 97–106.

Author

T. W. Yee

Note

The response must be a matrix with at least two columns. Apart from the first column, the differences between a column and its LHS adjacent column must all be positive. That is, each row must be strictly increasing.

See also

bigamma.mckay, gammaff.

Examples

if (FALSE) { # \dontrun{
data("mbflood", package = "VGAMdata")
mbflood <- transform(mbflood, VdivD = V / D)
fit <- vglm(cbind(Q, y2 = Q + VdivD) ~ 1,
            gammaff.mm, trace = TRUE, data = mbflood)
coef(fit, matrix = TRUE)
Coef(fit)
vcov(fit)
colMeans(depvar(fit))  # Check moments
head(fitted(fit), 1)
} # }