library(readr)
library(dplyr)
library(ggplot2)
library(gtsummary)
library(gt)
library(mice)

Dữ liệu: File Day3_Cardio.csv gồm 500 hồ sơ bệnh nhân.
Mục tiêu: Đánh giá nguy cơ tim mạch nền.

Bài 1

KHỞI TẠO VÀ KHÁM PHÁ DỮ LIỆU

Yêu cầu

Nhập dữ liệu và kiểm tra các giá trị khuyết (NA).

df <- read_csv("Day3_Cardio.csv")
head(df)

## # A tibble: 6 × 10
##   ID      Age Gender Height_cm Weight_kg Smoking Coffee_Cups Knowledge_Score
##   <chr> <dbl> <chr>      <dbl>     <dbl> <chr>         <dbl>           <dbl>
## 1 PT001    60 Male         171        61 Yes               4            53.4
## 2 PT002    44 Male         169        65 No                0            69.9
## 3 PT003    56 Male         170        68 Yes               4            74.4
## 4 PT004    54 Female       154        68 No                0            73.6
## 5 PT005    48 Female       154        56 No                1            57  
## 6 PT006    30 Male         157        77 No                2            54.8
## # ℹ 2 more variables: Triglycerides <dbl>, Systolic_BP <dbl>

summary(df)

##       ID                 Age           Gender            Height_cm    
##  Length:500         Min.   :30.00   Length:500         Min.   :134.0  
##  Class :character   1st Qu.:48.00   Class :character   1st Qu.:154.0  
##  Mode  :character   Median :56.00   Mode  :character   Median :161.0  
##                     Mean   :55.37                      Mean   :161.1  
##                     3rd Qu.:62.00                      3rd Qu.:168.0  
##                     Max.   :86.00                      Max.   :185.0  
##                                                                       
##    Weight_kg       Smoking           Coffee_Cups    Knowledge_Score
##  Min.   :31.00   Length:500         Min.   :0.000   Min.   :28.60  
##  1st Qu.:55.00   Class :character   1st Qu.:0.000   1st Qu.:55.90  
##  Median :63.00   Mode  :character   Median :1.000   Median :64.05  
##  Mean   :63.24                      Mean   :1.494   Mean   :64.23  
##  3rd Qu.:71.00                      3rd Qu.:2.000   3rd Qu.:71.90  
##  Max.   :96.00                      Max.   :8.000   Max.   :99.40  
##                                                                    
##  Triglycerides    Systolic_BP 
##  Min.   :0.470   Min.   :108  
##  1st Qu.:1.100   1st Qu.:141  
##  Median :1.440   Median :150  
##  Mean   :1.579   Mean   :150  
##  3rd Qu.:1.930   3rd Qu.:158  
##  Max.   :4.790   Max.   :194  
##  NA's   :25

[TƯ DUY]

Câu hỏi: Bạn phát hiện biến nào bị khuyết? Theo tiêu chuẩn Q1, bạn có được phép tự ý lấy số Trung bình để điền vào ô khuyết này không? Tại sao?

Trả lời:
Biến Triglycerides bị khuyết (có 25 giá trị NA, chiếm 5%)
Biến Triglycerides là MAR do là vd: những người NA thì hút thuốc nhiều hơn (28% > 16%) do là những người hút thuốc thì ít quan tâm tới sức khỏe hơn nên dễ bị không quan tâm tới làm xét nghiệm lipid máu

df$TG_missing <- ifelse(is.na(df$Triglycerides), "Missing", "Observed")
tab <- table(df$Smoking, df$TG_missing)
prop.table(tab, margin = 2) * 100

##      
##        Missing Observed
##   No  72.00000 83.78947
##   Yes 28.00000 16.21053

=> Không được tự ý lấy trung bình điền vô NA. Mình có thể dùng cách sau sẽ ok hơn.

# Điền giá trị khuyết bằng phương pháp Multiple Imputation (MI) với 50 lần lặp và phương pháp Predictive Mean Matching (PMM)
imp <- mice(df, m = 50, method = "pmm", seed = 123)
df <- complete(imp, 35) # Lấy 35 tại vì đó là sinh nhật em :)))
# df bây giờ sẽ không có NA nữa

Bài 2

CHẨN ĐOÁN “SỨC KHỎE” BIẾN SỐ

Yêu cầu

