On Fri, 16 May 2008, Rolf Turner wrote:
I am fitting a logistic binomial model of the form
glm(y ~ a*x,family=binomial)
where a is a factor (with 5 levels) and x is a continuous predictor.
To assess how much ``impact'' x has, I want to compare the fitted
success probability when x = its maximum value with the fitted
probability when x = its mean value. (The mean and the max are to be
taken by level of the factor ``a'', but that's not really an issue.)
I can of course easily calculate p.hat(x.max) - p.hat(x.mean) using
predict() (with type="response"). And I can get the standard error for
p.hat(x.max) and p.hat(x.mean) by specifiying se.fit=TRUE. No problem
there.
But how can I get a handle on the standard error of the difference?
The same way predict() does? It used a crude (delta-function)
transformation of SEs on linear predictor scale.
I am not sure that a difference in probabilites is interesting unless it
is small compared to the values, when linearization will be most likely be
satisfactory. (Why not ratios or log odds, for example?)
Linearization can also be used if the SEs on each probability are small
but the difference is not. Otherwise you may as well simulate from the
joint distribution of the coefficients, compute the difference for the
similations and summarize. (That might be as quick to do as to work out
the analytical approximations.)
In a linear model this would just be SE(beta_1.hat)*(x.max-x.mean) (where
beta_1.hat is specific to the particular level of `a' being considered).
If I am not mistaken. (Please correct me if I am!)
But in the logistic model, everything is entangled in the inverse link
function (the ``expit'' function as it is called by some), and I can see
no way of disentangling.
Is there any way of getting at this? I figure that simulation/Monte
Carlo inference/ parametric bootstrapping would provide a workaround,
but before I go that route, can anyone point me to a simpler method?
There wouldn't be anything built into R or an R package, would there?
(I did a fairly basic RSiteSearch() and came up empty handed.)
Thanks for any tips.
cheers,
Rolf Turner
--
Brian D. Ripley, [EMAIL PROTECTED]
Professor of Applied Statistics, http://www.stats.ox.ac.uk/~ripley/
University of Oxford, Tel: +44 1865 272861 (self)
1 South Parks Road, +44 1865 272866 (PA)
Oxford OX1 3TG, UK Fax: +44 1865 272595
______________________________________________
R-help@r-project.org mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
and provide commented, minimal, self-contained, reproducible code.
