HW2
# 先讀檔案，然後分IQ跟Size的水準
dta2 <- read.table("R/langMathDutch.txt", header = TRUE) %>%
mutate(Size = cut(size, quantile(size, probs = c(0, .33, .67, 1)),
c("Small", "Midium", "Large"),
include.lowest = TRUE, ordered = TRUE),
IQ = cut(IQV, quantile(IQV, probs = c(0, .33, .67, 1)),
c("Low", "Middle", "High"),
include.lowest = TRUE, ordered = TRUE),
Size.IQ = factor(paste(Size, IQ, sep = "."),
levels = c("Small.Low", "Small.Middle", "Small.High",
"Midium.Low", "Midium.Middle", "Midium.High",
"Large.Low", "Large.Middle", "Large.High")))
## Warning: package 'bindrcpp' was built under R version 3.4.4
# 繪圖
ggplot(dta2, aes(lang, arith))+
geom_point()+
stat_smooth(method = "lm")+
facet_wrap(~ Size.IQ)+
scale_y_continuous(breaks = seq(5, 35, 5))+
labs(x = "Language Socre", y = "Arithmetic Score")
HW3
# data
dta3 <- read.table("R/stateAnxiety.txt", h = T)
dta3 <- dta3 %>%
gather(key = "key", value = "Score") %>%
mutate(Gender = c(rep("Female", 250), rep("Male", 250)),
ID = rep(1:50, 10)) %>%
mutate(Time = rep(c(rep(-5, 50), rep(-4, 50), rep(-3, 50),
rep(-2, 50), rep(-1, 50)), 2))
# 看一下性別
pd <- position_dodge(.2)
ggplot(dta3, aes(Time, Score, group = Gender, color = Gender)) +
geom_point(position = pd) +
stat_smooth(aes(color = Gender), method = "lm") +
labs(x = "Weeks before exam", y = "Anxiety score") +
theme(legend.direction = "horizontal") +
theme(legend.position = "top")
# 再看不同時間的個體
ggplot(dta3, aes(Time, Score, group = ID)) +
geom_point() +
stat_smooth(method = "lm", se = F) +
labs(x = "Weeks before exam", y = "Anxiety score")
# 最後回答問題，看起來在性別方面，在焦慮程度上均存在差異；而在個體方面，似乎斜率沒有太大差異
HW4
# 讀檔案
dta4 <- read.table("R/math_attainment.txt", header = TRUE)
# 先畫基本散佈圖
ggplot(dta4, aes(math1, math2, color = cc))+
geom_point(size = rel(2))+
stat_smooth(method = "lm", color = "steelblue")+
theme(legend.position = c(.8, .2))+
labs(x = "Math score at Year 1",
y = "Math score at Year 2",
color = "Curriculum coverage")
# model
m <- lm(math2 ~ -1 + math1 + cc, dta4)
# standardized residual plot
dta4 %>% mutate(s.resi = scale(residuals(m)),
pred = fitted.values(m)) %>%
ggplot(., aes(pred, s.resi))+
geom_point()+
geom_hline(yintercept = 0, linetype = "dotted")+
labs(x = "Fitted values", y = "Standardized residuals")
# Q-Q plot
ggqqplot <- function(data)
{
y <- quantile(data[!is.na(data)], c(0.25, 0.75))
x <- qnorm(c(0.25, 0.75))
slope <- diff(y)/diff(x)
int <- y[1L] - slope * x[1L]
d <- data.frame(resids = data)
ggplot(d, aes(sample = resids))+
stat_qq()+
geom_abline(slope = slope, intercept = int)
}
ggqqplot(scale(resid(m)))
# Model fit
dta4 %>% mutate(pred = fitted.values(m)) %>%
ggplot(., aes(math1, math2, color = cc))+
geom_point(alpha = .7, pch = 21, size = rel(2))+
geom_point(aes(math1, pred, color = cc), size = rel(2))+
stat_smooth(aes(math1, pred, color = cc), method = "lm", color = "steelblue")+
theme(legend.position = c(.8, .2))+
labs(x = "Math score at Year 1",
y = "Fitted Values Math score at Year 2",
color = "Curriculum coverage")
# Parameter Coefficients
tidy(m, conf.int = TRUE) %>%
ggplot(., aes(term, estimate))+
geom_pointrange(aes(ymin = conf.low, ymax = conf.high))+
geom_hline(yintercept = 0, linetype = "dotted")+
scale_x_discrete(labels = c("Curriculum coverage", "Math score at Year 1"))+
labs(x = "term", y = "parameter coefficients")
HW5
dta5 <- read.table("R/hs0.txt", header = TRUE)
# 確立變項不同維度的特徵
cor <- cov(dta5$write, dta5$read)
v <- matrix(c(1,cor,cor,1), 2, 2)
slope <- eigen(v)$vector[1, 1]
intercept <- mean(dta5$write)-slope*mean(dta5$read)
# 繪圖
p <- ggplot(dta5, aes(read, write)) +
geom_point(shape = 21) +
stat_ellipse() +
geom_abline(intercept = 0, slope = 1, color = "gray") +
geom_abline(intercept = intercept, slope = slope, color = "tomato", size = rel(2)) +
geom_vline(xintercept = mean(dta5$read), color = "gray") +
geom_hline(yintercept = mean(dta5$write), color = "gray") +
stat_smooth(method = "lm", size = rel(1.1)) +
xlim(25, 80) +
ylim(25, 80) +
labs(x = "Reading score", y = "Writing score")
ggMarginal(p, type = "histogram",
xparams = list(binwidth = (IQR(dta5$read)*2)/(200)^(1/3),
fill = "gray"),
yparams = list(binwidth = (IQR(dta5$write)*2)/(200)^(1/3),
fill = "gray"))
HW7
# 處理檔案
dta7 <- USPersonalExpenditure %>% as.data.frame %>%
mutate(index = 1:5, category = c("Food\nTobacco", "Household\nOperation", "Medical\nHealth", "Personal\nCare", "Private\nEducation")) %>%
gather(year, exp, 1:5) %>%
group_by(category) %>%
mutate(exp_c = log(exp) - mean(log(exp)))
# 先畫第一部分
qplot(category, exp_c, data = dta7) +
geom_hline(yintercept = 0, colour = "grey50") +
geom_line(aes(group = 1)) +
facet_wrap(~ year) +
labs(x = "categories", y = "Centered personal expenditures")
# 再畫第二部分
qplot(exp_c, category, data = dta7) +
geom_segment(aes(xend = 0, yend = category)) +
geom_vline(xintercept = 0, colour = "grey50") +
facet_wrap(~ year, nrow = 1) +
labs(x = "Centered personal expenditures (log)", y = "categories")
HW8
dta8 <- read.table("R/hs0.txt", header = TRUE)
ggplot(dta8, aes(write, read, fill = female))+
stat_density2d(aes(color = ..level..))+
theme(legend.position = c(.05, .75))+
geom_point(pch = 21, alpha = .5)+
facet_wrap(~female)+
labs(x = "Writing score", y = "Reading score")
HW9
# 先處理檔案
dta9 <- data.frame(Race = c("White", "White", "Black", "Black"),
Gender = c("Male", "Female", "Male", "Female"),
Yes = c(43, 26, 29, 22),
No = c(134, 149, 23, 36)) %>%
mutate(P_yes = Yes/(Yes+No), P_no = No/(Yes+No))
# computate
new.dta <- dta9 %>% select(1:4) %>%
gather(Intercourse, Freq, 3:4) %>%
mutate(P = c(dta9$P_yes, dta9$P_no),
SE = sqrt(P*(1-P)/Freq))
# 繪圖
ggplot(new.dta, aes(Intercourse, P, fill = Gender))+
geom_bar(position = "dodge", stat = "identity")+
geom_errorbar(aes(ymin = P - SE, ymax = P + SE), position = "dodge")+
facet_wrap(~ Race)+
labs(y = "Proportion")
HW11
dta11 <- MASS::Cushings %>%
mutate(Type = factor(Type, levels = c("u", "b", "c", "a"),
labels = c("Unknown", "Bilateral Hyperplasia", "Carcinoma", "Adenoma")),
Label = "")
dta11$Label[c(1, 13, 21, 27)] = c("Adenoma", "Bilateral Hyperplasia", "Carcinoma", "Unknown")
ggplot(dta11, aes(Tetrahydrocortisone, Pregnanetriol, fill = Type))+
geom_point(pch = 21, size = rel(2))+
geom_text_repel(aes(label = Label, color = Type))+
scale_fill_discrete(guide = FALSE)+
scale_color_discrete(guide = FALSE)+
theme_hc()+
labs(x = "Tetrahydrocortisone (mg/24 hours)", y = "Pregnanetriol (mg/24 hours)",
title = "Cushings's syndrome")+
theme(plot.title = element_text(hjust = 1))
HW12
library(car)
## Warning: package 'car' was built under R version 3.4.2
##
## Attaching package: 'car'
## The following object is masked from 'package:dplyr':
##
## recode
## The following object is masked from 'package:purrr':
##
## some
dta12 <- Vocab
ggplot(dta12, aes(education))+
geom_histogram(aes(y = ..density..), breaks = c(0, 6, 9, 12, 16, 20),
fill = "gray", color = "black")+
scale_y_continuous(breaks = seq(0, 0.15, 0.05), limits = c(0, 0.15),
minor_breaks = seq(0, 0.15, 0.01)) +
scale_x_continuous(breaks = c(0, 6, 9, 12, 16, 20))+
theme(panel.grid.major.x = element_blank(),
panel.grid.minor.x = element_blank(),
panel.grid.major.y = element_line(linetype = "dotted", color = "gray50"),
panel.grid.minor.y = element_line(linetype = "dotted", color = "gray50"))+
labs(x = "Education level (years)", y = "Density")