Vẽ Histogram và Q-Q Plot, Q-Q lines cho Knowledge_Score và Triglycerides.
Thêm đường Mean (màu đỏ).

Knowledge_Score

ggplot(df, aes(x = Knowledge_Score)) +
  geom_histogram(bins = 30, fill = "grey", color = "black") +
  geom_vline(
    xintercept = mean(df$Knowledge_Score, na.rm = TRUE),
    color = "red",
    linewidth = 1
  ) +
  labs(
    x = "Knowledge Score",
    y = "Frequency"
  )

qqnorm(df$Knowledge_Score, main = "Q-Q Plot: Knowledge_Score",
       xlab = "Theoretical Quantiles",
       ylab = "Observed Knowledge Score")
qqline(df$Knowledge_Score, col = "red", lwd = 2)

Triglycerides

ggplot(df, aes(x = Triglycerides)) +
  geom_histogram(bins = 30, fill = "grey", color = "black") +
  geom_vline(
    xintercept = mean(df$Triglycerides, na.rm = TRUE),
    color = "red",
    linewidth = 1
  ) +
  labs(
    x = "Triglycerides",
    y = "Frequency"
  )

qqnorm(df$Triglycerides, 
       main = "Q-Q Plot: Triglycerides",
       xlab = "Theoretical Quantiles",
       ylab = "Observed Triglycerides")
qqline(df$Triglycerides, col = "red", lwd = 2)

[TƯ DUY]

Câu hỏi: Nhìn vào đồ thị Q-Q Plot, các chấm đen văng ra khỏi đường thẳng ở biến nào? Biến đó có tuân theo phân phối chuẩn không?

Trả lời:
Các điểm đen văng xa khỏi đường thẳng chủ yếu xuất hiện ở biến Triglycerides, cho thấy biến này không tuân theo phân phối chuẩn.
Mình có thể làm thêm kiểm định Shapiro-Wilk cho chắc

shapiro.test(df$Triglycerides)

## 
##  Shapiro-Wilk normality test
## 
## data:  df$Triglycerides
## W = 0.91082, p-value < 2.2e-16

p < 0.05 nên là Triglycerides chắc chắn là không phải phân phối chuẩn rồi.

Bài 3

PHẪU THUẬT DỮ LIỆU LỆCH

Yêu cầu

Dùng phép biến đổi Logarit cơ số 10 cho Triglycerides.

df <- df %>% mutate(log10_Triglycerides = log10(Triglycerides))
# Thêm cột log10_Triglycerides vào df

Vẽ lại Q-Q plot.

qqnorm(df$log10_Triglycerides, main = "Q-Q Plot: log10(Triglycerides)",
       xlab = "Theoretical Quantiles",
       ylab = "Observed log10(Triglycerides)")
qqline(df$log10_Triglycerides, col = "red", lwd = 2)

[TƯ DUY]

Câu hỏi: Sau khi “phẫu thuật”, dữ liệu đã đủ điều kiện để chạy các kiểm định tham số (như T-test) chưa? Tại sao?

Trả lời:
Rồi. - Tại vì là trên Q-Q plot thì biến log10_Triglycerides đã bám sát đường thẳng hơn với lại làm thêm kiểm định Shapiro-Wilk thì p > 0.05 rồi nên là nó là phân phối chuẩn

shapiro.test(df$log10_Triglycerides)

## 
##  Shapiro-Wilk normality test
## 
## data:  df$log10_Triglycerides
## W = 0.9963, p-value = 0.2996

log10_Triglycerides là biến liên tục
log10_Triglycerides là biến độc lập => Đủ điều kiện để chạy T-test rồi

Bài 4

SO SÁNH TRUNG BÌNH (T-TEST)

Yêu cầu

So sánh Điểm kiến thức (Knowledge_Score) giữa Nam và Nữ.

table(df$Gender)

## 
## Female   Male 
##    270    230

shapiro.test(df$Knowledge_Score[df$Gender == "Male"])

## 
##  Shapiro-Wilk normality test
## 
## data:  df$Knowledge_Score[df$Gender == "Male"]
## W = 0.99436, p-value = 0.5467

shapiro.test(df$Knowledge_Score[df$Gender == "Female"])

## 
##  Shapiro-Wilk normality test
## 
## data:  df$Knowledge_Score[df$Gender == "Female"]
## W = 0.99505, p-value = 0.5377

