Chương trình tập huấn phân tích dữ liệu bằng ngôn ngữ R - BV 108

Ngày 2: Hiển thị dữ liệu

Việc 1. Đọc dữ liệu vào R

ob = read.csv("C:\\Thach\\VN trips\\2024_2Aug\\Data Analysis workshop (Hospital 108)\\Datasets\\obesity data.csv")

Việc 2. Mô tả đặc điểm mẫu nghiên cứu:

library(table1)
## 
## Attaching package: 'table1'
## The following objects are masked from 'package:base':
## 
##     units, units<-
table1(~ age + gender + bmi + WBBMC + wbbmd + fat + lean + hypertension + diabetes, data = ob)
Overall
(N=1217)
age
Mean (SD) 47.2 (17.3)
Median [Min, Max] 48.0 [13.0, 88.0]
gender
F 862 (70.8%)
M 355 (29.2%)
bmi
Mean (SD) 22.4 (3.06)
Median [Min, Max] 22.2 [14.5, 37.1]
WBBMC
Mean (SD) 1720 (363)
Median [Min, Max] 1710 [695, 3040]
wbbmd
Mean (SD) 1.01 (0.113)
Median [Min, Max] 1.01 [0.650, 1.35]
fat
Mean (SD) 17300 (5210)
Median [Min, Max] 17000 [4280, 40800]
lean
Mean (SD) 35500 (7030)
Median [Min, Max] 33600 [19100, 63100]
hypertension
Mean (SD) 0.507 (0.500)
Median [Min, Max] 1.00 [0, 1.00]
diabetes
Mean (SD) 0.111 (0.314)
Median [Min, Max] 0 [0, 1.00]

Việc 3. Mô tả đặc điểm mẫu nghiên cứu theo giới tính:

table1(~ age + bmi + WBBMC + wbbmd + fat + lean + as.factor(hypertension) + as.factor(diabetes) | gender, data = ob)
F
(N=862)
M
(N=355)
Overall
(N=1217)
age
Mean (SD) 48.6 (16.4) 43.7 (18.8) 47.2 (17.3)
Median [Min, Max] 49.0 [14.0, 85.0] 44.0 [13.0, 88.0] 48.0 [13.0, 88.0]
bmi
Mean (SD) 22.3 (3.05) 22.7 (3.04) 22.4 (3.06)
Median [Min, Max] 22.1 [15.2, 37.1] 22.5 [14.5, 34.7] 22.2 [14.5, 37.1]
WBBMC
Mean (SD) 1600 (293) 2030 (336) 1720 (363)
Median [Min, Max] 1610 [695, 2660] 2030 [1190, 3040] 1710 [695, 3040]
wbbmd
Mean (SD) 0.988 (0.111) 1.06 (0.101) 1.01 (0.113)
Median [Min, Max] 0.990 [0.650, 1.35] 1.06 [0.780, 1.34] 1.01 [0.650, 1.35]
fat
Mean (SD) 18200 (4950) 15000 (5110) 17300 (5210)
Median [Min, Max] 17700 [6220, 40800] 15100 [4280, 29900] 17000 [4280, 40800]
lean
Mean (SD) 32000 (3970) 43800 (5820) 35500 (7030)
Median [Min, Max] 31500 [19100, 53400] 43400 [28600, 63100] 33600 [19100, 63100]
as.factor(hypertension)
0 430 (49.9%) 170 (47.9%) 600 (49.3%)
1 432 (50.1%) 185 (52.1%) 617 (50.7%)
as.factor(diabetes)
0 760 (88.2%) 322 (90.7%) 1082 (88.9%)
1 102 (11.8%) 33 (9.3%) 135 (11.1%)

Việc 4. So sánh đặc điểm mẫu nghiên cứu giữa nam và nữ:

