A Non-Iterative Alternative to Ordinal Log-Linear Models

Log-linear modeling is a popular statistical tool for analysing a contingency table. This presentation focuses on an alternative approach to modeling ordinal categorical data. The technique, based on orthogonal polynomials, provides a much simpler method of model fitting than the conventional approach of maximum likelihood estimation, as it does not require iterative calculations nor the fitting and re-fitting to search for the best model. Another advantage is that quadratic and higher order effects can readily be included, in contrast to conventional log-linear models which incorporate linear terms only. The focus of the discussion is the application of the new parameter estimation technique to multi-way contingency tables with at least one ordered variable. This will also be done by considering singly and doubly ordered two-way contingency tables. It will be shown by example that the resulting parameter estimates are numerically similar to corresponding maximum likelihood estimates for ordinal log-linear models.


Introduction
Log-linear modeling is a popular statistical tool for the analysis of a contingency table and is especially popular in English speaking countries; see Gower [18] and van der Heijden and de Leeuw [27].From the analysis of data using log-linear models, the researcher can determine important associations that exist in the data.The conventional method of estimating parameters from a log-linear model has been to use the maximum likelihood estimation technique, such as the Newton-Raphson iterative procedure or iterative proportional scaling (see Agresti [1]).The advantage of the technique lies in the ability to identify multi-level associations in multi-way contingency tables.However, a problem with this approach, as completely ordered tables are considered.Doubly ordered tables are also investigated where only the row and column variables are ordered.Similar results can be obtained for any combination of two ordered variables.Finally a singly ordered three-way contingency table is considered, where the ordered categories are associated with the row responses.Similar results can be obtained for the ordered column or ordered tube categories.The alternative procedure to parameter estimation is made in Sections 6 and 7 which are applicable to two-way and three-way contingency tables respectively and uses the association models described in Sections 4 and 5.The calculation of the orthogonal polynomials used in this paper require a set of scores that reflect the ordinal ordinal structure of the ordered variable.In many cases natural scores, i.e. 1, 2, . . ., are the easiest set of scores to use and give results that compare very well with those calculated using more computationally intensive scores.
A comparison of results using the technique discussed in this paper with those obtained using the conventional log-linear analysis approach is made in Section 8. Two examples are presented -a doubly ordered two-way contingency table and a completely ordered three-way table.It is shown through these examples that the non-iterative parameter estimation approach is an easier and more informative procedure to implement than the conventional maximum likelihood approach.

Notation
Consider an I × J two-way contingency table, N , where the (i, j) th cell entry is given by n ij for i = 1, 2, . . ., I and j = 1, 2, . . ., J. It is assumed that such a data set is randomly selected from a population.As the focus of this paper is to consider an analysis of ordinal categorical data, we will be considering situations when N has ordered rows and/or ordered columns.Let the grand total of N be n and the probability matrix be P so that the (i, j) th cell entry is p ij = n ij /n and For the analysis of multi-level data, suppose we have a three-way contingency table, N , consisting of I rows, J columns and, using the terminology of Kroonenberg [22], K tubes.Kendall and Stuart [21] used the term layer.
Let the grand total of N be n where the (i, j, k) th element is denoted as n ijk .We denote the same element of the matrix P by p ijk = n ijk /n.The i th row, j th column and k th tube row marginal proportions are defined as The distinction between population (parameter) values and their estimated (sample) values is made in the conventional way by the placement of a ∧ over the estimated sample value, while the parameter has no ∧ included.For example, the parameter value µ is estimated by the sample value μ.

