### SOURCE: https://rpubs.com/trjohns
library(SMPracticals)
## Loading required package: ellipse
##
## Attaching package: 'ellipse'
## The following object is masked from 'package:graphics':
##
## pairs
## Loading required package: msme
## Loading required package: MASS
##
## Attaching package: 'MASS'
## The following objects are masked from 'package:SMPracticals':
##
## cement, forbes, leuk, shuttle
## Loading required package: lattice
##
## Attaching package: 'lattice'
## The following object is masked from 'package:SMPracticals':
##
## barley
## Loading required package: sandwich
library(ggplot2)
data(fishing) # necessary as data are not automatically loaded for this package
p <- ggplot(fishing, aes(x = meandepth, y = totabund/sweptarea))
p <- p + geom_point() + facet_wrap(~ period) + theme_bw()
p <- p + xlab("Mean Trawl Depth") + ylab("Fish per Unit Trawl Area")
plot(p)
## Don't know how to automatically pick scale for object of type labelled/integer. Defaulting to continuous.
## Don't know how to automatically pick scale for object of type labelled/integer. Defaulting to continuous.
# The smoking variable is just coded as 0,1 but it would be nicer
# to create a factor with clearly labeled levels as follows.
smoking$smoker <- factor(smoking$smoker, labels = c("no", "yes"))
mymodel <- lm(dead/(alive + dead) ~ smoker + age, data = smoking)
summary(mymodel)$coefficients
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 0.006780201 0.02288455 0.29627851 7.770050e-01
## smokeryes 0.038932267 0.01618182 2.40592605 5.286637e-02
## age25-34 0.001774007 0.03027342 0.05859949 9.551740e-01
## age35-44 0.066899469 0.03027342 2.20984204 6.915833e-02
## age45-54 0.154522896 0.03027342 5.10424370 2.212161e-03
## age55-64 0.360782052 0.03027342 11.91745405 2.113521e-05
## age65-74 0.768004312 0.03027342 25.36893406 2.471213e-07
## age75+ 0.973753666 0.03027342 32.16530449 6.003231e-08
library(COUNT)
smoking$yhat <- predict(mymodel)
p <- ggplot(smoking, aes(x = age, y = dead/(alive + dead), fill = smoker))
p <- p + scale_fill_manual(values = c(grey(0.5), grey(0.95)))
p <- p + geom_point(aes(size = alive + dead), shape = 21, color = "black")
p <- p + geom_line(aes(y = yhat, group = smoker), color = "black")
p <- p + geom_point(aes(y = yhat, fill = smoker), shape = 21, color = "black")
p <- p + xlab("Age Group") + ylab("Proportion Deceased after 20 Years")
p <- p + labs(size = "Number of Women", fill = "Smoker?") + theme_bw()
plot(p)
