CESTI
Constructs an equivalent completely testable multivariate general linear hypothesis H BU = G from a partially testable hypothesis HpBU = Gp.
Required Arguments
HP — NHP by NCOEF matrix Hp with each row corresponding to a row in the hypothesis and containing the constants that specify a linear combination of the regression coefficients. (Input)
NDEP — Number of dependent (response) variables. (Input)
NU — U matrix option. (Input)
For positive NU, NU is the number of linear combinations of the dependent variables to be considered. If NU = 0, the hypothesis is HpB = Gp, and U is automatically taken to be the identity. NU must be less than or equal to NDEP .
GP — Matrix Gp containing the null hypothesis values. (Input)
If NU = 0, then GP is NHP by NDEP; otherwise, GP is NHP by NU.
R — NCOEF by NCOEF upper triangular matrix containing the R matrix. (Input)
The R matrix can come from a regression fit based on a QR decomposition of the matrix of regressors or based on a Cholesky factorization RTR of the matrix of sums of squares and crossproducts of the regressors. Elements to the right of a diagonal element of R that is zero must also be zero. A zero row indicates a nonfull rank model. For an R matrix that comes from a regression fit with linear equality restrictions on the parameters, each row of R corresponding to a restriction must have a corresponding diagonal element that is negative. The remaining rows of R must have positive diagonal elements. Only the upper triangle of R is referenced.
IRANKP — Rank of Hp. (Output)
NH — Number of rows in the completely testable hypothesis (also, the degrees of freedom for the hypothesis). (Output)
The degrees of freedom for the hypothesis (NH) classify the hypothesis Hp BU = Gp as nontestable (NH = 0), partially testable (0 < NH < IRANKP), or completely testable (0 < NH = IRANKP).
H — NH by NCOEF matrix H with each row corresponding to a row in the completely testable hypothesis and containing the constants that specify an estimable linear combination of the regression coefficients. (Output)
If HP is not needed, H and HP can occupy the same storage locations.
G — Matrix G containing the null hypothesis values for the completely testable hypothesis. (Output)
If NU = 0, then G is NH by NDEP, otherwise, G is NH by NU. If GP is not needed, G and GP can occupy the same storage locations.
Optional Arguments
NHP — Number of rows in the hypothesis. (Input)
Default: NHP = size (HP,1).
NCOEF — Number of regression coefficients in the model. (Input)
Default: NCOEF = size (HP,2).
LDHP — Leading dimension of HP exactly as specified in the dimension statement of the calling program. (Input)
Default: LDHP = size (HP,1).
LDGP — Leading dimension of GP exactly as specified in the dimension statement in the calling program. (Input)
Default: LDGP = size (GP,1).
LDR — Leading dimension of R exactly as specified in the dimension statement in the calling program. (Input)
Default: LDR = size (R,1).
LDH — Leading dimension of H exactly as specified in the dimension statement of the calling program. (Input)
Default: LDH = size (H,1).
LDG — Leading dimension of G exactly as specified in the dimension statement in the calling program. (Input)
Default: LDG = size (G,1).
FORTRAN 90 Interface
Generic: CALL CESTI (HP, NDEP, NU, GP, R, IRANKP, NH, H, G [, …])
Specific: The specific interface names are S_CESTI and D_CESTI.
FORTRAN 77 Interface
Single: CALL CESTI (NHP, NCOEF, HP, LDHP, NDEP, NU, GP, LDGP, R, LDR, IRANKP, NH, H, LDH, G, LDG)
Double: The double precision name is DCESTI.
Description
Once a general linear model y = X β + ɛ is fitted, particular hypothesis tests are frequently of interest. If the matrix of regressors X is not full rank (as evidenced by the fact that some diagonal elements of the R matrix output from the fit are equal to zero), methods that use the results of the fitted model to compute the hypothesis sum of squares (see routine RHPSS) require one to specify in the hypothesis only linear combinations of the regression parameters that are estimable. A linear combination of regression parameters cT β is estimable means that there exists some vector a such that cT = aTX, i.e., cT is in the space spanned by the rows of X. For a further discussion of estimable functions, see Maindonald (1984, pages 166‑168) and Searle (1971, pages 180 ‑ 188). Routine CESTI is only useful in the case of nonfull rank regression models, i.e., when the problem of estimability arises.
Peixoto (1986) noted that the customary definition of testable hypothesis in the context of a general linear hypothesis test H β = g is overly restrictive. He extended the notion of a testable hypothesis (a hypothesis composed of estimable functions of the regression parameters) to include partially testable and completely testable hypotheses. A hypothesis H β = g is partially testable means that the intersection of the row space of H (denoted by R(H)) and the row space of X(R(X)) is not essentially empty and is a proper subset of R(H), i.e., {0} ⊂ R(H) ∩ R(X) ⊂ R(H). A hypothesis H β = g is completely testable means that {0} ⊂ R(H) ⊆ R(X). Peixoto also demonstrated a method for converting a partially testable hypothesis to one that is completely testable so that the usual method for obtaining the sum of squares for the hypothesis from the results of the fitted model can be used. The method replaces Hp in the partially testable hypothesis Hp β = gp by a matrix H whose rows are a basis for the intersection of the row space of Hp and the row space of X. A corresponding conversion of the null hypothesis values from gp to g is also made. A sum of squares for the completely testable hypothesis can then be computed (see routine RHPSS). The sum of squares that is computed for the hypothesis H β = g equals the difference in the error sums of squares from two fitted models the restricted model with the partially testable hypothesis Hp β = gp adjoined to the model as linear equality restrictions (see routine RLEQU) and the unrestricted model.
Routines RGLM, RGIVN, RLEQU, and RCOV can be used to compute the fit of the general linear model prior to invoking CESTI. The R matrix is required for input to CESTI. After converting a partially testable hypothesis to a completely testable hypothesis, RHPSS can be invoked to compute the sum of squares for the hypothesis.
For the general case of the Multivariate General Linear Model Y = XB + E with possible linear equality restrictions on the regression parameters, CESTI converts the partially testable hypothesis Hp BU = Gp to a completely testable hypothesis H BU = G. For the case of the linear model with linear equality restrictions, the definitions of estimable functions, nontestable hypotheses, partially testable hypotheses, and completely testable hypothesis are similar to those previously given for the unrestricted model with the exception that R(X) is replaced by R(R) where R is the upper triangular matrix output from RLEQU. The nonzero rows of R form a basis for the rowspace of the matrix (XT, AT)T. The rows of H form an orthonormal basis for the intersection of two subspaces: the subspace spanned by the rows of Hp and the subspace spanned by the rows of R. The algorithm used by CESTI for computing the intersection of these two subspaces is based on an algorithm for computing angles between linear subspaces due to to Bjorck and Golub (1973). (See also Golub and Van Loan 1983, pages 429‑430). The method is closely related to a canonical correlation analysis discussed by Kennedy and Gentle (1980, 56‑565). The algorithm is as follows:
1. Compute a QR factorization of
with column permutations so that
Here, P1 is the associated permutation matrix that is also an orthogonal matrix. Determine the rank of Hp as the number of nonzero diagonal elements of R1, say n1. Partition Q1 = (Q11, Q12) so that Q11is the first n1columns of Q1. Set IRANKP = n1.
2. Compute a QR factorization of the transpose of the R matrix input to CESTI with column permutations so that
Determine the rank of R from the number of nonzero diagonal elements of R, say n2. Partition Q2 = (Q21, Q22) so that Q21 is the first n2 columns of Q2.
3. Form
4. Compute the singular values of A
and the left singular vectors W of the singular value decomposition of A so that
If σ1 < 1, then the dimension of the intersection of the two subspaces is s = 0. Otherwise, take the dimension of the intersection to be s if σs = 1 > σs+1. Set NH = s.
5. Let W1 be the first s columns of W. Set H = (Q1W1)T.
6. Take R11 to be a NHP by NHP matrix related to R1 as follows. If NHP ≤ NCOEF, R11 equals the first NHP rows of R1. Otherwise, R11 contains R1 in its first NCOEF rows and zeros in the remaining rows. Compute a solution Z to the linear system
using routine GIRTS (IMSL MATH/LIBRARY). If this linear system is declared inconsistent, an error message with error code equal to 2 is issued.
7. Partition
so that Z1 is the first n1 rows of Z. Set
The degrees of freedom (NH) classify the hypothesis Hp BU = Gp as nontestable (NH = 0), partially testable (0 < NH < IRANKP), or completely testable (0 < NH = IRANKP).
For further details concerning the algorithm, see Sallas and Lionti (1988).
Comments
1. Workspace may be explicitly provided, if desired, by use of C2STI/DC2STI. The reference is:
CALL C2STI (NCOEF, NHP, HP, LDHP, NDEP, NU, GP, LDGP, R, LDR, IRANKP, NH, H, LDH, G, LDG, IWK, WK)
The additional arguments are as follows:
IWK — Work vector of length max{NHP, NCOEF}.
WK — Work vector of length NCOEF * m + NCOEF2 + NHP2 + n * r + n2 + m + max{2 * m, n + r + max(n, r) ‑ 1}.
2. Informational errors
Type |
Code |
Description |
4 |
1 |
There is inadequate space to store the completely testable hypothesis. Increase LDH or LDG so that it is greater than or equal to NH. |
3 |
2 |
The hypothesis Hp BU = Gp is inconsistent. |
Example
A one-way analysis-of-variance model discussed by Peixoto (1986) is fitted to some data. The model is
yij= μ + α i + ɛ ij (i, j) = (1, 1), (2, 1), (2, 2)
The model is fitted using routine RGLM. Next, the partially testable hypothesis
is converted to a completely testable hypothesis using CESTI. Sum of squares associated with the hypothesis are computed using routine RHPSS. Finally, the F statistic is computed along with the associated p‑value using routine FDF (see Chapter 17, "Probability Distribution Functions and Inverses").
USE IMSL_LIBRARIES
IMPLICIT NONE
INTEGER LDB, LDG, LDGP, LDH, LDHP, LDR, LDSCPE, LDSCPH, &
LDX, LINDEF, MAXB, NCOL, NDEP, NEF,NHP, NROW, MAXCL, &
NCLVAR, J, NU
PARAMETER (LINDEF=1, MAXB=3, MAXCL=2, NCLVAR=1, NCOL=2, &
NDEP=1, NEF=1, NHP=2, NROW=3, LDB=MAXB, LDG=NHP, &
LDGP=NHP, LDH=NHP, LDHP=NHP, LDR=MAXB, LDSCPE=NDEP, &
LDSCPH=NDEP, LDX=NROW)
!
INTEGER INDCL(NCLVAR), INDDEP(NDEP),INDEF(LINDEF), INTCEP, &
IRANK, IRANKP, IRBEF(NEF+1),NCOEF, NH, NOUT, NVEF(NEF)
REAL B(LDB,NDEP), DFE, DFH, F, G(LDG,NDEP), GP(LDGP,NDEP),&
H(LDH,MAXB), HP(LDHP,MAXB), PVALUE, R(LDR,MAXB), &
SCPE(LDSCPE,NDEP), SCPH(LDSCPH,NDEP),X(LDX,NCOL)
!
DATA X/1.0, 2.0, 2.0, 17.3, 24.1, 26.3/
DATA INDCL/1/, NVEF/1/, INDEF/1/, INDDEP/2/
DATA (HP(1,J),J=1,MAXB)/0.0, 1.0, 0.0/
DATA (HP(2,J),J=1,MAXB)/0.0, 0.0, 1.0/
DATA GP/5.0, 3.0/
!
CALL RGLM (X, INDCL, NVEF, INDEF, NDEP, INDDEP, MAXCL, B, &
IRBEF=IRBEF, R=R, DFE=DFE, SCPE=SCPE)
NCOEF = IRBEF(NEF+1) - 1
!
NU = 0
CALL CESTI (HP, NDEP, NU, GP, R, IRANKP, NH, H, G, NCOEF=NCOEF)
!
CALL UMACH (2, NOUT)
IF (NH .EQ. 0) THEN
WRITE (NOUT,*) 'Nontestable hypothesis'
ELSE IF (NH .LT. IRANKP) THEN
WRITE (NOUT,*) 'Partially testable hypothesis'
ELSE
WRITE (NOUT,*) 'Completely testable hypothesis'
END IF
CALL WRRRN ('H', H, NH, NCOEF, LDH)
CALL WRRRN ('G', G, NH, NDEP, LDG)
CALL RHPSS (H, B, G, R, SCPH, DFH=DFH)
!
F = (SCPH(1,1)/DFH)/(SCPE(1,1)/DFE)
PVALUE = 1.0 - FDF(F,DFH,DFE)
WRITE (NOUT,*)
WRITE (NOUT,*) 'Degrees of Sum of Prob. of'
WRITE (NOUT,*) ' Freedom Squares F-statistic Larger F'
WRITE (NOUT,99999) DFH, SCPH(1,1), F, PVALUE
99999 FORMAT (F8.1, 3X, 1F10.3, F11.3, 2X, F10.4)
END
Output
Partially testable hypothesis
H
1 2 3
0.0000 0.7071 -0.7071
G
1.414
Degrees of Sum of Prob. of
Freedom Squares F-statistic Larger F
1.0 65.340 27.000 0.1210