#Thực hành hiển thị dữ liệu - chiều ngày 03/01/2024 ## Gọi package

library(ggplot2)
library (gridExtra)
library(ggthemes)

Việc 1. Đọc dữ liệu ob

ob = read.csv ("D:\\Deskop\\Học R tại CT\\DỮ LIỆU ĐÍNH KÈM BÀI TẬP\\Obesity data.csv")
head (ob)
##   id gender height weight  bmi age  bmc  bmd   fat  lean pcfat
## 1  1      F    150     49 21.8  53 1312 0.88 17802 28600  37.3
## 2  2      M    165     52 19.1  65 1309 0.84  8381 40229  16.8
## 3  3      F    157     57 23.1  64 1230 0.84 19221 36057  34.0
## 4  4      F    156     53 21.8  56 1171 0.80 17472 33094  33.8
## 5  5      M    160     51 19.9  54 1681 0.98  7336 40621  14.8
## 6  6      F    153     47 20.1  52 1358 0.91 14904 30068  32.2

Việc 2. Soạn biểu đồ phân bố histogram

grid.arrange: cho phép hiển thị 2 biểu đồ cùng lúc

p = ggplot(data=ob, aes(x=pcfat, fill = gender)) + geom_histogram( col="white")
p1 = p +
  labs (x="Percent body fat (%)", y=" Number of pepple", title = " Phân bố tỷ trọng mỡ") +
  theme_economist() 
grid.arrange(p, p1, ncol = 2)
## `stat_bin()` using `bins = 30`. Pick better value with `binwidth`.
## `stat_bin()` using `bins = 30`. Pick better value with `binwidth`.

Việc 3. Soạn biểu đồ thanh

###3.1. Phân nhóm BMI dựa trên BMI

ob$OB[ob$bmi< 18.5] = "Underweight"
ob$OB[ob$bmi>= 18.5 & ob$bmi< 25] = "Normal"
ob$OB[ob$bmi>= 25 & ob$bmi< 30] = "Overweight"
ob$OB[ob$bmi>= 30] = "Obese"

ob$OB = factor(ob$OB, levels = c("Underweight", "Normal", "Overweight", "Obese"))

Vẽ biểu đồ phân bố BMI

p = ggplot(data=ob, aes(x= OB, y=pcfat, fill = OB), width=0.5) + geom_bar(stat = "identity", just = 0.3, width = 0.5)
p + theme_pander() + theme(legend.position = "none") + labs (x= "Tình trạng cân nặng", y=" Tỷ trọng mỡ", title = "Tỷ trọng mỡ theo tình trạng cân nặng")

###3.2 Đếm số béo phì theo giới tính

library(tidyverse)
## ── Attaching core tidyverse packages ──────────────────────── tidyverse 2.0.0 ──
## ✔ dplyr     1.1.4     ✔ readr     2.1.5
## ✔ forcats   1.0.0     ✔ stringr   1.5.1
## ✔ lubridate 1.9.4     ✔ tibble    3.2.1
## ✔ purrr     1.0.2     ✔ tidyr     1.3.1
## ── Conflicts ────────────────────────────────────────── tidyverse_conflicts() ──
## ✖ dplyr::combine() masks gridExtra::combine()
## ✖ dplyr::filter()  masks stats::filter()
## ✖ dplyr::lag()     masks stats::lag()
## ℹ Use the conflicted package (<http://conflicted.r-lib.org/>) to force all conflicts to become errors
temp = ob %>% group_by(gender) %>% count(OB) %>% mutate(pct = n/sum(n))
temp$pct = round(temp$pct*100, 1)

###3.3 Vẽ biểu đồ ### position=“dodge”: các cột ở cạnh nhau chứ ko chồng lên nhau

p = ggplot(data = temp, aes(x = OB, y = pct, fill = gender, group = gender, width=0.5))
p1 = p + geom_bar(stat = "identity", position = "dodge", width=0.3) + geom_text(aes(x = OB, y = pct, label = pct, group = gender), position = position_dodge(width = 0.5), vjust = -0.25, col = "white") + labs(x = "Obesity status", y = "Percent (%)", title= "Phân bố BMI theo giới tính") + theme(legend.position = "top") +theme_economist() + scale_fill_brewer(palette = "Set1") 
p1

