Trực quan hoá dữ liệu bằng biểu đồ với ggplot2

Cài packages “ggplot2” Gọi các packages:

library(readr)
library(knitr)

– Load và view dữ liệu (Lưu ý: Nếu dùng dữ liệu trên mày thì dùng các lệnh khác để đọc dữ liệu hoặc thay thế đường dẫn đến folder lưu dữ liệu của bạn)

Litho <- "https://raw.githubusercontent.com/ngocdlu/Litho/master/Litho.csv"
data <- read_csv(Litho)
## Parsed with column specification:
## cols(
##   len = col_double(),
##   wid = col_double(),
##   rat = col_double(),
##   cir = col_double(),
##   pet = col_double(),
##   species = col_character()
## )
attach(data)

Gọi các packages “dplyr” và “ggplot2”

library(dplyr)
## 
## Attaching package: 'dplyr'
## The following objects are masked from 'package:stats':
## 
##     filter, lag
## The following objects are masked from 'package:base':
## 
##     intersect, setdiff, setequal, union
library(ggplot2)

Vẽ 1 biểu đồ boxplot đơn giản

data %>% 
  ggplot(aes(species, len, fill = species)) + 
  geom_boxplot() +
  labs(title = "Overview of Leaf_Length by Species", 
       x = "Species",
       y = "Length (cm)") +
  scale_y_continuous(breaks = seq(0,30, by = 5))

Tính giá trị trung bình Tạo thêm 1 dataframe tính độ dài trung bình của phiến lá

mean_length_by_species <- data %>% 
  group_by(species) %>% 
  summarise(mean_len = mean(len)) %>% 
  as.data.frame

Xem dữ liệu mới tổng hợp

mean_length_by_species

Biểu diễn độ dài trung bình lá bằng những điểm màu đỏ tương ứng trên biểu đồ. Dùng geom_point() để vẽ biểu đồ điểm

p1 <- data %>% 
  ggplot(aes(species, len, fill = species)) + 
 geom_boxplot() +
 labs(title = "Overview of Leaf Length by Species", 
            x = "Species",
            y = "Length (cm)") +
 scale_y_continuous(breaks = seq(0,30, by = 5)) +
 geom_point(aes(x = species, 
            y = round(mean_len,1)), 
            data = mean_length_by_species, 
            col = "red") + 
 geom_text(aes(label = round(mean_len,1), 
            x = species, 
            y = round(mean_len,1)), 
            data = mean_length_by_species, 
            check_overlap = TRUE, 
            vjust = -0.5)
p1

Thự hiện tương tự đối với các biến khác Chiều rộng

mean_width_by_species <- data %>% 
  group_by(species) %>% 
  summarise(mean_wid = mean(wid)) %>% 
  as.data.frame
mean_width_by_species

Boxplot cho chiều rộng phiển lá

p2 <- data %>% 
  ggplot(aes(species, wid, fill = species)) + 
  geom_boxplot() +
  labs(title = "Overview of Leaf Width by Species", 
       x = "Species",
       y = "Width (cm)") +
  scale_y_continuous(breaks = seq(0,30, by = 2)) +
  geom_point(aes(x = species, 
                 y = round(mean_wid,1)), 
             data = mean_width_by_species, 
             col = "red") + 
  geom_text(aes(label = round(mean_wid,1), 
                x = species, 
                y = round(mean_wid,1)), 
            data = mean_width_by_species, 
            check_overlap = TRUE, 
            vjust = -0.5)
p2

Giá trị trung bình cho dạng tròn phiến lá (Circularity)

mean_Cir_by_species <- data %>% 
  group_by(species) %>% 
  summarise(mean_cir = mean(cir)) %>% 
  as.data.frame
mean_Cir_by_species

Biểu đồ

p3 <- data%>% 
ggplot(aes(species, cir, fill = species)) + 
geom_boxplot() +
labs(title = "Overview of Leaf Circularity by Species", 
     x = "Species",
     y = "Circularity") +
scale_y_continuous(breaks = seq(0,1, by = 0.1)) +
geom_point(aes(x = species, 
     y = round(mean_cir,2)), 
     data = mean_Cir_by_species, 
     col = "red") + 
geom_text(aes(label = round(mean_cir,2), 
     x = species, 
     y = round(mean_cir,2)), 
     data = mean_Cir_by_species, 
     check_overlap = TRUE, 
     vjust = -0.5)

p3

Tỉ lệ lá (Aspect_ratio)

mean_aspect_by_species <- data %>% 
  group_by(species) %>% 
  summarise(mean_aspect = mean(rat)) %>% 
  as.data.frame
mean_aspect_by_species

Boxplot