library(compareGroups)
createTable(compareGroups(gender ~ age + bmi + WBBMC + wbbmd + fat + lean + hypertension + diabetes, data = ob))
## 
## --------Summary descriptives table by 'gender'---------
## 
## ________________________________________________ 
##                   F            M       p.overall 
##                 N=862        N=355               
## ¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯ 
## age          48.6 (16.4)  43.7 (18.8)   <0.001   
## bmi          22.3 (3.05)  22.7 (3.04)    0.013   
## WBBMC         1599 (293)   2030 (336)   <0.001   
## wbbmd        0.99 (0.11)  1.06 (0.10)   <0.001   
## fat          18240 (4954) 14978 (5113)  <0.001   
## lean         32045 (3966) 43762 (5819)  <0.001   
## hypertension 0.50 (0.50)  0.52 (0.50)    0.527   
## diabetes     0.12 (0.32)  0.09 (0.29)    0.181   
## ¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯
ob$hypert = as.factor(ob$hypertension)
ob$dm = as.factor(ob$diabetes)
createTable(compareGroups(gender ~ age + bmi + WBBMC + wbbmd + fat + lean + hypert + dm, data = ob))
## 
## --------Summary descriptives table by 'gender'---------
## 
## ___________________________________________ 
##              F            M       p.overall 
##            N=862        N=355               
## ¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯ 
## age     48.6 (16.4)  43.7 (18.8)   <0.001   
## bmi     22.3 (3.05)  22.7 (3.04)    0.013   
## WBBMC    1599 (293)   2030 (336)   <0.001   
## wbbmd   0.99 (0.11)  1.06 (0.10)   <0.001   
## fat     18240 (4954) 14978 (5113)  <0.001   
## lean    32045 (3966) 43762 (5819)  <0.001   
## hypert:                             0.569   
##     0   430 (49.9%)  170 (47.9%)            
##     1   432 (50.1%)  185 (52.1%)            
## dm:                                 0.238   
##     0   760 (88.2%)  322 (90.7%)            
##     1   102 (11.8%)   33 (9.30%)            
## ¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯

Việc 4. Biểu đồ histogram

4.1 Phân bố tỉ trọng mỡ

library(ggplot2)
library(gridExtra) 