Việc 4. Soạn biểu đồ hộp

p2 = ggplot (data=ob, aes(x=gender, y = pcfat, col = gender)) + geom_boxplot(width=0.3) + geom_jitter(alpha= 0.2, width=0.2)
p2 + theme_tufte() + theme(legend.position = "none") + labs(title = "Phân bố phần trăm mỡ cơ thể theo giới tính", 
       x = "Giới tính", 
       y = "Tỷ trọng mỡ cơ thể (%)", 
       color = "Giới tính")

Việc 5. Soạn biểu đồ tương quan

p3 = ggplot (data=ob, aes(x= bmi, y= pcfat, col=gender)) + geom_point() + geom_smooth() 
p3 + theme_economist() + labs( x= "BMI", y=" Tỷ trọng mỡ (%)", title=" Tương quan tỷ trọng mỡ và BMI") 
## `geom_smooth()` using method = 'loess' and formula = 'y ~ x'

Thêm phân bố

library(ggExtra)
p3 = ggplot (data=ob, aes(x= bmi, y= pcfat, col=gender)) + geom_point() + geom_smooth(method="lm", formula=y~x+ I(x^2))  
p3 + theme_economist() + labs( x= "BMI", y=" Tỷ trọng mỡ (%)", title=" Tương quan tỷ trọng mỡ và BMI") 
ggMarginal(p3, groupColour = T, groupFill = T)

Ridge plot

library (ggridges)
p4 = ggplot(data =ob, aes(x= pcfat, y=OB, fill=gender))+ geom_density_ridges() + theme_ridges() + theme (legend.position = "none")+ labs (x="Percent body fat")
p4
## Picking joint bandwidth of 1.59

Logistic regression plot

library (broom)

Plot of logistic regression

dat = data.frame(y = sample(0:1,1000, replace=TRUE))
dat$x1 = rnorm(1000,5,2) * (dat$y+1)
m=glm (y~x1, family=binomial, data = dat)
plot.dat = augment(m, type.predict = "response")
p5 = ggplot(plot.dat, aes (x=x1)) + geom_line (aes(y=.fitted),color="blue")+ labs (x="Size", y="Survial")+ geom_point(aes(y=y), alpha=0.2)
p5 + stat_summary_bin(geom= "point", fun=mean, aes(y=y), bins=60) 

p5 + geom_histogram(aes(x=1, y= stat(count)/1000), bins=30, na.rm = TRUE, col="white", fill="blue")+ geom_histogram(aes(x= dat$x1, y= -1*stat(count)/1000),bins=30,na.rm= TRUE,position= position_nudge(y=1),col="white",fill="blue")
## Warning: `stat(count)` was deprecated in ggplot2 3.4.0.
## ℹ Please use `after_stat(count)` instead.
## This warning is displayed once every 8 hours.
## Call `lifecycle::last_lifecycle_warnings()` to see where this warning was
## generated.
## Warning in geom_histogram(aes(x = 1, y = stat(count)/1000), bins = 30, na.rm = TRUE, : All aesthetics have length 1, but the data has 1000 rows.
## ℹ Please consider using `annotate()` or provide this layer with data containing
##   a single row.

Phần 2: So sánh 2 nhóm - Biến liên tục

df = read.csv("D:\\Deskop\\Học R tại CT\\DỮ LIỆU ĐÍNH KÈM BÀI TẬP\\Bone data.csv")
head(df)
##   id    sex age weight height prior.fx fnbmd smoking fx
## 1  1   Male  73     98    175        0  1.08       1  0
## 2  2 Female  68     72    166        0  0.97       0  0
## 3  3   Male  68     87    184        0  1.01       0  0
## 4  4 Female  62     72    173        0  0.84       1  0
## 5  5   Male  61     72    173        0  0.81       1  0
## 6  6 Female  76     57    156        0  0.74       0  0

So sánh 2 nhóm