# p-value của nam và nữ đều > 0.05 nên đều chuẩn hết; n của cả nam và nữ đều > 30 nên xài t-test 
var.test(Knowledge_Score ~ Gender, data = df)

## 
##  F test to compare two variances
## 
## data:  Knowledge_Score by Gender
## F = 1.0762, num df = 269, denom df = 229, p-value = 0.5669
## alternative hypothesis: true ratio of variances is not equal to 1
## 95 percent confidence interval:
##  0.8372276 1.3797988
## sample estimates:
## ratio of variances 
##            1.07617

# p-value của var.test > 0.05 nên để var.equal = TRUE được
test_gender_knowledge <- t.test(Knowledge_Score ~ Gender, data = df, var.equal = TRUE)
test_gender_knowledge

## 
##  Two Sample t-test
## 
## data:  Knowledge_Score by Gender
## t = 0.32049, df = 498, p-value = 0.7487
## alternative hypothesis: true difference in means between group Female and group Male is not equal to 0
## 95 percent confidence interval:
##  -1.795663  2.495663
## sample estimates:
## mean in group Female   mean in group Male 
##                64.39                64.04

[TƯ DUY]

Câu hỏi: Nhìn vào Khoảng tin cậy 95% (95% CI) của sự chênh lệch. Nếu CI chứa giá trị 0, điều đó có nghĩa là gì? Sự khác biệt này có đủ lớn để bạn đưa ra một chương trình giáo dục sức khỏe riêng biệt cho từng giới tính không?

Trả lời:

test_gender_knowledge$conf.int

## [1] -1.795663  2.495663
## attr(,"conf.level")
## [1] 0.95

95% CI có chứa 0 và p-value > 0.05 => Chỉ kết luận được là không bằng chứng cho thấy 2 giới khác nhau thôi chứ ko kết luận chắc được là 2 giới khác nhau.
Mean của nam và nữ chỉ hơn nhau 0.35 thôi => Quá nhỏ để dù khác biệt thật thì cũng ko có đủ ý nghĩa trong lâm sàng.
=> Sự khác biệt này không đủ lớn để đưa ra một chương trình giáo dục sức khỏe riêng biệt cho từng giới tính

Bài 5

KIỂM ĐỊNH TỶ LỆ (CHI-SQUARE)

Yêu cầu

Tìm mối liên quan giữa Giới tính và Hành vi Hút thuốc (Smoking).

tab_gender_smoking <- table(df$Gender, df$Smoking)
tab_gender_smoking

##         
##           No Yes
##   Female 264   6
##   Male   152  78

chi_gender_smoking <- chisq.test(tab_gender_smoking)
chi_gender_smoking

## 
##  Pearson's Chi-squared test with Yates' continuity correction
## 
## data:  tab_gender_smoking
## X-squared = 86.986, df = 1, p-value < 2.2e-16

chi_gender_smoking$expected

##         
##              No   Yes
##   Female 224.64 45.36
##   Male   191.36 38.64

# Tất cả các ô kỳ vọng đều >5 nên ko cần dùng Fisher

[TƯ DUY]

Câu hỏi: Nếu giá trị P < 0.05, bạn kết luận giới tính có liên quan đến hút thuốc. Nhưng hãy nhìn vào bảng tần suất: Tỷ lệ nam hút thuốc cao hơn nữ bao nhiêu lần? Ý nghĩa thực tế ở đây là gì?

Trả lời:

prop_gender_smoking <- prop.table(tab_gender_smoking, margin = 1)
prop_gender_smoking

##         
##                  No        Yes
##   Female 0.97777778 0.02222222
##   Male   0.66086957 0.33913043

prop_gender_smoking["Male", "Yes"]/prop_gender_smoking["Female", "Yes"]

## [1] 15.26087

Tỷ lệ nam hút thuốc là 33.9%, trong khi tỷ lệ nữ hút thuốc là 2.22% => tỷ lệ hút thuốc ở nam cao hơn khoảng 15.26 lần so với nữ.
Về ý nghĩa thực tế, đây không chỉ là khác biệt có ý nghĩa thống kê mà còn là khác biệt rất lớn về mặt lâm sàng và y tế công cộng => gợi ý rằng các chương trình phòng chống tác hại thuốc lá nên đặc biệt ưu tiên nhóm nam giới. Nhưng mà số lượng nữ hút thuốc cũng khá ít (6 người) nên đây chỉ nên xem là gợi ý thôi, muốn chắc hơn thì em nghĩ phải có nhiều mẫu hơn.

