Electronic Journal: Southwest Journal of Pure and Applied Mathematics Internet: http://rattler.cameron.edu/swjpam.html ISSN 1083-0464 Vol. No. 02, 1996, pp. 1–5. Submitted: February 3, 1996. Published: July 8, 1996.

PREDICTIVE DISTRIBUTION FOR POLYTOMOUS LOGISTIC MODELS

WITH NONLINEAR LINK FUNCTIONS

M.A. Tabatabai and K.P. Singh

Summary. This article is concerned with predictive distribution for polytomous logistic models with nonlinear link functions. Within a Bayesian framework, an approximation of the predictive distribution for determining probabilities is presented. We describe this method using nonlinear link functions. The generalized linear link function is simply a special case of nonlinear models. All of the generalized logistic models require nonlinear estimation procedures to evaluate maximum likelihood estimates.

(1991) A.M.S. Subject Classification Codes. 34.

Key Words and Phrases. polytomous logistic model,link function,predictive distribution

1. INTRODUCTION.

Statistical methods based on logistic models are routinely used in experiments with binary outcome data. It can also be used to model the association between a polytomous response variable and some explanatory variable of interest. This article’s emphasis is on obtaining a predictive distribution for future observation based on previous data. Prediction analysis plays the key role in the real decision-making process implicit in most medical survival studies. Discussions of prediction analysis can be found in Geisser (1971, 1980, 1982) and Aitchison and Dunsmore (1975). Data with polytomous response variables are frequently encountered in social and biomedical sciences. The response may be recorded on an ordinal or nominal scale and the link function can be linear or nonlinear functions. Here we use a Bayesian formulation to make predictions for future individuals. Given data

D = {(x_{i}, y_{i}), i = 1, …, n}

and x_{n+1} for future individual, we can predict the distribution of unobserved response y_{n+1}. The future response is polytomous response variable. This usage has previously been applied by Percy (1993) to find predictive distribution for linear models, including its application to survival data. Heitjan (1991) applied nonlinear modeling to serial immunologic data from a multiple sclerosis clinical trial. Rosner (1989) presented a polychotomous logistic model to analyze ophthalmologic data. We focus on polytomous logistic models with nonlinear links.

2. THE POLYTOMOUS LOGISTIC MODEL.

Let Y represent the response variable and assume Y ∈ {0, 1, …, k}, and let a sample of size n be taken. When the response is categorical, the following logit model may be used to relate the response

ln(P_{ij}/P_{io}) = g_{j}(x_{i}; β_{j}), j = 1, …, k,  (2.1)

Typeset by AMS-TEX

Cameron University, Lawton, Oklahoma

E-Mail Address: mohammat@cameron.edu

University of Alabama, Birmingham, Alabama

\begin{matrix}c \\ ◯ \end{matrix}1997 Cameron University

where g_{j}(x_{i}; β_{j}) is a function of a parameter vector β_{j} and covariate vector x_{i} of dimension P, and

P_{ij} = P (y_{i} = j|x_{i}), i = 1, …, n, j = 0, 1, …, k  (2.2)

the conditional probabilities P_{ij} of each outcome category given the covariate vector are

P_{ij} = exp[g_{j}(x_{i}; β_{j})]/ {∑_{l=0}^{k} exp[g_{l}(x_{i}; β_{l})]} , j = 0, 1, …, k, i = 1, …, n  (2.3)

where g_{0}(x_{i}; β_{0}) ≡ 0, ∑_{j=0}^{k} P_{ij} = 1, and β_{j} = (β_{j1}, β_{j2}, …, β_{jP_{j}})′. To construct the likelihood function, we

formulate k + 1 binary variables coded as zero or 1 to indicate group membership of an observation. The variables are coded as follows:

If y_{i} = j, then the coded vector y_{ij} = (y_{i0}, y_{i1}, … y_{ih}, … y_{ik}) has y_{ih} = 1 if h = j and 0 if h ≠ j.

3. THE PROPOSED TECHNIQUE.

Let y_{1}, y_{2}, …, y_{n} be n independent random variable with joint probability mass function

f(y_{1}, y_{2}, …, y_{n}| P_{1j}, P_{2j}, … P_{nj}, j = 0, 1, …, k) =
∏_{i=1}^{n} {(P_{i0})^{y_{io}}(P_{i1})^{y_{i1}} … (P_{ik})^{y_{ik}}}  (3.1)

where for each i, ∑_{j=0}^{k} y_{ij} = 1.

If we use (2.3), we obtain the joint distribution

f(y_{1}, y_{2}, …, y_{n}|x_{i}, i = 1, …, n; β_{j}, j = 1, …, k) =
∏_{i=1}^{n} {exp[∑_{j=0}^{k} y_{ij}g_{j}(x_{i}, β_{j})]/ ∑_{j=0}^{k} exp[g_{j}(x_{i}; β_{j})]}  (3.2)

One major application of linear as well as nonlinear polytomous regression is to make forecasts. With these models, one can incorporate recent data

{(y_{i}, x_{i}), i = 1, …, n} and x_{n+1} into