p = ggplot(data = ob, aes(x = pcfat))
p1 = p + geom_histogram()
p2 = p + geom_histogram(fill = "blue", col = "white") + labs(x = "Tỉ trọng mỡ (%)", y = "Số người", title = "Phân bố tỉ trọng mỡ")

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`.

4.2 Phân bố tỉ trọng mỡ theo giới tính

p = ggplot(data = ob, aes(x = pcfat, fill = gender))
p1 = p + geom_histogram(col="white") + labs(x = "Tỉ trọng mỡ", y = "Số người", title = "Phân bố tỉ trọng mỡ")
p2 = p + geom_density(alpha = 0.5) + labs(x = "Tỉ trọng mỡ", y = "Số người", title = "Phân bố tỉ trọng mỡ")

grid.arrange(p1, p2, ncol = 2)
## `stat_bin()` using `bins = 30`. Pick better value with `binwidth`.

Việc 5. Biểu đồ thanh

5.1 Tạo biến số OB từ biế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"

5.2 Phân bố của tình trạng béo phì

p = ggplot(data = ob, aes(x = OB, fill = OB, col = OB))
p + geom_bar(position = "dodge")

5.3 Phân bố của tình trạng béo phì theo giới tính

p = ggplot(data = ob, aes(x = OB, y = pcfat, fill = gender, group = gender))
p + geom_bar(stat = "identity", position = "dodge")

5.4 Thêm tỉ lệ %

Tính tỉ lệ %:

library(tidyverse)
## ── Attaching core tidyverse packages ──────────────────────── tidyverse 2.0.0 ──
## ✔ dplyr     1.1.3     ✔ readr     2.1.4
## ✔ forcats   1.0.0     ✔ stringr   1.5.0
## ✔ lubridate 1.9.3     ✔ tibble    3.2.1
## ✔ purrr     1.0.2     ✔ tidyr     1.3.0
## ── 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)

Thêm % vào biểu đồ

p = ggplot(data = temp, aes(x = OB, y = pct, fill = gender, group = gender))
p1 = p + geom_bar(stat = "identity", position = "dodge") + geom_text(aes(x = OB, y = pct, label = pct, group = gender), position = position_dodge(width = 1), vjust = -0.5, col = "blue")
p1

p = ggplot(data = temp, aes(x = OB, y = pct, fill = gender, label = pct))
p2 = p + geom_bar(stat = "identity") + geom_text(size = 3, position = position_stack(vjust = 0.5))
p2

grid.arrange(p1, p2, ncol = 2)

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

6.1 So sánh phân bố của tỉ trọng mỡ giữa các nhóm béo phì ở nữ

Tạo tập dữ liệu chỉ gồm nữ

women = subset(ob, gender == "F")
dim(women)
## [1] 862  16

Vẽ biểu đồ hộp

p = ggplot(data = women, aes(x = OB, y = pcfat, fill = OB, col = OB))
p1 = p + geom_boxplot(col = "black") + geom_jitter(alpha = 0.05) + labs(x = "Tình trạng béo phì", y = "Tỉ trọng mỡ (%)") + ggtitle("Tỉ trọng mỡ theo tình trạng béo phì ở nữ")
p1

women$OB.n = factor(women$OB, levels = c("Underweight", "Normal", "Overweight", "Obese"))
p = ggplot(data = women, aes(x = OB.n, y = pcfat, fill = OB.n, col = OB.n))
p2 = p + geom_boxplot(col = "black") + geom_jitter(alpha = 0.05) + labs(x = "Tình trạng béo phì", y = "Tỉ trọng mỡ (%)") + ggtitle("Tỉ trọng mỡ theo tình trạng béo phì ở nữ")
p2

grid.arrange(p1, p2, ncol = 2)

6.2 So sánh phân bố của tỉ trọng mỡ giữa các nhóm béo phì ở nam

Tạo tập dữ liệu chỉ gồm nam

men = subset(ob, gender == "M")
dim(men)
## [1] 355  16

Vẽ biểu đồ hộp

p = ggplot(data = men, aes(x = OB, y = pcfat, fill = OB, col = OB))
p1 = p + geom_boxplot(col = "black") + geom_jitter(alpha = 0.05) + labs(x = "Tình trạng béo phì", y = "Tỉ trọng mỡ (%)") + ggtitle("Tỉ trọng mỡ theo tình trạng béo phì ở nam")
p1

men$OB.n = factor(men$OB, levels = c("Underweight", "Normal", "Overweight", "Obese"))
p = ggplot(data = men, aes(x = OB.n, y = pcfat, fill = OB.n, col = OB.n))
p2 = p + geom_boxplot(col = "black") + geom_jitter(alpha = 0.05) + labs(x = "Tình trạng béo phì", y = "Tỉ trọng mỡ (%)") + ggtitle("Tỉ trọng mỡ theo tình trạng béo phì ở nam")
p2

grid.arrange(p1, p2, ncol = 2)

6.3 Trình bày trong cùng 1 biểu đồ

ob$OB.n = factor(ob$OB, levels = c("Underweight", "Normal", "Overweight", "Obese"))
p = ggplot(data = ob, aes(x = OB.n, y = pcfat, fill = gender, col = gender))
p1 = p + geom_boxplot(col = "black") + geom_jitter(alpha = 0.05) + labs(x = "Tình trạng béo phì", y = "Tỉ trọng mỡ (%)") + ggtitle("Tỉ trọng mỡ theo tình trạng béo phì và giới tính")
p1

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

7.1 Mối liên quan giữa chỉ số khối cơ thể và tỉ trọng mỡ

p = ggplot(data = ob, aes(x = bmi, y = pcfat))
p1 = p + geom_point()
p2 = p + geom_point() + geom_smooth()

grid.arrange(p1, p2, ncol = 2)
## `geom_smooth()` using method = 'gam' and formula = 'y ~ s(x, bs = "cs")'

7.2 Mối liên quan giữa chỉ số khối cơ thể và tỉ trọng mỡ theo giới tính

p = ggplot(data = ob, aes(x = bmi, y = pcfat, fill = gender, col = gender))
p1 = p + geom_point() + geom_smooth() + labs(x = "Chỉ số khối cơ thể (kg/m2)", y = "Tỉ trọng mỡ (%)") + ggtitle("Liên quan giữa chỉ số khối cơ thể và tỉ trọng mỡ theo giới tính")
p1
## `geom_smooth()` using method = 'loess' and formula = 'y ~ x'

p = ggplot(data = ob, aes(x = bmi, y = pcfat, fill = gender, col = gender))
p2 = p + geom_point() + geom_smooth(method = "lm", formula = y ~ x + I(x^2)) + labs(x = "Chỉ số khối cơ thể (kg/m2)", y = "Tỉ trọng mỡ (%)") + ggtitle("Liên quan giữa chỉ số khối cơ thể và tỉ trọng mỡ theo giới tính")
p2

grid.arrange(p1, p2, ncol = 2)
## `geom_smooth()` using method = 'loess' and formula = 'y ~ x'

Việc 8. Ghi lại tất cả các hàm/lệnh trên và chia sẻ lên tài khoản rpubs (https://rpubs.com/ThachTran/1213967)