Bài 6

BẪY TƯƠNG QUAN VS SO SÁNH NHÓM (PDTR)

Yêu cầu

Kiểm định mối liên hệ giữa BMI và Triglycerides.

Cách 1: Dùng tương quan

df <- df %>% mutate(BMI = Weight_kg / ((Height_cm/100)^2)) # Thêm cột BMI vào df
cor.test(df$BMI, df$Triglycerides, method = "spearman")

## 
##  Spearman's rank correlation rho
## 
## data:  df$BMI and df$Triglycerides
## S = 19666692, p-value = 0.2113
## alternative hypothesis: true rho is not equal to 0
## sample estimates:
##      rho 
## 0.055995

cor.test(df$BMI, df$log10_Triglycerides, method = "pearson")

## 
##  Pearson's product-moment correlation
## 
## data:  df$BMI and df$log10_Triglycerides
## t = 1.5029, df = 498, p-value = 0.1335
## alternative hypothesis: true correlation is not equal to 0
## 95 percent confidence interval:
##  -0.02061785  0.15397753
## sample estimates:
##        cor 
## 0.06719423

Cách 2: Dùng t-test

df <- df %>%
  mutate(BMI_group = ifelse(BMI < 25, "Normal/Underweight", "Overweight/Obese"))
# do là Age từ 30-86 nên là em dùng mốc 25 của WHO để phân chia giữa bth với thừa cân ở người lớn
table(df$BMI_group)

## 
## Normal/Underweight   Overweight/Obese 
##                286                214

# Triggerides không phân phối chuẩn nên dùng Wilcoxon test
wilcox.test(Triglycerides ~ BMI_group, data = df)

## 
##  Wilcoxon rank sum test with continuity correction
## 
## data:  Triglycerides by BMI_group
## W = 29618, p-value = 0.5384
## alternative hypothesis: true location shift is not equal to 0

# Log10_Triglycerides phân phối chuẩn nên dùng t-test
t.test(log10_Triglycerides ~ BMI_group, data = df)

## 
##  Welch Two Sample t-test
## 
## data:  log10_Triglycerides by BMI_group
## t = -0.77088, df = 460.59, p-value = 0.4412
## alternative hypothesis: true difference in means between group Normal/Underweight and group Overweight/Obese is not equal to 0
## 95 percent confidence interval:
##  -0.04330887  0.01890393
## sample estimates:
## mean in group Normal/Underweight   mean in group Overweight/Obese 
##                        0.1567417                        0.1689442

[TƯ DUY]

Câu hỏi: Nếu bài này muốn biết mối liên hệ thì dùng tương quan hay t-test? Tại sao? Nêu giải pháp bằng 2 cách.

Trả lời:
Nếu mục tiêu là tìm mối liên hệ giữa BMI và Triglycerides, thì về bản chất nên ưu tiên phân tích tương quan hơn là t-test, vì tương quan đánh giá mối liên hệ giữa hai biến liên tục, còn t-test chỉ so sánh trung bình giữa hai nhóm.
Có thể giải quyết theo 2 cách:
- (1) giữ BMI là biến liên tục và dùng Spearman với biến Triglycerides (do ko phải phân phối chuẩn); dùng Pearson với biến log10_Triglycerides (do đây là phân phối chuẩn)
- (2) chia BMI thành 2 nhóm (Normal/Underweight” và “Overweight/Obese) rồi so sánh Triglycerides giữa hai nhóm bằng Wilcoxon (đối với Triglycerides) và t-test (đối với log10_Triglycerides).

Bài 7

ĐÓNG GÓI TABLE 1 (CHUẨN QUỐC TẾ)

Yêu cầu

Sử dụng gtsummary tạo Bảng 1 phân tầng theo Smoking.

num_vars <- df %>% select(where(is.numeric))
lapply(num_vars, shapiro.test)

