Finger Exercises To Be Discussed Fri Nov 1, 2002

1. Logistic regression / theory. One formulation of the logistic regression model is as follows. Given data y_{1}, …, y_{n} with corresponding covariate (row) vectors x_{1}^{T} , …, x_{n}^{T}, the model is that

y_{i}|x_{i} ∼^{indep} Bin (1, \frac{1}{1 + exp[−x_{i}^{T} θ]}) , i = 1, …, n

Do the following, starting from our discussion of scoring in generalized linear models in class. Compare with Monahan, p. 218.

(a) Show that the gradient of the log-likelihood (the score function) for logistic regression

is

∇ℓ_{n}(θ) = X^{T} [Y − μ(θ)]

where

X = [ x_{1}^{T} ; x_{2}^{T} ; ⋮ ; x_{n}^{T} ]_{n×p}

is the matrix of covariates (model matrix or design matrix as in linear regression), Y_{n×1} = (y_{1}, …, y_{n})^{T} is the data (a column vector of 0’s and 1’s),

μ(θ) = [ p_{1} ; p_{2} ; ⋮ ; p_{n} ]_{n×1} = [ \frac{1}{1+exp[−x_{1}^{T} θ]} \frac{1}{1+exp[−x_{2}^{T} θ]} ⋮ \frac{1}{1+exp[−x_{n}^{T} θ]} ]_{n×1}

is the vector of means (p_{i} = E[y_{i}|θ]), and θ = (θ_{1}, …, θ_{p})^{T} is the vector of “regression coefficients”.

(b) Show that the negative Hessian and the Fisher Information for logistic regression are the same value, X^{T} ΣX:

−∇^{2}ℓ_{n}(θ) = I(θ) = X^{T} ΣX

where X is as before and now

Σ = [ \begin{matrix}p_{1}(1 − p_{1}) & 0 & ⋯ & 0 \\ 0 & p_{2}(1 − p_{2}) & ⋯ & 0 \\ ⋮ & ... & ⋮ & 0 \\ 0 & ⋯ & 0 & p_{n}(1 − p_{n}) \end{matrix} ]

2. Logistic regression / example. Cox (1970) presented the following data on ingot manufacture.

Condition c Number of trials n_{c} Number of successes y_{c} Covariate x_{c1} Covariate x_{c2}
1  55  0  1  7
2  157  2  1  14
3  159  7  1  27
4  16  3  1  57

(a) Formulate the problem as a logistic regression problem (trials are independent of each other within and between conditions; covariates are the same for all trials within a single condition. What is the total number of trials?).

(b) Find the MLE’s of the regression coefficeints, and their standard errors, using either the methods of problem 2 or using a general optimizer like ms() in Splus or R.

(c) [Optional.] Compare your results to what you would get from glm() in Splus or R.

3. Delta Method. The multivariate delta method is a way of obtaining approximate standard errors for transformed MLE’s. Suppose θ̂ is the MLE of the p-dimensional parameter θ, and suppose τ = g(θ) is some r-dimensional transformation of θ. Then by Taylor’s theorem,

τ̂ − τ = g(θ̂) − g(θ) ≈ ∇g(θ)(θ̂ − θ)

and so, approximately, the variance-covariance matrix of τ̂ may be related to the variancecovariance matrix of θ̂:

Var_{τ}(τ̂) = ∇g(θ)Var_{θ}(θ̂)∇g(θ)^{T} .

Multiplying the Taylor expansion by \sqrt{n}, we get a similar result relating asymptotic normality of \sqrt{n}(τ̂ − τ) to \sqrt{n}(θ̂− θ). (Note: If g(θ) is an r × 1 vector valued function of the p × 1 vector variable θ, then ∇g(θ) = [\frac{∂g_{i}}{∂g_{j}} ]_{r×p}.)

(a) The fitted model P̂[Y = 1|x_{1}, x_{2}] = e^{θ̂_{1}x_{1}+θ̂_{2}x_{2}} can be used to estimate the probability of success at various levels (x_{1}, x_{2}) of the covariates. Compute this probability, for x_{1} = 1, x_{2} = 42 in problem 3. Use the multivariate delta method to compute an approximate 95% confidence interval for this probability.

(b) Plot the fitted response curve P̂(x) = P̂[Y = 1|x_{1} = 1, x_{2} = x], over the range x = 0 to x = 100. On the same graph, plot an approximate 95% confidence envelope (that is, plot the curves corresponding to the upper and lower edpoints of approximate 95% confidence intervals for P (x) at each value of x, using the multivariate delta method. (This can be done with a single matplot(...,type="l") command, or with several plot and lines commands).

4. Exponential families / theory. The following items were needed in the discussion of the scoring algorithm for exponential family models.

(a) Review assumptions and argument needed to assert that, when θ = (θ_{1}, …, θ_{p}),

I(θ) = E_{θ}[−∇^{2}ℓ_{n}(θ)] = Var_{θ}[∇ℓ_{n}(θ)] = E_{θ}[∇_{θ}ℓ_{n}(θ)^{2}]

where ℓ_{n}(θ) = log ∏_{i} f(y_{i}|θ) is the log-likelihood of iid data y_{1}, …, y_{n}. (Notation: E_{θ}[·] = E[·|θ], i.e. it is the expectation over the data y_{1}, …, y_{n} for fixed parameter value θ.)

(b) Let L(θ) = g(Y) exp[B(θ) + K(Y)^{T} γ(θ)] be an exponential family likelihood, with μ(θ) = E_{θ}[K(Y)] and Σ(θ) = Var_{θ}[K(Y)]. Complete the argument (assuming that integration and differentiation can be interchanged) that shows

∇μ(θ) = ∇ ∫ K(Y)L(θ)dν(Y)
= ⋯
= Σ(θ)∇γ(θ)

Note: If μ(θ) is an r × 1 vector valued function of the p × 1 vector variable θ, then

∇μ(θ) = [\frac{∂μ_{i}}{∂θ_{j}} ]_{r×p}