p4 <- data %>% 
ggplot(aes(species, rat, fill = species)) + 
geom_boxplot() +
labs(title = "Overview of Leaf_Aspect_ratio by Species", 
       x = "Species",
       y = "Aspect_ratio") +
scale_y_continuous(breaks = seq(0,5, by = 1)) +
geom_point(aes(x = species, 
       y = round(mean_aspect,1)), 
       data = mean_aspect_by_species, 
       col = "red") + 
geom_text(aes(label = round(mean_aspect,1), 
       x = species, 
       y = round(mean_aspect,1)), 
       data = mean_aspect_by_species, 
       check_overlap = TRUE, 
       vjust = -0.5)

p4

Chiều dài cuống lá (pet)

mean_pet_length_by_species <- data %>% 
  group_by(species) %>% 
  summarise(mean_pet = mean(pet)) %>% 
  as.data.frame
mean_pet_length_by_species

Boxplot

p5 <- data %>% 
  ggplot(aes(species, pet, fill = species)) + 
  geom_boxplot() +
  labs(title = "Overview of pet length by Species", 
       x = "Species",
       y = "Petiole length (cm)") +
  scale_y_continuous(breaks = seq(0,4, by = 1)) +
  geom_point(aes(x = species, 
                 y = round(mean_pet,1)), 
             data = mean_pet_length_by_species, 
             col = "red") + 
  geom_text(aes(label = round(mean_pet,1), 
                x = species, 
                y = round(mean_pet,1)), 
            data = mean_pet_length_by_species, 
            check_overlap = TRUE, 
            vjust = -0.5)

p5

Vẽ biểu đồ Scatter plot mối tương quan giữa Apspect ratio và Circularity

qplot(rat, cir, data=data)

Fill màu theo Species

p6 <- qplot(rat, cir, data=data, facets=species~., col=species)
p6

Sắp xếp các biểu đồ trên cùng một trang install.packages(“ggpubr”)

library(ggpubr)
## Loading required package: magrittr

Thực hiện sắp xếp các biểu đồ trên cùng một trang bằng lệnh:

ggarrange(p1, p2, p3, p4, labels = c("A", "B", "C", "D", "E", "F"),
            ncol = 2, nrow = 2)

Phân tích phương sai và hậu định

VD: Phân tích phương sai của của biến chiều dài giữa các loài

av <-aov(len~species)
summary(av)
##             Df Sum Sq Mean Sq F value Pr(>F)    
## species      2 1299.8   649.9   97.95 <2e-16 ***
## Residuals   54  358.3     6.6                   
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1

Phân tích hậu định bằng Tukey’s HSD test

t<-TukeyHSD(av)
t
##   Tukey multiple comparisons of means
##     95% family-wise confidence level
## 
## Fit: aov(formula = len ~ species)
## 
## $species
##                                   diff        lwr        upr     p adj
## L_campylolepis-L_cadamensis  -8.638421 -10.652467 -6.6243755 0.0000000
## L_gougerotae-L_cadamensis   -11.148947 -13.162993 -9.1349018 0.0000000
## L_gougerotae-L_campylolepis  -2.510526  -4.524572 -0.4964808 0.0110701

Trực quan hoá phân tích hậu định

plot(t, xlim =c(-15,0))

Kruskal test (kiểm định phi tham số) Cài đặt packages “pgirmess”

library(pgirmess)

VD: Kiểm định kruskal cho biến aspect ratio (rat)

kruskal.test(rat~species)
## 
##  Kruskal-Wallis rank sum test
## 
## data:  rat by species
## Kruskal-Wallis chi-squared = 38.826, df = 2, p-value = 3.707e-09

Phân tích hậu định

kruskalmc(rat, species)
## Multiple comparison test after Kruskal-Wallis 
## p.value: 0.05 
## Comparisons
##                              obs.dif critical.dif difference
## L_cadamensis-L_campylolepis 14.47368     12.89198       TRUE
## L_cadamensis-L_gougerotae   18.97368     12.89198       TRUE
## L_campylolepis-L_gougerotae 33.44737     12.89198       TRUE

VD: Kiểm định kruskal cho biến aspect ratio (cir)

kruskal.test(cir~species)
## 
##  Kruskal-Wallis rank sum test
## 
## data:  cir by species
## Kruskal-Wallis chi-squared = 39.482, df = 2, p-value = 2.671e-09

Phân tích hậu định

kruskalmc(cir, species)
## Multiple comparison test after Kruskal-Wallis 
## p.value: 0.05 
## Comparisons
##                              obs.dif critical.dif difference
## L_cadamensis-L_campylolepis 13.71053     12.89198       TRUE
## L_cadamensis-L_gougerotae   19.86842     12.89198       TRUE
## L_campylolepis-L_gougerotae 33.57895     12.89198       TRUE

——- end———–

Ngoc Nguyen