## $Age
## 
##  Shapiro-Wilk normality test
## 
## data:  X[[i]]
## W = 0.99612, p-value = 0.2609
## 
## 
## $Height_cm
## 
##  Shapiro-Wilk normality test
## 
## data:  X[[i]]
## W = 0.99427, p-value = 0.05741
## 
## 
## $Weight_kg
## 
##  Shapiro-Wilk normality test
## 
## data:  X[[i]]
## W = 0.99605, p-value = 0.2471
## 
## 
## $Coffee_Cups
## 
##  Shapiro-Wilk normality test
## 
## data:  X[[i]]
## W = 0.79008, p-value < 2.2e-16
## 
## 
## $Knowledge_Score
## 
##  Shapiro-Wilk normality test
## 
## data:  X[[i]]
## W = 0.99823, p-value = 0.8943
## 
## 
## $Triglycerides
## 
##  Shapiro-Wilk normality test
## 
## data:  X[[i]]
## W = 0.91082, p-value < 2.2e-16
## 
## 
## $Systolic_BP
## 
##  Shapiro-Wilk normality test
## 
## data:  X[[i]]
## W = 0.99482, p-value = 0.09113
## 
## 
## $log10_Triglycerides
## 
##  Shapiro-Wilk normality test
## 
## data:  X[[i]]
## W = 0.9963, p-value = 0.2996
## 
## 
## $BMI
## 
##  Shapiro-Wilk normality test
## 
## data:  X[[i]]
## W = 0.99466, p-value = 0.07929

# Ta thấy là có mỗi Triglycerides với Coffee_Cups là không chuẩn thôi nên trong bảng sẽ thể hiện bằng median với tứ phân vị; còn lại là mean với sd

df %>%
  select(Age, Gender, BMI, Coffee_Cups, Triglycerides, Smoking, Systolic_BP) %>%
  tbl_summary(
    by = Smoking,
    type = list(Coffee_Cups ~ "continuous"),
    statistic = list(
      Age ~ "{mean} ({sd})",
      Gender ~ "{n} ({p}%)",
      BMI ~ "{mean} ({sd})",
      Coffee_Cups ~ "{median} ({p25}, {p75})",
      Triglycerides ~ "{median} ({p25}, {p75})",
      Systolic_BP ~ "{mean} ({sd})"
    ),
    label = list(
      Age ~ "Age (years)",
      Gender ~ "Gender",
      BMI ~ "BMI (kg/m\u00B2)",
      Coffee_Cups ~ "Coffee intake (cups/day)",
      Triglycerides ~ "Triglycerides (mmol/L)",
      Systolic_BP ~ "Systolic blood pressure (mmHg)"
    ),
    missing_text = "Missing"
  ) %>%
  add_p() %>%
  bold_labels() %>%
  modify_spanning_header(all_stat_cols() ~ "**Smoking status**") %>%
  modify_caption("**Table 1. Baseline characteristics stratified by smoking status**")

**Table 1. Baseline characteristics stratified by smoking status**
Characteristic	Smoking status		p-value²
Characteristic	No N = 416¹	Yes N = 84¹	p-value²
Age (years)	55 (10)	57 (11)	0.2
Gender			<0.001
Female	264 (63%)	6 (7.1%)
Male	152 (37%)	78 (93%)
BMI (kg/m²)	24.4 (4.7)	24.7 (4.5)	0.5
Coffee intake (cups/day)	1.00 (0.00, 2.00)	4.00 (3.00, 6.00)	<0.001
Triglycerides (mmol/L)	1.43 (1.12, 1.89)	1.49 (1.03, 1.99)	>0.9
Systolic blood pressure (mmHg)	147 (11)	165 (11)	<0.001
¹ Mean (SD); n (%); Median (Q1, Q3)
² Wilcoxon rank sum test; Pearson’s Chi-squared test

[TƯ DUY]

Câu hỏi: Viết 3 câu diễn giải kết quả Bảng 1 này như cách bạn sẽ viết trong phần “Results” của một bài báo quốc tế.

Trả lời:
Of the 500 participants, 84 (16.8%) were current smokers and 416 (83.2%) were non-smokers. Compared with non-smokers, smokers were more often male (92.9% vs 36.5%), consumed more coffee (median 4.0 [IQR 3.0–6.0] vs 1.0 [IQR 0.0–2.0] cups/day), and had higher systolic blood pressure (165 ± 11 vs 147 ± 11 mmHg; all p < 0.001). Age, BMI, and triglyceride levels did not differ between the two groups (p = 0.2, 0.5, and >0.9, respectively).

Homework Day 3

Trần Anh Vũ

2026-04-04

Bài 1

Bài 2

Bài 3

Bài 4

Bài 5

Bài 6

Bài 7