forecast for unobserved response y_{n+1}. To solve this problem and more importantly to estimate the predictive distribution we use Bayes theorem.

According to Bayes theorem, the predictive distribution

f(y_{n+1}|x_{1}, x_{2}, …, x_{n+1}; y_{1}, …, y_{n})
= ∫_{β} f(y_{n+1}|x_{n+1}; β) f(β|x_{1}, x_{2}, …, x_{n})dβ  (3.3)

where β = (β_{1}, β_{2}, …, β_{k}) and dim(β_{j}) = P_{j}, j = 1, …, k and ∑_{j=1}^{k} P_{j} < n.

To evaluate this integral we use an approximation similar to Percy’s (1993) for binary linear logistic regression. Before proceeding, we must decide on a suitable prior distribution for the parameter matrix β = (β_{1}, β_{2}, …, β_{k}). We assume β_{1}, β_{2}, …, β_{k} are k independent random vectors with each having a prior distribution as multivariate normal. This allows considerable flexibility for practical applications. The parameter values for these priors are easily equated to the first two moments of β based on hyperprior knowledge.

Let

f(β_{j}) = (det A_{j})^{\frac{1}{2}} /(2π)^{\frac{P_{j}}{2}} exp \{−\frac{1}{2}(β_{j} − b_{j})′ A_{j}(β_{j} − b_{j})\} , j = 1, 2, …, k  (3.4)

where A_{j} is symmetric positive definite matrix. Then the likelihood function is given by