Orthogonal Polynomials
For the parameter estimation technique discussed in this paper, we make use of orthogonal polynomials based on the ordinal nature of the categories within a variable.The polynomials used as those considered by Beh [3,4], Beh and Davy [6,7], and Rayner and Best [24].
For ordered column categories of a two-way contingency table, say, let {b v (j)} be the set of orthogonal polynomials on {p •j }, for v = 1, 2, . . ., (J − 1) and j = 1, 2, . . ., J.These polynomials have the property and are defined so that where and Note that b −1 (j) = 0 and b 0 (j) = 1 for all j.For the calculation of (1), (2), (3) and (4), a set of scores is required to represent the ordinal nature of the categories.Scores may chosen subjectively, such as natural scores where s J (j) = j, or by more objective means, such as Nishisato scores (Nishisato and Arri [23]), or midrank scores.Beh [4] shows that often natural scores give very similar results to the more computationally intensive scores.For the assignement of scores to ordinal log-linear models, Ishii-Kuntz [20] views the use of natural scores as very simple and are therefore an advantage of the procedure.Row orthogonal polynomials can be calculated in a similar manner to the column polynomials.For ordered row categories, let {a u (i)} denote the set of orthogonal polynomials on {p i• }, for u = 1, 2, . . ., (I − 1) and i = 1, 2, . . ., I.These polynomials have the property When a third variable is included, consisting of K tubes, the set of tube orthogonal polynomials of order w, {c w (k)}, are defined so that they are orthogonal across the tube marginal probabilities {p ••k } subject to the constraint

Two-Way Contingency Tables
The RC canonical correlation model can be used to measure the departure from independence of two ordered variables from a two-way crossclassification and has been extensively investigated.See for example Gilula [14], Gilula and Haberman [15,16], and Gilula and Ritov [17].In this section we consider the correlation models for two-way ordinal cross -classifications as seen in Beh [3], Beh and Davy [6,7], Rayner and Best [24] and Best and Rayner [8].

Doubly ordered two-way contingency tables
The RC bivariate correlation model for a doubly ordinal two-way crossclassifications has been discussed by Rayner and Best [24] and Beh [3].These authors demonstrate that p ij can be reconstituted by the saturated model where The √ n θuv values are asymptotically standard normally distributed random variables and are asymptotically independent, and are calculated using the set of row scores s I (i) and the set of column scores s J (j).To determine significant association values we can compare each of these terms with what is expected under the standard normal distribution.Rayner and Best [25] showed that when natural scores are used, θ11 is the Pearson product moment correlation.They also showed that when midrank scores are used to calculate this quantity, then it is equivalent to Spearman's rank correlation.
When model ( 5) is unsaturated, it can sometimes give negative p ij values.In this case, Rayner and Best [24] suggested the following alternative exponential model It is with this and similar exponential approximations that we establish a parallel between our non-iterative approach and conventional log-linear models.

Singly ordered two-way contingency tables
Suppose now we consider a singly ordered two-way contingency table, where the columns are of an ordinal nature.The saturated bivariate correlation model discussed in Beh [3] and Rayner and Best [25] is where The exponential approximation of model ( 7) is ORDINAL LOG-LINEAR MODELS 73

Three-Way Contingency Tables
For the analysis of three-way ordinal contingency tables, Beh and Davy [6,7] presented some models which are extensions of the bivariate models presented in the previous section.These trivariate models can be used to estimate the parameter values from a log-linear model of the data.

Completely ordered three-way tables
Suppose that a three-way contingency table consists of three variables which have an ordered categorical responses.Beh and Davy [6] show that the saturated RC model can be defined by where Beh and Davy [6] showed that θuv0 , θu0w and θ0vw are bivariate rowcolumn, row-tube and column-tube correlation values respectively.The θuvw are extensions of the Pearson product moment correlation when natural scores are used for trivariate data.When midrank scores are used, they are extensions to the Spearman rank correlation for the same data structure.
The term θuvw is a measure of the deviation of the rows, columns and tubes up to the (u, v, w) th trivariate moment in the data from might be expected under complete independence.For example, θ111 is the linear-bylinear-by-linear association and can be used to assess the trivariate location of the three variables.Refer to Davy, Rayner and Beh [10] for information on the interpretation of θuvw .
The unsaturated form of model of ( 9) can sometimes calculate negative values in much the same way as the bivariate models can.To avoid this, the value of p ijk can be approximately reconstituted by the exponential model

