check.HLP: HLP Link Status Check
In cta: Contingency Table Analysis Based on ML Fitting of MPH Models

Description Usage Arguments Details Value Author(s) References See Also Examples

Checks whether the link function L(\cdot) is a candidate HLP link function.

Specifically, this program checks whether L(\cdot) satisfies certain necessary conditions that follow from a sufficient condition for HLP link status.

If the necessary conditions are satisfied then there is corroborating evidence that L(\cdot) has HLP link status. If the necessary conditions are not satisfied, then the sufficient condition for HLP link status is not satisfied, so L(\cdot) may or may not have HLP link status.

1	check.HLP(L.fct, Z, tol = NULL)

`L.fct`	An R function object, indicating the link L(\cdot) for HLP link status check.
`Z`	Population (aka strata) matrix Z.
`tol`	The pre-set tolerance with which `norm(diff)` is to be compared with.

The main idea:

The model L(m) = Xβ is an HLP model if L(\cdot) is a smooth link function that satisfies the HLP conditions with respect to Z (i.e. strata s) and X. That is,

(1) L(\cdot) has HLP link status with respect to Z, and
(2) The implied constraint function h(m) = U'L(m) is Z homogeneous. Here, null(U') = span(X).

Here, (1) L(\cdot) has HLP link status with respect to Z if, for m = Diag(Zγ)p, equivalently, for γ = Z'm and p = Diag^{-1}(ZZ'm)m,

(1)(a) L(m) = a(γ) + L(p), where a(γ_{1}/γ_{2}) - a(1) = a(γ_{1}) - a(γ_{2}), i.e. a(γ) has the form C \logγ + \texttt{constant}; or
(1)(b) L(m) = G(γ) L(p), where G(γ) is a diagonal matrix with diagonal elements that are powers of the γ elements, i.e. L(\cdot) is Z homogeneous (see Lang (2004)); or
(1)(c) The components of L(\cdot) are a mixture of types (a) and (b): L_{j}(m) = a_{j}(γ) + L_{j}(p) or L_{j}(m) = G_{j}(γ) L_{j}(p), j = 1, …, l.

N.B. Lang (2005) defined HLP models as those satisfying (1)(a) and (2). mph.fit uses a broader definition of HLP model. Specifically, models satisfying (1)(b) and (2) or (1)(c) and (2) are also considered HLP models.

Conditions (1)(b) and (2) can be checked using the check.homog function. Condition (1)(c) is not checked.

This function, check.HLP, is concerned with sufficient condition (1)(a) only. If L(\cdot) satisfies (1)(a) then

(i) \texttt{diff1} = [L(Diag(Zγ_{1})p_{1}) - L(Diag(Zγ_{2})p_{1})] - [L(Diag(Z γ_{1}/γ_{2})p_{1}) - L(p_{1})] = 0, and
(ii) \texttt{diff2} = [L(Diag(Zγ_{1})p_{1}) - L(Diag(Zγ_{1})p_{2})] - [L(p_{1}) - L(p_{2})] = 0.

Here p_{i} = Diag^{-1}(ZZ'm_{i})m_{i}, where m_{i} = Diag(Zγ_{i})p_{i}, i = 1, 2.

This program randomly generates g1 (γ_{1}), g2 (γ_{2}), p1, p2, and computes norm(diff) = sqrt(norm(diff1)^2 + norm(diff2)^2). It returns a warning if norm(diff) is too far from 0.

check.HLP returns a character string chk. If chk = "", then there is corroborating evidence that L(\cdot) has HLP link status. If chk = paste("L(m) may not be an HLP link [based on tol=",tol,"]!"), then the sufficient condition for HLP link status is not satisfied, so L(\cdot) may or may not have HLP link status.

Joseph B. Lang

Lang, J. B. (2004) Multinomial-Poisson homogeneous models for contingency tables, Annals of Statistics, 32, 340–383.

Lang, J. B. (2005) Homogeneous linear predictor models for contingency tables, Journal of the American Statistical Association, 100, 121–134.

mph.fit, check.homog, check.zero.order.homog

# 3-by-3-by-3 Table.
# For a description of the model, see Michael Haber's Example 2,
# p. 433, in  Biometrics (in Shorter Communications), Vol. 42,
# No. 2. (Jun., 1986), pp. 429-435.
A <- gl(3, 9, 27)
B <- gl(3, 3, 27)
C <- gl(3, 1, 27)
MAB <- kronecker(diag(9), matrix(1, 1, 3))
MAC <- kronecker(diag(3), kronecker(matrix(1, 1, 3), diag(3)))
MBC <- kronecker(matrix(1, 1, 3), diag(9))
M <- rbind(MAB, MAC, MBC)
Mr <- M[-c(3, 6, 7, 8, 9, 12, 15, 16, 17, 18, 21, 24,
           25, 26, 27), ]
C <- c(1, -1, -1, 1, 0, 0, 0, 0, 0, 0, 0, 0,
       0, 0, 0, 0, 1, -1, -1, 1, 0, 0, 0, 0,
       0, 0, 0, 0, 0, 0, 0, 0, 1, -1, -1, 1)
C <- matrix(C, 3, 12, byrow = TRUE)
L.fct <- function(m) {
  p <- m / sum(m)
  C %*% log(Mr %*% p)
}
Z <- matrix(rep(1, 27), ncol = 1)
check.HLP(L.fct, Z)