The code below uses ggplot with stat_smooth(method="glm",
family=binomial, ...)
to plot the data on survival of passengers on the Titanic, with the
logistic regression
curves for each sex on the scale of Pr(survived). This works (quite
nicely!) because
I've explicitly transformed the factor survived to 0/1 in the ggplot call.
Some questions:
- Is it possible, and if so, how, to plot the same data and fitted
smooths on the logit
scale, i.e., the linear predictor for the binomial glm?
- the response, survived, is a factor. Is it possible to avoid using
as.numeric(survived)-1
in the call to ggplot()? This is cosmetic, but requires an extra bit of
explanation to use
in teaching or writing.
i.e., glm() is quite happy to fit the model survived ~ age+sex
in the binomial family, and gives the same predicted probabilities and
logits.
install.packages("vcdExtra")# data from the most recent version,
vcdExtra_0.5-11
data(Titanicp, package="vcdExtra")
str(Titanicp)
'data.frame': 1309 obs. of 6 variables:
$ pclass : Factor w/ 3 levels "1st","2nd","3rd": 1 1 1 1 1 1 1 1 1 1 ...
$ survived: Factor w/ 2 levels "died","survived": 2 2 1 1 1 2 2 1 2 1 ...
$ sex : Factor w/ 2 levels "female","male": 1 2 1 2 1 2 1 2 1 2 ...
$ age : num 29 0.917 2 30 25 ...
$ sibsp : num 0 1 1 1 1 0 1 0 2 0 ...
$ parch : num 0 2 2 2 2 0 0 0 0 0 ...
>
require(ggplot2)
# remove missings on age
Titanicp <- Titanicp[!is.na(Titanicp$age),]
ggplot(Titanicp, aes(age, as.numeric(survived)-1, color=sex)) +
stat_smooth(method="glm", family=binomial, formula=y~x, alpha=0.2,
size=2, aes(fill=sex)) +
geom_point(position=position_jitter(height=0.02, width=0), size=1.5)
# equivalent logistic regression model, survived as a factor
mod <- glm(survived ~ age+sex, family=binomial, data=Titanicp)
summary(mod)
--
Michael Friendly Email: friendly AT yorku DOT ca
Professor, Psychology Dept. & Chair, Quantitative Methods
York University Voice: 416 736-2100 x66249 Fax: 416 736-5814
4700 Keele Street Web: http://www.datavis.ca
Toronto, ONT M3J 1P3 CANADA
______________________________________________
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.
