# freund61: Freund's (1961) Bivariate Extension of the Exponential... In VGAM: Vector Generalized Linear and Additive Models

## Description

Estimate the four parameters of the Freund (1961) bivariate extension of the exponential distribution by maximum likelihood estimation.

## Usage

 ```1 2 3``` ```freund61(la = "loglink", lap = "loglink", lb = "loglink", lbp = "loglink", ia = NULL, iap = NULL, ib = NULL, ibp = NULL, independent = FALSE, zero = NULL) ```

## Arguments

 `la, lap, lb, lbp` Link functions applied to the (positive) parameters alpha, alpha', beta and beta', respectively (the “`p`” stands for “prime”). See `Links` for more choices. `ia, iap, ib, ibp` Initial value for the four parameters respectively. The default is to estimate them all internally. `independent` Logical. If `TRUE` then the parameters are constrained to satisfy alpha=alpha' and beta=beta', which implies that y1 and y2 are independent and each have an ordinary exponential distribution. `zero` A vector specifying which linear/additive predictors are modelled as intercepts only. The values can be from the set {1,2,3,4}. The default is none of them. See `CommonVGAMffArguments` for more information.

## Details

This model represents one type of bivariate extension of the exponential distribution that is applicable to certain problems, in particular, to two-component systems which can function if one of the components has failed. For example, engine failures in two-engine planes, paired organs such as peoples' eyes, ears and kidneys. Suppose y1 and y2 are random variables representing the lifetimes of two components A and B in a two component system. The dependence between y1 and y2 is essentially such that the failure of the B component changes the parameter of the exponential life distribution of the A component from alpha to alpha', while the failure of the A component changes the parameter of the exponential life distribution of the B component from beta to beta'.

The joint probability density function is given by

f(y1,y2) = alpha * beta' * exp(-beta' * y2 - (alpha+beta-beta') * y1)

for 0 < y1 < y2, and

f(y1,y2) = beta * alpha' * exp(-alpha' * y1 - (alpha+beta-alpha') * y2)

for 0 < y2 < y1. Here, all four parameters are positive, as well as the responses y1 and y2. Under this model, the probability that component A is the first to fail is alpha/(alpha+beta). The time to the first failure is distributed as an exponential distribution with rate alpha+beta. Furthermore, the distribution of the time from first failure to failure of the other component is a mixture of Exponential(alpha') and Exponential(beta') with proportions beta/(alpha+beta) and alpha/(alpha+beta) respectively.

The marginal distributions are, in general, not exponential. By default, the linear/additive predictors are eta1=log(alpha), eta2=log(alpha'), eta3=log(beta), eta4=log(beta').

A special case is when alpha=alpha' and beta'=beta', which means that y1 and y2 are independent, and both have an ordinary exponential distribution with means 1/alpha and 1/beta respectively.

Fisher scoring is used, and the initial values correspond to the MLEs of an intercept model. Consequently, convergence may take only one iteration.

## Value

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

## Note

To estimate all four parameters, it is necessary to have some data where y1<y2 and y2<y1.

The response must be a two-column matrix, with columns y1 and y2. Currently, the fitted value is a matrix with two columns; the first column has values (alpha'+beta)/(alpha' * (alpha+beta)) for the mean of y1, while the second column has values (beta'+alpha)/(beta' * (alpha+beta)) for the mean of y2. The variance of y1 is

[(alpha')^2 + 2 * alpha * beta + beta^2]/ [(alpha')^2 * (alpha + beta)^2],

the variance of y2 is

[(beta')^2 + 2 * alpha * beta + alpha^2]/ [(beta')^2 * (alpha + beta)^2],

the covariance of y1 and y2 is

[alpha' * beta' - alpha * beta]/ [alpha' * beta' * (alpha + beta)^2].

T. W. Yee

## References

Freund, J. E. (1961). A bivariate extension of the exponential distribution. Journal of the American Statistical Association, 56, 971–977.

`exponential`.
 ``` 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17``` ```fdata <- data.frame(y1 = rexp(nn <- 1000, rate = exp(1))) fdata <- transform(fdata, y2 = rexp(nn, rate = exp(2))) fit1 <- vglm(cbind(y1, y2) ~ 1, fam = freund61, data = fdata, trace = TRUE) coef(fit1, matrix = TRUE) Coef(fit1) vcov(fit1) head(fitted(fit1)) summary(fit1) # y1 and y2 are independent, so fit an independence model fit2 <- vglm(cbind(y1, y2) ~ 1, freund61(indep = TRUE), data = fdata, trace = TRUE) coef(fit2, matrix = TRUE) constraints(fit2) pchisq(2 * (logLik(fit1) - logLik(fit2)), # p-value df = df.residual(fit2) - df.residual(fit1), lower.tail = FALSE) lrtest(fit1, fit2) # Better alternative ```