Phân tích biểu đồ đơn giản

Làm quen với những package sau đây: lattice, ggplot2, gridExtra

Task 1: Biểu đồ phân bố với hist()

Đọc dữ liệu PISA Data Vietnam 2015.csv và gọi là đối tượng pisa

pisa = read.csv("C:\\Users\\Nguyen\\Desktop\\TDT Workshop 2020\\TDTU Datasets for 2020 Workshop\\PISA Data Vietnam 2015.csv")

Thử nghiệm với hist()

hist(pisa$Math, col="blue")

hist(pisa$Math, col="blue", border="white")

hist(pisa$Science, col="blue" , border="white", xlab="Science Score", ylab="Frequency", main="Distribution of Science Scores")

Thêm đường biểu diễn

hist(pisa$Science, col="blue", border="white", prob=T)
lines(density(pisa$Science), col="red")

Vẽ hai phân bố cho nam và nữ

p1 = hist(pisa$Science[pisa$Gender=="Boys"] , plot=F)
p2 = hist(pisa$Science[pisa$Gender=="Girls"] , plot=F)
plot(p1, col="skyblue", border="white")
plot(p2, add=T, col=scales::alpha("green", 0.5), border="white")

Task 2: Biểu đồ phân bố với lattice()

Vẽ biểu đồ phân bố cho nhiều nhóm dùng hàm densityplot

library(lattice)
densityplot(~Science, groups=Gender, data=pisa)

densityplot(~Science | Gender, data=pisa)

densityplot(~Science, groups=Area, data=pisa, auto.key=list(space="top"))

Chia 3 cửa sổ (dùng package lattice và gridExtra)

p1 = densityplot(~Science, groups=Area, data=pisa, auto.key=list(space="top"))
p2 = densityplot(~Math, groups=Area, data=pisa)
p3 = densityplot(~Read, groups=Area, data=pisa)
library(gridExtra)
grid.arrange(p1, p2, p3, ncol=3)

Task 3: Biểu đồ hộp với boxplot()

Chia màn hình ra 3 cửa sổ

par(mfrow=c(1, 3))
boxplot(pisa$Science, col="purple")
boxplot(pisa$Science ~ pisa$Gender, col=c("blue", "red"))
boxplot(pisa$Science ~ pisa$Region, col=c("blue", "red", "purple"))

Task 4: Biểu đồ tương quan với plot()

Tìm hiểu mối liên quan giữa PARED và điểm môn khoa học

plot(pisa$Science ~ pisa$Math, col="blue")

Vẽ thêm đường hồi quy:

plot(pisa$Science ~ pisa$Math, col="blue")
abline(lm(pisa$Science ~ pisa$Math), col="red")

Phân tích biểu đồ chất lượng cao (ggplot2)

Task 5: Biểu đồ tương quan

Đọc dữ liệu obesity data vào R và gọi đối tượng là ob

ob = read.csv("C:\\Users\\Nguyen\\Desktop\\TDT Workshop 2020\\TDTU Datasets for 2020 Workshop (new)\\Obesity data.csv")

Tạo ra biến mới gọi là OB dựa vào biến bmi: nếu bmi<18.5 thì OB=“Underweight”; bmi từ 18.5 đến 24.9 thì OB=“Normal”; bmi từ 25.0 đến 29.9 thì OB=“Overweight”; bmi từ 30.0 trở lên, OB=“Obese”

ob$OB [ob$bmi < 18.5] = "Underweight"
ob$OB [ob$bmi >= 18.5 & ob$bmi < 24.9] = "Normal"
ob$OB [ob$bmi >= 25.0 & ob$bmi < 29.9] = "Overweight"
ob$OB [ob$bmi >= 30] = "Obese"
ob$OB = factor(ob$OB, levels=c("Underweight", "Normal", "Overweight", "Obese"))
table(ob$OB)
## 
## Underweight      Normal  Overweight       Obese 
##         107         857         230          15

Tìm hiểu phân bố của OB

library(ggplot2)
p = ggplot(data=ob, aes(OB, fill=OB)) + geom_bar()
p = p + xlab("Obesity group") + ylab("Frequency")
p + theme(legend.position="none")

Task 6: Vẽ biểu đồ tương quan giữa weight và pcfat dùng ggplot2

Biểu đồ đơn giản

library(ggplot2) ; library(ggthemes); library(gridExtra)
p = ggplot(data=ob, aes(x=weight, y=pcfat))
p = p + geom_point() + geom_smooth()
p = p + xlab("Weight") + ylab("Percent body fat") + ggtitle("Weight and Percent Body Fat")
p = p + theme(plot.title=element_text(hjust=0.5))
p
## `geom_smooth()` using method = 'gam' and formula 'y ~ s(x, bs = "cs")'

Biểu đồ theo nhóm nam nữ

p = ggplot(data=ob, aes(x=weight, y=pcfat, fill=gender,
col=gender))
p = p + geom_point() + geom_smooth()
p = p + xlab("Weight") + ylab("Percent body fat") + ggtitle("Weight and Percent Body Fat for men and women separately")
p = p + theme(plot.title=element_text(hjust=0.5))
p + theme_economist()
## `geom_smooth()` using method = 'loess' and formula 'y ~ x'