Doubly ordered three-way tables
For a three-way contingency table with ordered row and column categories, Beh and Davy [7] showed that the saturated model associated with such a data structure is where The value of θuv(k) is the (u, v) th bivariate association between the rows and columns at the k th tube.For example, θ11(2) is the linear-by-linear association between the rows and column at the second tube category.Beh and Davy [7] defined the (u, v) th bivariate association for the ordered rows and ordered columns of the three-way contingency table as The value of the (i, j, k) th cell proportion can be approximately reconstituted using the exponential form of (11) ORDINAL LOG-LINEAR MODELS 75 5.3.Singly ordered three-way tables For the analysis of a three-way table with ordered row categories only (singly ordered), Beh and Davy [7] used where The value of the (i, j) th cell proportion can be approximately reconstituted using the exponential model

Two-Way Log-linear Models
Log-linear models are a popular way to determine if there is any relationship between any sets of variables of a cross-classification.In this section we look at a new method of calculating the parameters of a log-linear model and show that it is a far more versatile approach than the conventional maximum likelihood approach discussed in Agresti [1, p80] and Fienberg [13, p54].

Doubly ordered two-way contingency tables
Suppose that we wish to conduct a log-linear analysis on a two-way contingency table with ordered rows and ordered columns.Agresti [1] showed that the most applicable theoretical log-linear model is where µ I and µ J are the mean row and column scores respectively and The value of m ij is the expected value of the (i, j) th cell entry in N .
The value of φ in (17) describes the association between the two ordinal categorical variables.Note that, just as described by Agresti [1, p77], when φ = 0, this leads to the independence model, while the quantity φ [s I (i) − µ I ] [s J (j) − µ J ] measures the deviation from independence for ln m ij .Haberman [19] proposed a method of estimating this parameter value using maximum-likelihood estimation, which Agresti [1] describes.Now, consider the exponential form of the RC canonical correlation model given by (6).Multiplying this model by n and taking the natural logarithm of both sides yields  18) is the ordinal form of the log-linear model presented in van der Heijden et al. [28].

Linear Terms
If only the linear terms of the row and column categories are considered, that is when u = v = 1, then (18) simplifies to where μI = I i=1 s I (i) p i• and σ2

Quadratic Terms
Suppose now we wish to include the quadratic terms with the linear terms in (18) so that we may determine the quadratic nature of the variables and ORDINAL LOG-LINEAR MODELS 77 be more precise with our estimate of φ.Then the conventional log-linear model is Therefore, considering the linear and quadratic terms in (18), we obtain, After a little algebra, the parameters φ, φ 12 and φ 21 can be estimated by where A r 2 and A c 2 are the values of ( 4) for the row and column categories respectively and γI and γJ denote the (standardised) skewness coefficients for the row and column scores respectively.
Observing the above equations, if the quadratic terms do not contribute to the association in the contingency table, then (21) reduces to (20).
Alternatively, the orthogonal polynomials are constructed to include a (standardised) linear, quadratic and higher order vector and could be included instead.In this case, the log-linear model with the quadratic terms included becomes If the researcher was interested in expressing the model in the conventional form, the row (and column) orthogonal polynomials can be regressed such that where A ug are unknown constants.For example, a 2 (i) can be expressed as Thus, we can model the saturated version of the log-linear model by considering all terms in (18).But for the sake of approximation, it is practical to only consider the first few terms; the linear and quadratic terms, say.The investigation of quadratic terms will be left for future work.

Singly ordered two-way contingency tables
The approximation of parameters for singly ordered two-way contingency tables can be made in much the same way as discussed in the previous section.
Suppose that we have a two-way table with ordered column categories.The log-linear model presented in Agresti [1] which is applicable to this table is ln where The value of τ J(i) , such that I i=1 τ J(i) = 0, can be calculated using maximum likelihood estimation, and is the measure of linearity if the ordered column categories on the j th row category.This quantity is shown to be related to θ(i)1 as seen in (8).
Consider the exponential two-way model of (7).Multiplying it by n and taking the natural logarithm of both sides of the expression yields Consider the linear component of the ordered column categories.Then (23) reduces to Comparing (24) with the classical log-linear model of ( 22), shows that the parameter τ J(i) can be approximated by Three-Way Log-Linear Models Beh and Davy [6,7] described a method of parameter estimation for threeway contingency tables with any number of the variables ordered.In this section, we provide some of the derivations in their papers and extend them.