L(β; x_{1}, x_{2}, …, x_{n}; y_{1}, y_{2}, …, y_{n}) = ∏_{i=1}^{n} {exp[∑_{j=0}^{k} y_{ij}g_{j}(x_{i}; β_{j})/ ∑_{j=0}^{k} exp[g_{j}(x_{i}; β_{j})]}  (3.5)

Now, let D(β) denote the posterior density. Then

D(β) = f(β|x_{1}, x_{2}, …, x_{n}; y_{1}, y_{2}, …, y_{n})

αL(β; x_{1}, x_{2}, …, x_{n}; y_{1}, y_{2}, …, y_{n}). ∏_{j=1}^{k} f(β_{j})

α ∏_{i=1}^{n} {exp[∑_{j=0}^{k} y_{ij}g_{j}(x_{i}; β_{j})/ ∑_{j=0}^{k} exp[g_{j}(x_{i}; β_{j})]} .

exp \{−\frac{1}{2} ∑_{j=1}^{k}(β_{j} − b_{j})′A_{j}(β_{j} − b_{j})\}  (3.6)

By taking natural logarithm we obtain

ln D(β) = ∑_{j=0}^{k} {y_{ij}g_{j}(x_{i}; β_{j}) − ln(∑_{j=0}^{k} exp[g_{j}(x_{i}; β_{j})])}

− \frac{1}{2} ∑_{j=1}^{k}(β_{j} − b_{j})′A_{j}(β_{j} − b_{j})  (3.7)

Now, consider the vector w(β) defined as

w(β) = \frac{∂}{∂β} ln D(β) = [ w_{1}(β) ; ⋮ ; w_{k}(β) ]

= [ \frac{∂}{∂β_{1}} ln D(β) ; ⋮ ; \frac{∂}{∂β_{k}} ln D(β) ]

where

w_{j}(β) = ∑_{i=1}^{n}[y_{ij} − exp(g_{j}(x_{i}; β_{j})/ ∑_{l=0}^{k} exp[g_{l}(x_{i}; β_{l})] · \frac{∂g_{j}}{∂β_{j}} (x_{i}; β_{j})

−A_{j}(β_{j} − b_{j})

and

\frac{∂g_{j}}{∂β_{j}} (x_{i}; β_{j}) = [ \frac{∂g_{j}}{∂β_{j1}} (x_{i}; β_{j}) ; ⋮ ; \frac{∂g_{j}}{∂β_{jP_{j}}} (x_{i}, β_{j}) ]  (3.8)

j = 1, …, k. Define

l(β) = \frac{∂}{∂β} w(β) = [ \frac{∂}{∂β_{1}} w_{1}(β) … \frac{∂}{∂β_{k}} w_{1}(β) ; ⋮ ; \frac{∂}{∂β_{1}} w_{k}(β) … \frac{∂}{∂β_{k}} w_{k}(β)  ]  (3.9)

where l(β) is the information matrix of order (∑_{j=1}^{k} P_{j}) by (∑_{j=1}^{k} P_{j}) which is adjusted for prior information

and if j = j′, j = 1, …, k, then

\frac{∂}{∂β_{j}} w_{j}(β) = ∑_{i=1}^{n}\{\frac{exp[2g_{j}(x_{i}; β_{j})] − exp[g_{j}(x_{i}; β_{j})][∑_{l=0}^{k} exp[g_{l}(x_{i}; β_{l})]}{[∑_{l=0}^{k} exp[g_{l}(x_{i}; β_{l})]]^{2}}.

[\frac{∂}{∂β_{j}} g_{j}(x_{i}; β_{j})] · [\frac{∂}{∂β_{j}} g_{j}(x_{i}; β_{j})]′

+ [y_{ij} − \frac{exp[g_{j}(x_{i}; β_{j})}{∑_{l=0}^{k} exp[g_{l}(x_{i}; β_{l})]}] · \frac{∂^{2}g_{j}(x_{i}; β_{j})}{∂β_{j}∂β′_{j′}}\} − A_{j}  (3.10)

But if j ≠ j′, j = 1, …, k, j′ = 1, …, k, then we have

\frac{∂}{∂β_{j′}} w_{j}(β) = ∑_{i=1}^{n} \frac{exp[g_{j}(x_{i}; β_{j}) + g_{j}, (x_{i}; β_{j}_{′})]}{[∑_{l=0}^{k} exp[g_{l}(x_{i}; β_{l})]]^{2}} · (\frac{∂}{∂β_{j}_{′}} g_{j}_{′}(x_{i}; β_{j}_{′})) · (\frac{∂}{∂β_{j}} g_{j}(x_{i}; β_{j}))′  (3.11)

The first order approximation is defined by

f(y_{n+1}|x_{1}, x_{2}, …, x_{n+1}; y_{1}, y_{2}, …, y_{n}) = f(y_{n+1}|x_{n+1}, β^{∗})

where β^{∗} denotes the posterior mode, which can be found by solving the vector equation

w(β) = 0  (3.12)

For this, an initial guess β^{0} is made for β^{∗} and the iterative algorithm is given by

β^{(l)} = β^{(l−1)} + [l(β^{(l−1)})]^{−1}w(β^{(l−1)})  (3.13)

for l = 1, 2, … until convergence is achieved according to satisfaction of the criterion ∥ β^{(l)} − β^{(l−1)} ∥≤∈ for a chosen norm ∥ · ∥ and small positive value ∈.

Having achieved convergence, β^{∗} is taken to be β^{(l)} having estimated β = (β_{1}, β_{2}, …, β_{k}), our apparent predictive distribution for polytomous nonlinear logistic model is

f(y_{n+1}|x_{1}, x_{2}, …, x_{n+1}; y_{1}, y_{2}, …, y_{n})

= exp{∑_{j=0}^{k}[y_{n+1,j} g_{j}(x_{n+1}; β_{j})]}/ ∑_{j=0}^{k} exp{g_{j}(x_{n+1}; β_{j})}

It should be mentioned that the above method takes into account full prior knowledge which is available in most studies. If no prior information is available, set A_{j} = 0, and b_{j} = 0 for j = 1, …, k in equation (2.6).

4. AN EXAMPLE.

As an example, consider the data appearing in Percy (1993) consisting of a sample of n = 100 patients who were admitted to intensive care units in the UK with acute respiratory and renal failure. The respose y is died (0) or survived (1) and covariates are age, shunt (a measure of lung damage), cieatinine (mg/100ml), albumin (g/100ml) and sepsis (0 = absent or 1 = present).

Percy used linear logistic model with linear link function (g(x, β) = β_{0} + β_{1}x_{1} + … + β_{5}x_{5}) with prior information A = diag (1, 100, .01, 10000, 100, 1) and β = (−1, −0.1, −10, 0.01, 0.1, 1)′. For this data he obtained β^{∗} = (−1.2574, −0.052341, −7.8132, 0.0023812, 0.084193, 1.4057)′ if the vector of covariates for the next patient admitted to an intensive care unit with acute respiratory and renal failure is

x_{101} = (1, 49, 0.21, 706, 38, 1)′

then Percy approximated the predictive distribution of y_{101} by

f(y_{101}|x − 1, x − 2, …, x_{101}; y_{1}, y_{2}, …, y_{100}) = \{ \begin{matrix}0.305 & if y_{101} = 0 \\ 0.695 & if y_{101} = 1. \end{matrix}

References

J. Aitchison and I.R. Dunsmore, Statistical Prediction Analysis, Cambridge University Press, 1975. S. Geisser, The inferential use of predictive distributions, Foundations of Statistical Inference (V.P. Godambe and D.A. Spratt,

eds.), Holt, Rinehart and Winston, 1971, pp. 456–469. S. Geisser, The estimation of distribution functions and prediction of future observations, Proceedings of the conference on

Recent Developments in Statistical Methods and Applications, Taiwan: Academia Sinica (1980), 193–208. S. Geisser, Aspects of the predictive and estimative approaches in the determination of probabilities, A Supplement to Biometrics

38, March (1982). D.F. Heitjan, nonlinear modeling of serial immunologic data: a case study, Journal of the American Statistical Association 20

(1991), no. 416. D.F. Percy, Prediction for generalized linear models, Journal of Applied Statistics 20 (1993), no. 2. B. Rosner, Multivariate methods for clustered binary data with more than one level of nesting, Journal of the American

Statistical Association 84 (1989), no. 406.

Department of Mathematical Sciences, Cameron University, Lawton, Oklahoma

Department of Biostatistics, University of Alabama at Birmingham, Birmingham, Alabama