Re: [R] Firths bias correction for log-linear models

From: David Firth <>
Date: Fri 13 Jan 2006 - 09:00:23 EST

On 12 Jan 2006, at 20:54, wrote:

> Dear R-Help List,
> I'm trying to implement Firth's (1993) bias correction for log-linear
> models.
> Firth (1993) states that such a correction can be implemented by
> supplementing
> the data with a function of h_i, the diagonals from the hat matrix,
> but doesn't
> provide further details. I can see that for a saturated log-linear
> model, h_i=1
> for all i, hence one just adds 1/2 to each count, which is equivalent
> to the
> Jeffrey's prior, but I'd also like to get bias corrected estimates for
> other log
> linear models. It appears that I need to iterate using GLM, with the
> weights
> option and h_i, which I can get from the function hatvalues. For
> logistic
> regression, this can be performed by splitting up each observation
> into response
> and nonresponse, and using weights as described in Heinze, G. and
> Schemper, M.
> (2002), but I'm unsure of how to implement the analogue for log-linear
> models. A
> procedure using IWLS is described by Firth (1992) in Dodge and
> Whittaker (1992),
> but this book isn't in the local library, and its $141+ on Amazon.
> I've tried
> looking at the code in the logistf and brlr libraries, but I haven't
> had any
> (successful) ideas. Can anyone help me in describing how to implement
> this in R?

I don't recommend the adjusted IWLS approach in practice, because that algorithm is only first-order convergent. It is mainly of theoretical interest.

The brlr function (in the brlr package) provides a template for a more direct approach in practice. The central operation there is an application of optim(), with objective function

I hope that helps. Please feel free to contact me off the list if anything is unclear.

Kind regards,

Professor David Firth

______________________________________________ mailing list
PLEASE do read the posting guide!
Received on Fri Jan 13 09:07:09 2006

This archive was generated by hypermail 2.1.8 : Fri 13 Jan 2006 - 14:19:48 EST