CTWLS
Performs a generalized linear least‑squares analysis of transformed probabilities in a two‑dimensional contingency table.
Required Arguments
TABLE — NRESP by NPOP matrix containing the frequency count in each cell of each population. (Input)
The i‑th column of TABLE contains the counts for the i‑th population.
ITRAN — Vector of length NTRAN containing the transformation code for each of the NTRAN transformations to be applied. (Input)
ITRAN is not referenced and can be a vector of length 1 in the calling program if NTRAN = 0. Let a “response” denote a transformed cell probability. Then, ITRAN(1) contains the first transformation to be applied to the cell probabilities, ITRAN(2) contains the second transformation, which is to be applied to the responses obtained after ITRAN(1) is performed, etc. Note that the k‑th transformation takes the ISIZE(k ‑ 1) responses at step k into ISIZE(k) responses, where ISIZE(0) is taken to be NPOP * NRESP. Let y denote the vector result of a transformation, x denote the responses before the transformation is applied, A denote a matrix of constants, and v denote a vector of constants. Then, the possible transformations are
ITRAN |
Transformation |
1 |
Linear, defined over all populations (y = Ax) |
2 |
Logarithmic (y(i, j) = ln(x(i, j)) |
3 |
Exponential (y(i, j) = exp(x(i, j))) |
4 |
Additive (y(i, j) = y(i, j) + v(i, j)) |
5 |
Linear, defined for one population and, identically, applied over all populations (y(i) = Ax(i)) |
where y(i) and x(i) are the subvectors for the i‑th population, y(i, j) and x(i, j) denote the j‑th response in the i‑th population, and v(i, j) denotes the corresponding element of the vector “v”. Transformation type 5 is the same as transformation type 1 when the same linear transformation is applied in each population (i.e., the type 1 matrix is block diagonal with identical blocks). Because the size of the type 5 transformation matrix A is NPOP2 times smaller than the type 1 transformation matrix, the type 5 transformation is usually preferred where it can be used.
ISIZE — Vector of length NTRAN containing the number of response functions defined by the k‑th transformation. (Input)
Transformation types 2, 3, and 4 have the same number of output responses as are input, and elements of ISIZE corresponding to transformations of these types should reflect this fact. Transformation types 1 and 5 can either increase or, more commonly, decrease the number of responses. For transformation type 5, if m linear transformations are defined for each population, the corresponding element of ISIZE should be m * NPOP.
AMATS — Vector containing the transformation constants. (Input)
AMATS contains the transformation matrices and vectors needed in the type 1, 4 and 5 transformations. While AMATS is a vector, its elements may be treated as a number of matrices or vectors where the number of structures depends upon the transformation types as follows:
ITRAN |
Type |
Dimension |
Length |
1 |
Matrix |
m by n |
M * n |
2, 3 |
Not referenced |
|
0 |
4 |
Vector |
M |
M |
5 |
Matrix |
m/NPOP by n/NPOP |
M * n/(NPOP * NPOP) |
Here, m = ISIZE(i) and n = ISIZE(i ‑ 1), and ISIZE(0) is not input (in ISIZE) but is taken to be NPOP * NRESP. Matrices and vectors are stored consecutively in AMATS with column elements for matrices stored consecutively as is standard in FORTRAN. Thus, if ITRAN(1) = 5 and ITRAN(2) = 4, with NREP = 3, NPOP= 2, and ISIZE(1) = ISIZE(2) = 2, then the vector AMATS would contain in consecutive positions A(1, 1), A(2, 1), A(1, 2), A(2, 2), A(1, 3), A(2, 3), v(1), v(2), v(3), v(4) where A is the matrix for transformation type 5 and v is the vector for transformation type 4.
X — Design matrix of size ISIZE(NTRAN) by NCOEF. (Input, if NCOEF > 0)
X contains the design matrix for predicting the transformed cell probabilities F from the covariates stored in X. If NCOEF = 0, X is not referenced and can be a 1 by 1 matrix in the calling program.
NH — Vector of length NUMH. (Input, if NCOEF > 0)
NH(i) contains the number of consecutive rows in H used to specify hypothesis i. If NCOEF = 0, NH is not referenced and can be a vector of length 1 in the calling program.
H — Matrix of size m by NCOEF containing the constants to be used in the multivariate hypothesis tests. (Input, if NCOEF > 0)
Here, m is the sum of the elements in NH. Each hypothesis is of the form
H0 : C * COEF = 0, where C for the i‑th hypothesis is NH(i) rows of H, and COEF is estimated in the linear model. The first NH(1) rows of H make up the first hypothesis, the next NH(2) rows make up the second hypothesis, etc. If NCOEF = 0, H is not referenced and can be a 1 by 1 matrix in the calling program.
CHSQ — NUMH + 1 by 3 matrix containing the results of the hypothesis tests. (Output, if NCOEF > 0)
The first row of CHSQ contains the results for test 1, the next row contains the results for test 2, etc. The last row of CHSQ contains a test of the adequacy of the model. Within each row, the first column contains the chi‑squared statistic, the second column contains its degrees of freedom, and the last column contains the probability of a larger chi‑squared. If NCOEF = 0, CHSQ is not referenced and can be a 1 by 1 matrix in the calling program.
COEF — NCOEF by 4 matrix containing the coefficient estimates and related statistics. (Output, if NCOEF > 0)
The columns of coefficient are as follows:
Col. |
Statistic |
1 |
Coefficient estimate |
2 |
Estimated standard error of the coefficient |
3 |
z‑statistic for a test that the coefficient equals 0 versus the Two‑sided alternative |
4 |
p‑value corresponding to the z‑statistic |
If NCOEF = 0, COEF is not referenced and can be a 1 by 1 matrix in the calling program.
COVCF — NCOEF by NCOEF matrix containing the estimated variances and covariances of COEF. (Output, if NCOEF > 0)
If NCOEF = 0, COVCF is not referenced and can be a 1 by 1 matrix in the calling program.
F — Vector of length ISIZE(NTRAN) containing the transformed probabilities, the responses. (Output)
COVF — Matrix of size ISIZE(NTRAN) by ISIZE(NTRAN) containing the estimated variances and covariances of F. (Output)
RESID — ISIZE(NTRAN) by 4 matrix containing a case analysis for the transformed probabilities as estimated by the linear model. (Output, if NCOEF > 0)
The linear model gives F = X * BETA. The columns of RESID are as follows:
Col. |
Description |
1 |
Residual |
2 |
Standard error |
3 |
Leverage |
4 |
Standardized residual |
If NCOEF = 0, RESID is not referenced and can be a 1 by 1 matrix in the calling program.
Optional Arguments
NRESP — Number of cells in each population. (Input)
Default: NRESP = size (TABLE,1).
NPOP — Number of populations. (Input)
Default: NPOP = size (TABLE,2).
LDTABL — Leading dimension of TABLE exactly as specified in the dimension statement in the calling program. (Input)
Default: LDTABL = size (TABLE,1).
NTRAN — Number of transformations to be applied to the cell probabilities. (Input)
Cell probabilities are computed as the frequency count for the cell divided by the population sample size. Set NTRAN = 0 if a linear model predicting the cell probabilities is to be used.
Default: NTRAN = size (ITRAN,1).
NCOEF — Number of coefficients in the linear model relating the transformed probabilities F to the design matrix X. (Input)
Let F denote the vector result of applying the NTRAN transformations, and assume that the model gives F = X * COEF. Then, NCOEF is the length of COEF.
Default: NCOEF = size (X,2).
LDX — Leading dimension of X exactly as specified in the dimension statement in the calling program. (Input)
Default: LDX = size (X,1).
NUMH — Number of multivariate hypotheses to be tested on the coefficients in COEF. (Input, if NCOEF > 0)
If NCOEF = 0, NUMH is not referenced.
Default: NUMH = size (NH,1).
LDH — Leading dimension of H exactly as specified in the dimension statement in the calling program. (Input)
Default: LDH = size(H,1).
IPRINT — Printing option. (Input)
Default: IPRINT = 0.
IPRINT |
Action |
0 |
No printing is performed. |
1 |
Print all output arrays and vectors. |
2 |
Print all output arrays and vectors as well as the matrices and vectors in AMATS. |
LDCHSQ — Leading dimension of CHSQ exactly as specified in the dimension statement in the calling program. (Input)
Default: LDCHSQ = size (CHSQ,1).
LDCOEF — Leading dimension of COEF exactly as specified in the dimension statement in the calling program. (Input)
Default: LDCOEF = size (COEF,1).
LDCOVC — Leading dimension of COVCF exactly as specified in the dimension statement in the calling program. (Input)
Default: LDCOVC = size (COVCF,1)
LDCOVF — Leading dimension of COVF exactly as specified in the dimension statement in the calling program. (Input)
Default: LDCOVF = size (COVF,1).
LDRESI — Leading dimension of RESID exactly as specified in the dimension statement in the calling program. (Input)
Default: LDRESI = size (RESID,1).
FORTRAN 90 Interface
Generic: CALL CTWLS (TABLE, ITRAN, ISIZE, AMATS, X, NH, H, CHSQ, COEF, COVCF, F, COVF, RESID [, …])
Specific: The specific interface names are S_CTWLS and D_CTWLS.
FORTRAN 77 Interface
Single: CALL CTWLS (NRESP, NPOP, TABLE, LDTABL, NTRAN, ITRAN, ISIZE, AMATS, NCOEF, X, LDX, NUMH, NH, H, LDH, IPRINT, CHSQ, LDCHSQ, COEF, LDCOEF, COVCF, LDCOVC, F, COVF, LDCOVF, RESID, LDRESI)
Double: The double precision name is DCTWLS.
Description
Routine CTWLS performs weighted least‑squares analysis of a general p = NPOP population by r = NRESP response categories per population contingency table. After division by the sample size, there are n = pr cell probabilities.
Define s = ISIZE(NTRAN) responses fi such that each response is obtained from the cell probabilities as fi = gi(p1, p2, …, pn), for i = 1, …, s. Call the functions gi the response functions. Then, if
is the asymptotic covariance matrix of the responses, and X is a design matrix for a linear model predicting f = X β with q = NCOEF coefficients β = COEF, then CTWLS performs a weighted least‑squares analysis of the model f = X β where the generalized weights are given by
Estimates obtained in this way are best asymptotic normal estimates of β.
Let
denote the estimated variance‑covariance matrix of the estimated cell probabilities, and let (∂gi/∂pj) denote the matrix of partial derivatives of gi with respect to pj. Then,
is given by
where the (i, j)-th element in
is computed as
pi(δij − pj)
Here, δij = 1 if i = j and is zero otherwise.
In CTWLS, the transformations gi are defined by successive application of one of five types of simpler transformations. Let pi = h0,j for j = 1, …, n denote the n cell probabilities, and let hi,j denote the ISIZE(i) responses obtained after i simple transformations have been performed with hi denoting the corresponding vector of estimates. Then, the simple transformations are defined by:
1. Linear: hi+1 = Aihi where Ai is a matrix of coefficients specified via the vector AMATS in CTWLS.
2. Logarithmic: hi+1,j = ln(hi,j) where j = 1, …, ISIZE(i). That is, take the logarithm of each of the responses.
3. Exponential: hi+1,j = exp(hi,j) where j = 1, …, ISIZE(i). That is, take the exponential of each of the responses.
4. Additive: hi+1,j = hi,j + vj , where j = 1, …, ISIZE(i), and vj is specified via the vector AMATS in CTWLS. Additive transformations are generally used to adjust for zero cells or to apply a continuity correction to the cell probabilities.
5. Linear (by population):
is the vector of responses at stage i in the j‑th population, and Ai is a matrix of coefficients specified via AMATS.
Given the responses fi and their covariances
estimates for β are computed via generalized least squares as
Let Σ β denote the asymptotic covariance matrix of β. Then, Σ β is estimated by
Hypothesis tests of the form Ho : Ci β = 0 are performed when requested. Here, Ci is a matrix of coefficients specified via a submatrix of the matrix H. Results are returned in the vector CHSQ. The asymptotic chi‑squared test for testing the null hypothesis is given by
This test has qi = rank(Ci) degrees of freedom. If zero degrees of freedom are returned, the hypothesis cannot be tested in the original parameterization.
A test of the model checks that the residuals obtained from the model f = X β are not too large. This test, which has s ‑ q degrees of freedom, is an asymptotic chi‑squared test and is computed as
Residuals from the generalized linear model are easily computed as
where xi is the row of the design matrix X corresponding to the i‑th observation. This residual has the asymptotic variance
where (A)ii denotes the i‑th diagonal element of matrix A. A standardized residual is then computed as
which has an asymptotic standard normal distribution if the model is correct.
The leverage of observation i, vi, is computed as
It is a measure of the importance of the observation in the predicted values. Values greater than 2q/s are large.
Because the tests performed by CTWLS are asymptotic ones, the user should treat the results with caution. The reported asymptotic p‑values are most likely to be exact when the number of counts in each cell is large (say 5 or more), and less exact for smaller cell counts. Care should also be taken to avoid illegal operations. For example, the routine returns an error message when the log of a negative or zero value is attempted. When this occurs, the user should either use a continuity correction (i.e. modify the transformations used by adding a constant to all cells or to the cell resulting in the illegal operation) or abandon the model.
Comments
1. Workspace may be explicitly provided, if desired, by use of C2WLS/DC2WLS. The reference is:
CALL C2WLS (NRESP, NPOP, TABLE, LDTABL, NTRAN, ITRAN, ISIZE, AMATS, NCOEF, X, LDX, NUMH, NH, H, LDH, IPRINT, CHSQ, LDCHSQ, COEF, LDCOEF, COVCF, LDCOVC, F, COVF, LDCOVF, RESID, LDRESI, PDER, FRQ, EST, XX, WK, IWK, WWK)
The additional arguments are as follows:
PDER — Work vector of length ISIZE(NTRAN) * max(NPOP * NRESP, ISIZE(i)) if NTRAN is greater than zero. PDER is not used and can be dimensioned of length 1 if NTRAN = 0.
FRQ — Work vector of length NPOP.
EST — Work vector of length NPOP * NRESP + ISIZE(1) + … + ISIZE(NTRAN).
XX — Work vector of length (NCOEF + 1) * ISIZE(NTRAN) if NCOEF is greater than zero. If NCOEF = 0, XX is not referenced and can be a vector of length 1 in the calling program.
WK — Work vector of length 3(max(NPOP * NRESP, ISIZE(i))) + NCOEF + 1.
IWK — Work vector of length max(NH(i)) if NUMH is greater than 0. If NCOEF = 0, IWK is not referenced and can be a vector of length 1 in the calling program.
WWK — Work vector of length max(NH(i)) * (4 + NCOEF + max(NCOEF, max(NH(i))) if NUMH is greater than 0. If NUMH = 0, WWK is not referenced and can be a vector of length 1 in the calling program.
2. Informational error
Type |
Code |
Description |
4 |
1 |
A negative response occurred while performing a logarithmic transformation. The logarithm of a negative number is not allowed. |
Examples
Example 1
This example is taken from Landis, Stanish, Freeman, and Koch (1976), pages 213‑217. Generalized kappa statistics are computed via vector functions of the form:
F(p) = exp(A4 ln(A3 exp(A2 ln(A1p))))
where p is the cell probabilities. The raw frequencies are given as two 4 × 4 contingency tables. These tables are reorganized as a single 16 × 2 table for input into CTWLS. The input tables are
Two generalized kappa statistics using two different sets of weights are computed for each population. Hypothesis tests are then performed on the four resulting generalized kappa statistics. In this example, the matrix of covariates is an identity matrix so that tests on the responses are performed.
USE CTWLS_INT
IMPLICIT NONE
INTEGER IPRINT, LDCHSQ, LDCOEF, LDCOVC, LDCOVF, LDH, LDRESI, &
LDTABL, LDX, NCOEF, NPOP, NTRAN
PARAMETER (IPRINT=2, LDCHSQ=10, LDCOEF=4, LDCOVC=4, LDCOVF=4, &
LDH=10, LDRESI=4, LDTABL=16, LDX=4, NCOEF=4, NPOP=2, &
NTRAN=8)
!
INTEGER ISIZE(NTRAN), ITRAN(NTRAN), NH(9)
REAL A1(10,16), A2(18,10), A3(4,18), A4(2,4), AMATS(420), &
CHSQ(LDCHSQ,3), COEF(LDCOEF,4), COVCF(LDCOVC,NCOEF), &
COVF(LDCOVF,LDCOVF), F(LDX), H(LDH,4), &
RESID(LDRESI,4), TABLE(LDTABL,NPOP), X(LDX,NCOEF)
!
EQUIVALENCE (A1, AMATS(1)), (A2, AMATS(161)), (A3, AMATS(341)), &
(A4, AMATS(413))
!
DATA TABLE/38, 5, 0, 1, 33, 11, 3, 0, 10, 14, 5, 6, 3, 7, 3, 10, &
5, 3, 0, 0, 3, 11, 4, 0, 2, 13, 3, 4, 1, 2, 4, 14/
DATA X/1, 0, 0, 0, 0, 1, 0, 0, 0, 0, 1, 0, 0, 0, 0, 1/
DATA NH/1, 1, 1, 1, 1, 1, 2, 1, 1/
DATA H/1, 0, 1, 0, 0, 0, 1, 0, 1, 0, 0, 1, -1, 0, 0, 0, 0, 1, 0, &
1, 0, 0, 0, 1, 0, 1, -1, 0, -1, 0, 0, 0, 0, 0, 1, -1, 0, &
-1, 0, -1/
DATA ITRAN/5, 2, 5, 3, 5, 2, 5, 3/
DATA ISIZE/20, 20, 36, 36, 8, 8, 4, 4/
DATA A1/1, 0, 0, 0, 1, 0, 0, 0, 1, 1, 1, 0, 0, 0, 0, 1, 0, 0, 0, &
.5, 1, 0, 0, 0, 0, 0, 1, 0, 0, .25, 1, 0, 0, 0, 0, 0, 0, 1,&
0, 0, 0, 1, 0, 0, 1, 0, 0, 0, 0, .5, 0, 1, 0, 0, 0, 1, 0, &
0, 1, 1, 0, 1, 0, 0, 0, 0, 1, 0, 0, .5, 0, 1, 0, 0, 0, 0, &
0, 1, 0, .25, 0, 0, 1, 0, 1, 0, 0, 0, 0, .25, 0, 0, 1, 0, &
0, 1, 0, 0, 0, .5, 0, 0, 1, 0, 0, 0, 1, 0, 1, 1, 0, 0, 1, &
0, 0, 0, 0, 1, 0, .5, 0, 0, 0, 1, 1, 0, 0, 0, 0, 0, 0, 0, &
0, 1, 0, 1, 0, 0, 0, .25, 0, 0, 0, 1, 0, 0, 1, 0, 0, .5, 0, &
0, 0, 1, 0, 0, 0, 1, 1, 1/
DATA A2/1, 1, 1, 1, 0, 0, 0, 0, 0, 0, 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, 0, &
0, 0, 0, 0, 0, 1, 1, 1, 1, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, &
0, 0, 0, 0, 0, 0, 0, 1, 1, 1, 1, 0, 0, 1, 0, 0, 0, 1, 0, 0, &
0, 1, 0, 0, 0, 1, 0, 0, 0, 0, 0, 0, 1, 0, 0, 0, 1, 0, 0, 0, &
1, 0, 0, 0, 1, 0, 0, 0, 0, 0, 0, 1, 0, 0, 0, 1, 0, 0, 0, 1, &
0, 0, 0, 1, 0, 0, 0, 0, 0, 0, 1, 0, 0, 0, 1, 0, 0, 0, 1, 0, &
0, 0, 1, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, &
0, 1, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, &
1/
DATA A3/-1, -1, 0, 0, 0, -.5, 1, .5, 0, -.25, 1, .75, 0, 0, 1, &
1, 0, -.5, 1, .5, -1, -1, 0, 0, 0, -.5, 1, .5, 0, -.25, 1, &
.75, 0, -.25, 1, .75, 0, -.5, 1, .5, -1, -1, 0, 0, 0, -.5, &
1, .5, 0, 0, 1, 1, 0, -.25, 1, .75, 0, -.5, 1, .5, -1, -1, &
0, 0, 1, 0, 0, 0, 0, 1, 0, 0/
DATA A4/1, 0, 0, 1, -1, 0, 0, -1/
!
CALL CTWLS (TABLE, ITRAN, ISIZE, AMATS, X, NH, H, &
CHSQ, COEF, COVCF, F, COVF, RESID, IPRINT=IPRINT)
!
END
Output
Hypothesis Tests on Coefficients
H-1
1 0 0 0
H-2
0 1 0 0
H-3
1 -1 0 0
H-4
0 0 1 0
H-5
0 0 0 1
H-6
0 0 1 -1
H-7
1 0 -1 0
0 1 0 -1
H-8
1 0 -1 0
H-9
0 1 0 -1
Hypothesis Chi-Squared Statistics
Degrees of
Hypothesis Chi-Squared freedom p-value
1 16.99 1 0.0000
2 39.70 1 0.0000
3 39.54 1 0.0000
4 14.27 1 0.0002
5 30.07 1 0.0000
6 28.76 1 0.0000
7 1.07 2 0.5850
8 0.90 1 0.3425
9 1.06 1 0.3040
Degrees of
Chi-Squared freedom p-value
Model Test 0.00 0 NaN
Coefficient Statistics
Coefficient Standard Error Statistic p-value
1 0.2079 0.05 4.12 0.0000
2 0.3150 0.05 6.30 0.0000
3 0.2965 0.08 3.78 0.0002
4 0.4069 0.07 5.48 0.0000
Asymptotic Coefficient Covariance
1 2 3 4
1 2.5457E-03 2.3774E-03 0. 0.
2 2.4988E-03 0. 0.
3 6.1629E-03 5.6229E-03
4 5.5069E-03
Residual Analysis
Standard Standardized
Residual Error Leverage Residual
1 0.0000 0.0000 1.0000 NaN
2 0.0000 0.0000 1.0000 NaN
3 0.0000 0.0000 1.0000 NaN
4 0.0000 0.0000 1.0000 NaN
Transformed Probabilities
1 0.2079
2 0.3150
3 0.2965
4 0.4069
Asymptotic Covariance of the Transformed Probabilities
1 2 3 4
1 2.5457E-03 2.3774E-03 0. 0.
2 2.4988E-03 0. 0.
3 6.1629E-03 5.6229E-03
4 5.5069E-03
Linear transformation matrix, by population, for transformation 5
1 2 3 4 5 6 7 8 9
1 1.000 1.000 1.000 1.000 0.000 0.000 0.000 0.000 0.000
2 0.000 0.000 0.000 0.000 1.000 1.000 1.000 1.000 0.000
3 0.000 0.000 0.000 0.000 0.000 0.000 0.000 0.000 1.000
4 0.000 0.000 0.000 0.000 0.000 0.000 0.000 0.000 0.000
5 1.000 0.000 0.000 0.000 1.000 0.000 0.000 0.000 1.000
6 0.000 1.000 0.000 0.000 0.000 1.000 0.000 0.000 0.000
7 0.000 0.000 1.000 0.000 0.000 0.000 1.000 0.000 0.000
8 0.000 0.000 0.000 1.000 0.000 0.000 0.000 1.000 0.000
9 1.000 0.000 0.000 0.000 0.000 1.000 0.000 0.000 0.000
10 1.000 0.500 0.250 0.000 0.500 1.000 0.500 0.250 0.250
10 11 12 13 14 15 16
1 0.000 0.000 0.000 0.000 0.000 0.000 0.000
2 0.000 0.000 0.000 0.000 0.000 0.000 0.000
3 1.000 1.000 1.000 0.000 0.000 0.000 0.000
4 0.000 0.000 0.000 1.000 1.000 1.000 1.000
5 0.000 0.000 0.000 1.000 0.000 0.000 0.000
6 1.000 0.000 0.000 0.000 1.000 0.000 0.000
7 0.000 1.000 0.000 0.000 0.000 1.000 0.000
8 0.000 0.000 1.000 0.000 0.000 0.000 1.000
9 0.000 1.000 0.000 0.000 0.000 0.000 1.000
10 0.500 1.000 0.500 0.000 0.250 0.500 1.000
Linear transformation matrix, by population, for transformation 5
1 2 3 4 5 6 7 8 9
1 1.000 0.000 0.000 0.000 1.000 0.000 0.000 0.000 0.000
2 1.000 0.000 0.000 0.000 0.000 1.000 0.000 0.000 0.000
3 1.000 0.000 0.000 0.000 0.000 0.000 1.000 0.000 0.000
4 1.000 0.000 0.000 0.000 0.000 0.000 0.000 1.000 0.000
5 0.000 1.000 0.000 0.000 1.000 0.000 0.000 0.000 0.000
6 0.000 1.000 0.000 0.000 0.000 1.000 0.000 0.000 0.000
7 0.000 1.000 0.000 0.000 0.000 0.000 1.000 0.000 0.000
8 0.000 1.000 0.000 0.000 0.000 0.000 0.000 1.000 0.000
9 0.000 0.000 1.000 0.000 1.000 0.000 0.000 0.000 0.000
10 0.000 0.000 1.000 0.000 0.000 1.000 0.000 0.000 0.000
11 0.000 0.000 1.000 0.000 0.000 0.000 1.000 0.000 0.000
12 0.000 0.000 1.000 0.000 0.000 0.000 0.000 1.000 0.000
13 0.000 0.000 0.000 1.000 1.000 0.000 0.000 0.000 0.000
14 0.000 0.000 0.000 1.000 0.000 1.000 0.000 0.000 0.000
15 0.000 0.000 0.000 1.000 0.000 0.000 1.000 0.000 0.000
16 0.000 0.000 0.000 1.000 0.000 0.000 0.000 1.000 0.000
17 0.000 0.000 0.000 0.000 0.000 0.000 0.000 0.000 1.000
18 0.000 0.000 0.000 0.000 0.000 0.000 0.000 0.000 0.000
10
1 0.000
2 0.000
3 0.000
4 0.000
5 0.000
6 0.000
7 0.000
8 0.000
9 0.000
10 0.000
11 0.000
12 0.000
13 0.000
14 0.000
15 0.000
16 0.000
17 0.000
18 1.000
Linear transformation matrix, by population, for transformation 5
1 2 3 4 5 6 7 8 9
1 -1.000 0.000 0.000 0.000 0.000 -1.000 0.000 0.000 0.000
2 -1.000 -0.500 -0.250 0.000 -0.500 -1.000 -0.500 -0.250 -0.250
3 0.000 1.000 1.000 1.000 1.000 0.000 1.000 1.000 1.000
4 0.000 0.500 0.750 1.000 0.500 0.000 0.500 0.750 0.750
10 11 12 13 14 15 16 17 18
1 0.000 -1.000 0.000 0.000 0.000 0.000 -1.000 1.000 0.000
2 -0.500 -1.000 -0.500 0.000 -0.250 -0.500 -1.000 0.000 1.000
3 1.000 0.000 1.000 1.000 1.000 1.000 0.000 0.000 0.000
4 0.500 0.000 0.500 1.000 0.750 0.500 0.000 0.000 0.000
Linear transformation matrix, by population, for transformation 5
1 2 3 4
1 1.000 0.000 -1.000 0.000
2 0.000 1.000 0.000 -1.000
Example 2
The second example is taken from Prentice (1976) and involves a logistic fit to the mortality of beetles after exposure to various concentrations of carbon disulphide. Because one of the cells on input has a count of zero and it is not possible to take the logarithm of zero, a constant 0.5 is added to each cell prior to calling CTWLS. The model can be expressed as
where i indexes the 8 populations. The data is given as:
X |
fi1 |
fi2 |
1.690 |
6 |
53 |
1.724 |
13 |
47 |
1.755 |
18 |
44 |
1.784 |
28 |
28 |
1.811 |
52 |
11 |
1.836 |
53 |
6 |
1.861 |
61 |
1 |
1.883 |
60 |
0 |
For comparison, a maximum fit yields
(see routine CTGLM).
USE CTWLS_INT
IMPLICIT NONE
INTEGER IPRINT, LDCHSQ, LDCOEF, LDCOVC, LDCOVF, LDH, LDRESI, &
LDTABL, LDX, NCOEF, NPOP, NRESP, NTRAN, NUMH
PARAMETER (IPRINT=2, LDCOVF=8, LDH=1, LDX=8, NCOEF=2, NPOP=8, &
NRESP=2, NTRAN=2, NUMH=0, LDCHSQ=NUMH+1, &
LDCOEF=NCOEF, LDCOVC=NCOEF, LDRESI=LDX, LDTABL=NRESP)
!
INTEGER ISIZE(NTRAN), ITRAN(NTRAN), NH(1)
REAL AMATS(2), CHSQ(LDCHSQ,3), COEF(LDCOEF,4), &
COVCF(LDCOVC,NCOEF), COVF(LDCOVF,LDCOVF), F(LDX), &
H(LDH,4), RESID(LDRESI,4), TABLE(LDTABL,NPOP), &
X(LDX,NCOEF)
!
DATA TABLE/6, 53, 13, 47, 18, 44, 28, 28, 52, 11, 53, 6, 61, 1, &
60, 0/, ITRAN/2, 5/, ISIZE/16, 8/, AMATS/1, -1/
DATA X/8*1, 1.690, 1.724, 1.755, 1.784, 1.811, 1.836, 1.861, &
1.883/
!
TABLE = TABLE + 0.5
!
CALL CTWLS (TABLE, ITRAN, ISIZE, AMATS, X, NH, H, &
CHSQ, COEF, COVCF, F, COVF, RESID, NUMH=NUMH, &
IPRINT=IPRINT)
!
END
Output
Test of the Model
Degrees of
Chi-Squared freedom p-value
8.43 6 0.2081
Coefficient Statistics
Coefficient Standard Error Statistic p-value
1 -55.6590 5.02 -11.10 0.0000
2 31.4177 2.83 11.09 0.0000
Asymptotic Coefficient Covariance
1 2
1 25.16 -14.20
2 8.024
Residual Analysis
Standard Standardized
Residual Error Leverage Residual
1 0.4552 0.3232 0.6052 1.4086
2 0.2368 0.2480 0.6468 0.9548
3 -0.3568 0.2413 0.7608 -1.4787
4 -0.3902 0.2285 0.7440 -1.7076
5 0.2800 0.2761 0.7192 1.0141
6 0.0840 0.3484 0.7036 0.2410
7 0.9042 0.7749 0.8791 1.1670
8 1.2953 1.3777 0.9413 0.9402
Transformed Probabilities
1 -2.108
2 -1.258
3 -0.878
4 0.000
5 1.518
6 2.108
7 3.714
8 4.796
Asymptotic Covariance of the Transformed Probabilities
1 2 3 4 5
1 0.1725 0. 0. 0. 0.
2 9.5127E-02 0. 0. 0.
3 7.6526E-02 0. 0.
4 7.0175E-02 0.
5 0.1060
6 7 8
1 0. 0. 0.
2 0. 0. 0.
3 0. 0. 0.
4 0. 0. 0.
5 0. 0. 0.
6 0.1725 0. 0.
7 0.6829 0.
8 2.017
Linear transformation matrix, by population, for transformation 5
1 2
1.000 -1.000