Completely ordered three-way contingency tables
For a three-way contingency table, where all three variables are ordered, Fienberg [13] and Agresti [1] offer the following log-linear model: where γ k = 0.The value of s K (k) is the score associated with the k th ordered tube category.The values φ IJ , φ IK and φ JK describe the bivariate association for each pair of variables.These values correspond to the linear-by-linear associations for each pair of variables, while φ IJK is the association term and corresponds to the linear-by-linearby-linear association.To show this, consider the exponential model (10).
Taking the natural logarithm of both sides after multiplying by n yields where If only the linear components for each variable are considered then (26) simplifies to ln n ijk = μ + αi + βj + γk + θ110 [s Comparing ( 25) and ( 26) using natural scores, the parameters from the log-linear model of ( 25 In model (28), τ IK(k) is the effect the linearity of the row categories has on the k th tube profile while τ JK(k) is the effect of the linearity of the columns on the k th tube category.These values correspond to θ10(k) and θ01(k) respectively from (12).The value of φ IJ(k) is related to the linear-bylinear association of the ordered rows and columns at the k th non-ordered tube category.This value differs slightly from the parameter in the model presented in Agresti [1], φ IJ , which is the overall linear-by-linear association value.To verify these results, consider the exponential model of association (14).Taking the natural logarithm of this model after multiplying by n yields ln n ijk = μ + αi + βj + γk + Consider the linear terms in model (29).Then the linear model is Therefore, comparing models ( 28) and (30), the parameters τ IK(k) , τ JK(k) and φ IJ(k) can be approximated by The overall linear-by-linear association between the rows and column of the doubly ordered three-way contingency table, using ( 13) and (31), is while the parameter φ IJ from the model of Agresti [1] and Fienberg [13] can be approximated by γ k = 0.The value τ IJ(j) is a measure of the linearity of the rows and its effect on the j th column category.Similarly,  20) can be made in terms of the fitted values.A goodness-of-fit test is performed using only the above first order parameter estimates by calculating the chi-squared statistic for each of the three methods.Doing so yields the goodness-of-fit chi-squared statistic using the Haberman and Agresti estimates, X 2 HAB = 11.51938 and X 2 AGR = 11.44213respectively.Using the estimate obtained from the non-iterative procedure outlined here, X 2 N I = 11.40662.The p-values associated with these statistics, and (I − 1) (J − 1) − 1 = 14 degrees of freedom are 0.645, 0.651 and 0.654 respectively.This suggests that the non-iterative approach for Table 1 produces estimates which are slightly better than those obtained from the iterative methods.This single parameter also explains nearly three quarters of the total variation from independence for each of the three estimates, since the Pearson chi-squared statistic is X 2 = 45.98526.

Example 2 -A completely ordered three-way table
Consider Table 2 which is a 4 × 5 × 3 contingency table classifying 1517 people according to the number years of completed schooling, their number of siblings and their general level of happiness.Clogg [9] and Beh and Davy [6,7] considered this contingency table for their analyses.The three variables can be seen to be of an ordinal nature, and so the log-linear analysis of sub-section 7.1 can be applied.Beh and Davy [6,7] showed that the the three bivariate relationships were highly significant and that the trivariate was not.Therefore, we only need to concern ourselves with estimating the parameters φ IJ , φ IK and φ JK .Table 3 gives a comparison of each of the parameters calculated using maximum likelihood estimation, denoted by φ, and those estimates obtained using (27), denoted by φ.

p
ij = 1.Define the i th row marginal proportion as p i• = J j=1 p ij and define the j th column marginal probability as p •j = I i=1 p ij so that I i=1 − σJ γJ φ12 − σI γI φ21(21) .3.Singly ordered three-way contingency tables Consider now a three-way contingency table with ordered row categories.The theoretical log-linear model for such a table is ln m ijk

Table 1 .
Cross-classification of 1660 Patients from Mid-town Manhattan According to Mental Health Status and Parental Socio-Economic Status.

Table 2 .
Cross-classification of 1517 People According to their Happiness, Schooling and Number of Siblings.

Table 3 .
MLE and Alternative Values of Parameters for Table2.