Có thể thử theme_tufte(), theme_few(), theme_wsj(), theme_clean(), theme_hc()

Biểu đồ tương quan + biểu đồ phân bố dùng package ggExtra

library(ggExtra)
p = ggplot(data=ob, aes(x=bmi, y=pcfat, fill=gender, col=gender))
p = p + geom_point() + geom_smooth()
ggMarginal(p, type="histogram", groupColour=T, groupFill=T)
## `geom_smooth()` using method = 'loess' and formula 'y ~ x'
## `geom_smooth()` using method = 'loess' and formula 'y ~ x'

Có thể thử nghiệm thêm với type ("density", "violin", "boxplot"):

ggMarginal(p, type="density", groupColour=T, groupFill=T)
## `geom_smooth()` using method = 'loess' and formula 'y ~ x'
## `geom_smooth()` using method = 'loess' and formula 'y ~ x'

Biểu đồ tương quan đa biến (package GGally) Chọn biến quan tâm

library(GGally)
## Registered S3 method overwritten by 'GGally':
##   method from   
##   +.gg   ggplot2
dat = ob[, c("gender", "age", "bmi", "weight", "height", "pcfat")]
ggpairs(dat)
## `stat_bin()` using `bins = 30`. Pick better value with `binwidth`.
## `stat_bin()` using `bins = 30`. Pick better value with `binwidth`.
## `stat_bin()` using `bins = 30`. Pick better value with `binwidth`.
## `stat_bin()` using `bins = 30`. Pick better value with `binwidth`.
## `stat_bin()` using `bins = 30`. Pick better value with `binwidth`.

Thêm màu theo nhóm

ggpairs(data=dat, mapping = aes(color = gender))
## `stat_bin()` using `bins = 30`. Pick better value with `binwidth`.
## `stat_bin()` using `bins = 30`. Pick better value with `binwidth`.
## `stat_bin()` using `bins = 30`. Pick better value with `binwidth`.
## `stat_bin()` using `bins = 30`. Pick better value with `binwidth`.
## `stat_bin()` using `bins = 30`. Pick better value with `binwidth`.

So sánh với:

ggpairs(data=ob, mapping = aes(color = gender), columns = c("age", "weight", "bmi", "pcfat"))

Task 7: Vẽ biểu đồ phân bố (histogram)

Đọc dữ liệu PISA Data Vietnam 2015.csv và gọi đối tượng là pisa

Vẽ biểu đồ phân bố điểm môn khoa học

library(gridExtra)
p = ggplot(data=pisa, aes(x=Science))
p1 = p + geom_histogram(color="white", fill="blue")
p = ggplot(data=pisa, aes(x=Science))
p = p + geom_histogram(aes( y=..density..), color="white", fill="blue")
p2 = p + geom_density(col="red")
grid.arrange(p1, p2, ncol=2)
## `stat_bin()` using `bins = 30`. Pick better value with `binwidth`.
## `stat_bin()` using `bins = 30`. Pick better value with `binwidth`.

Phân bố điểm môn khoa học theo vùng (Area)

p = ggplot(data=pisa, aes(x=Science, fill=Area))
p1 = p + geom_histogram(position="dodge")
p2 = ggplot(data=pisa, aes(x=Science, fill=Area, color=Area)) + geom_density(alpha = 0.1)
grid.arrange(p1, p2, nrow=2)
## `stat_bin()` using `bins = 30`. Pick better value with `binwidth`.

Histogram và xác suất tích lũy

p = ggplot(data=pisa, aes(HISCED))
p = p + stat_ecdf(color="red", lwd=1)
p = p + geom_bar(aes(y = (..count..)/sum(..count..)), fill="blue", colour="yellow")
p
## Warning: Removed 14 rows containing non-finite values (stat_ecdf).
## Warning: Removed 14 rows containing non-finite values (stat_count).

Task 8: Vẽ biểu đồ hộp (box plot)

Biểu đồ hộp theo vùng

p = ggplot(data=pisa, aes(x=Area, y=Science, col=Area, fill=Area))
p1 = p + geom_boxplot(col="black")
p2 = p + geom_boxplot(col="black") + geom_jitter(alpha=0.05)
grid.arrange(p1, p2, ncol=2)

Biểu đồ hộp theo kinh tế

p = ggplot(data=pisa, aes(x=PARED, y=Science, fill=PARED))
p1 = p + geom_boxplot(col="black") + geom_jitter(alpha=0.02)
p = ggplot(data=na.omit(pisa) , aes(x=factor(PARED), y=Science,
fill=factor(PARED), color=factor(PARED)))
p2 = p + geom_boxplot(col="black") + geom_jitter(alpha=0.05)
grid.arrange(p1, p2, ncol=2)
## Warning: Continuous x aesthetic -- did you forget aes(group=...)?
## Warning: Removed 14 rows containing missing values (stat_boxplot).
## Warning: Removed 14 rows containing missing values (geom_point).