Bai tap ngay 2
Tri Nguyen
2026-01-07
Ngày 2: So sánh 2 nhóm
Việc 1. Đoc dữ liệu “Bone data.csv” và gọi dữ liệu này là “df”
df = read.csv("Bone data.csv")
dim(df)## [1] 2162 9head(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 0summary(df)## id sex age weight
## Min. : 1.0 Length:2162 Min. :57.00 Min. : 34.00
## 1st Qu.: 545.2 Class :character 1st Qu.:65.00 1st Qu.: 60.00
## Median :1093.5 Mode :character Median :69.00 Median : 69.00
## Mean :1094.7 Mean :70.65 Mean : 70.14
## 3rd Qu.:1639.8 3rd Qu.:75.00 3rd Qu.: 79.00
## Max. :2216.0 Max. :94.00 Max. :133.00
## NA's :1
## height prior.fx fnbmd smoking
## Min. :136.0 Min. :0.0000 Min. :0.2800 Min. :0.0000
## 1st Qu.:158.0 1st Qu.:0.0000 1st Qu.:0.7300 1st Qu.:0.0000
## Median :164.0 Median :0.0000 Median :0.8200 Median :0.0000
## Mean :164.9 Mean :0.1554 Mean :0.8287 Mean :0.4251
## 3rd Qu.:171.0 3rd Qu.:0.0000 3rd Qu.:0.9300 3rd Qu.:1.0000
## Max. :196.0 Max. :1.0000 Max. :1.5100 Max. :1.0000
## NA's :40
## fx
## Min. :0.0000
## 1st Qu.:0.0000
## Median :0.0000
## Mean :0.2521
## 3rd Qu.:1.0000
## Max. :1.0000
## Việc 2. So sánh mật độ xương giữa nam và nữ
2.1 Histogram đánh giá phân bổ
# Kiểm tra nhanh
str(df)## 'data.frame': 2162 obs. of 9 variables:
## $ id : int 1 2 3 4 5 6 7 8 9 10 ...
## $ sex : chr "Male" "Female" "Male" "Female" ...
## $ age : int 73 68 68 62 61 76 63 64 76 64 ...
## $ weight : int 98 72 87 72 72 57 97 85 48 89 ...
## $ height : int 175 166 184 173 173 156 173 167 153 166 ...
## $ prior.fx: int 0 0 0 0 0 0 0 0 0 0 ...
## $ fnbmd : num 1.08 0.97 1.01 0.84 0.81 0.74 1.01 0.86 0.65 0.98 ...
## $ smoking : int 1 0 0 1 1 0 1 0 0 0 ...
## $ fx : int 0 0 0 0 0 0 0 0 0 0 ...# Vẽ histogram
hist(df$fnbmd,
main = "Histogram of Femoral Neck BMD",
xlab = "Femoral neck BMD (g/cm2)",
col = "blue",
border = "white")
2.2 So sánh mật độ cổ xương đùi giữa nam và nữ
library(lessR); ##
## lessR 4.5 feedback: gerbing@pdx.edu
## --------------------------------------------------------------
## > d <- Read("") Read data file, many formats available, e.g., Excel
## d is the default data frame, data= in analysis routines optional
##
## Many examples of reading, writing, and manipulating data, graphics,
## testing means and proportions, regression, factor analysis,
## customization, forecasting, and aggregation to pivot tables.
## Enter: browseVignettes("lessR")
##
## View lessR updates, now including modern time series forecasting
## and many, new Plotly interactive visualizations output. Most
## visualization functions are now reorganized to three functions:
## Chart(): type="bar", "pie", "radar", "bubble", "treemap", "icicle"
## X(): type="histogram", "density", "vbs" and more
## XY(): type="scatter" for a scatterplot, or "contour", "smooth"
## Most previous function calls still work, such as:
## BarChart(), Histogram, and Plot().
## Enter: news(package="lessR"), or ?Chart, ?X, or ?XY
## There is also Flows() for Sankey flow diagrams, see ?Flows
##
## Interactive data analysis for constructing visualizations.
## Enter: interact()library(table1) ##
## Attaching package: 'table1'## The following object is masked from 'package:lessR':
##
## label## The following objects are masked from 'package:base':
##
## units, units<-ttest(fnbmd ~ sex, data = df)##
## Compare fnbmd across sex with levels Male and Female
## Grouping Variable: sex
## Response Variable: fnbmd
##
##
## ------ Describe ------
##
## fnbmd for sex Male: n.miss = 23, n = 822, mean = 0.910, sd = 0.153
## fnbmd for sex Female: n.miss = 17, n = 1300, mean = 0.778, sd = 0.132
##
## Mean Difference of fnbmd: 0.132
##
## Weighted Average Standard Deviation: 0.141
##
##
## ------ Assumptions ------
##
## Note: These hypothesis tests can perform poorly, and the
## t-test is typically robust to violations of assumptions.
## Use as heuristic guides instead of interpreting literally.
##
## Null hypothesis, for each group, is a normal distribution of fnbmd.
## Group Male: Sample mean assumed normal because n > 30, so no test needed.
## Group Female: Sample mean assumed normal because n > 30, so no test needed.
##
## Null hypothesis is equal variances of fnbmd, homogeneous.
## Variance Ratio test: F = 0.023/0.018 = 1.336, df = 821;1299, p-value = 0.000
## Levene's test, Brown-Forsythe: t = 3.449, df = 2120, p-value = 0.001
##
##
## ------ Infer ------
##
## --- Assume equal population variances of fnbmd for each sex
##
## t-cutoff for 95% range of variation: tcut = 1.961
## Standard Error of Mean Difference: SE = 0.006
##
## Hypothesis Test of 0 Mean Diff: t-value = 21.080, df = 2120, p-value = 0.000
##
## Margin of Error for 95% Confidence Level: 0.012
## 95% Confidence Interval for Mean Difference: 0.120 to 0.144
##
##
## --- Do not assume equal population variances of fnbmd for each sex
##
## t-cutoff: tcut = 1.961
## Standard Error of Mean Difference: SE = 0.006
##
## Hypothesis Test of 0 Mean Diff: t = 20.407, df = 1560.981, p-value = 0.000
##
## Margin of Error for 95% Confidence Level: 0.013
## 95% Confidence Interval for Mean Difference: 0.119 to 0.145
##
##
## ------ Effect Size ------
##
## --- Assume equal population variances of fnbmd for each sex
##
## Standardized Mean Difference of fnbmd, Cohen's d: 0.939
##
##
## ------ Practical Importance ------
##
## Minimum Mean Difference of practical importance: mmd
## Minimum Standardized Mean Difference of practical importance: msmd
## Neither value specified, so no analysis
##
##
## ------ Graphics Smoothing Parameter ------
##
## Density bandwidth for sex Male: 0.044
## Density bandwidth for sex Female: 0.034
# Trung bình và độ lệch chuẩn theo giới
tapply(df$fnbmd, df$sex, mean, na.rm = TRUE)## Female Male
## 0.7775231 0.9096959tapply(df$fnbmd, df$sex, sd, na.rm = TRUE)## Female Male
## 0.1323548 0.1529898boxplot(fnbmd ~ sex,
data = df,
main = "Femoral Neck BMD by Sex",
xlab = "Sex",
ylab = "Femoral neck BMD (g/cm2)",
col = c("pink", "lightblue"))
Việc 3. So sánh khác biệt về chỉ số RER giữa giữa 2 nhóm
3.1 Giả sử bạn có dữ liệu về chỉ số RER của 2 nhóm chứng (placebo) và café như sau:
Nhóm chứng (n= 9): 105, 119, 100, 97, 96, 101, 94, 95, 98 Nhóm café (n= 9): 96, 99, 94, 89, 96, 93, 88, 105, 88 Nhập dữ liệu trên vào R.
placebo <- c(105, 119, 100, 97, 96, 101, 94, 95, 98)
coffee <- c(96, 99, 94, 89, 96, 93, 88, 105, 88)
summary(placebo)## Min. 1st Qu. Median Mean 3rd Qu. Max.
## 94.0 96.0 98.0 100.6 101.0 119.0summary(coffee)## Min. 1st Qu. Median Mean 3rd Qu. Max.
## 88.00 89.00 94.00 94.22 96.00 105.003.2 Thực hiện phép kiểm t để đánh giá khác biệt về RER giữa 2 nhóm. Nhận xét kết quả.
# Kiểm t hai mẫu độc lập
t_test_result <- t.test(placebo, coffee,
alternative = "two.sided",
var.equal = TRUE)
t_test_result##
## Two Sample t-test
##
## data: placebo and coffee
## t = 1.9948, df = 16, p-value = 0.06339
## alternative hypothesis: true difference in means is not equal to 0
## 95 percent confidence interval:
## -0.3972398 13.0639064
## sample estimates:
## mean of x mean of y
## 100.55556 94.22222Nhận xét kết quả phép kiểm t so sánh RER giữa 2 nhóm
Phép kiểm t hai mẫu độc lập được sử dụng để so sánh chỉ số RER trung bình giữa nhóm chứng (placebo) và nhóm sử dụng café. Kết quả cho thấy giá trị trung bình RER của hai nhóm có sự khác biệt. Tuy nhiên, sự khác biệt này (có/không) đạt ý nghĩa thống kê với mức ý nghĩa 5% (p = …).
Điều này cho thấy việc sử dụng café (có/không) làm thay đổi đáng kể chỉ số RER so với nhóm chứng trong mẫu nghiên cứu.
Nếu p < 0.05 → “có ý nghĩa thống kê”
Nếu p ≥ 0.05 → “không có ý nghĩa thống kê”
3.3 Thực hiện bootstrap để đánh giá khác biệt về RER giữa 2 nhóm. Nhận xét kết quả.
library(simpleboot)## Simple Bootstrap Routines (1.1-8)library(boot)
set.seed(123) # Để kết quả tái lập được
boot_result <- two.boot(placebo, coffee, FUN = mean, R = 10000)
boot.ci(boot_result, type = c("perc","bca","norm", "basic"))## BOOTSTRAP CONFIDENCE INTERVAL CALCULATIONS
## Based on 10000 bootstrap replicates
##
## CALL :
## boot.ci(boot.out = boot_result, type = c("perc", "bca", "norm",
## "basic"))
##
## Intervals :
## Level Normal Basic
## 95% ( 0.383, 12.188 ) ( 0.111, 11.778 )
##
## Level Percentile BCa
## 95% ( 0.889, 12.556 ) ( 0.778, 12.444 )
## Calculations and Intervals on Original ScaleNhận xét kết quả bootstrap so sánh RER giữa 2 nhóm
Phương pháp bootstrap với 10.000 lần lặp được sử dụng để ước lượng khoảng tin cậy 95% của chênh lệch trung bình RER giữa nhóm placebo và nhóm caffee. Kết quả cho thấy khoảng tin cậy bootstrap chênh lệch trung bình là [0.889 ; 12.556].
Do khoảng tin cậy không chứa giá trị 0, nên có thể kết luận rằng sự khác biệt RER giữa hai nhóm có ý nghĩa thống kê. Kết quả bootstrap phù hợp với kết quả thu được từ phép kiểm t, qua đó làm tăng độ tin cậy của kết luận.
Bootstrap không yêu cầu giả định phân phối chuẩn, do đó kết quả có tính vững (robust) hơn so với phương pháp kiểm định tham số truyền thống.
set.seed(123) # Để kết quả tái lập được
# Số lần bootstrap
B <- 10000
# Chênh lệch trung bình quan sát
obs_diff <- mean(placebo) - mean(coffee)
# Bootstrap
boot_diff <- replicate(B, {
mean(sample(placebo, replace = TRUE)) -
mean(sample(coffee, replace = TRUE))
})
# Khoảng tin cậy 95% bootstrap
boot_ci <- quantile(boot_diff, probs = c(0.025, 0.975))
obs_diff## [1] 6.333333boot_ci## 2.5% 97.5%
## 0.8888889 12.4444444Việc 4. So sánh khác biệt về tỷ lệ gãy xương giữa nam và nữ:
4.1 So sánh tỉ lệ gãy xương (fx) giữa nam và nữ (sex).
Đặc điểm Nam (n= xxx) Nữ (n= xxx) Giá trị P Gãy xương Ca gãy xương (%) Ca gãy xương (%)
library(table1)
table1(~ fx + as.factor(fx) | sex, data = df)| Female (N=1317) |
Male (N=845) |
Overall (N=2162) |
|
|---|---|---|---|
| fx | |||
| Mean (SD) | 0.304 (0.460) | 0.170 (0.376) | 0.252 (0.434) |
| Median [Min, Max] | 0 [0, 1.00] | 0 [0, 1.00] | 0 [0, 1.00] |
| as.factor(fx) | |||
| 0 | 916 (69.6%) | 701 (83.0%) | 1617 (74.8%) |
| 1 | 401 (30.4%) | 144 (17.0%) | 545 (25.2%) |
chisq.test(df$fx, df$sex)##
## Pearson's Chi-squared test with Yates' continuity correction
##
## data: df$fx and df$sex
## X-squared = 48.363, df = 1, p-value = 3.542e-12Nhận xét kết quả.
Trong tổng số 2.162 đối tượng nghiên cứu, tỉ lệ gãy xương chung là 25,2%. Ở nhóm nữ (N = 1.317), có 401 trường hợp gãy xương, chiếm 30,4%, trong khi ở nhóm nam (N = 845), có 144 trường hợp gãy xương, tương ứng 17,0%. Như vậy, tỉ lệ gãy xương ở nữ cao hơn rõ rệt so với nam.
Giá trị trung bình của biến gãy xương (fx) ở nữ là 0,304 (SD = 0,460), cao hơn so với nam là 0,170 (SD = 0,376). Điều này cho thấy nguy cơ gãy xương trung bình ở nữ cao hơn so với nam trong quần thể nghiên cứu.
📌 (Vì fx là biến nhị phân 0/1 nên mean chính là tỉ lệ)
Kết quả cho thấy nữ giới có tỉ lệ gãy xương cao hơn đáng kể so với nam giới trong thời gian theo dõi. Điều này gợi ý rằng giới tính có thể là một yếu tố liên quan đến nguy cơ gãy xương.
# Bảng chéo --> phuong an cua Chat GPT
tab_fx_sex <- table(df$sex, df$fx)
tab_fx_sex##
## 0 1
## Female 916 401
## Male 701 144# Tỉ lệ theo hàng (%)
prop.table(tab_fx_sex, margin = 1) * 100##
## 0 1
## Female 69.55201 30.44799
## Male 82.95858 17.04142chisq_test <- chisq.test(tab_fx_sex)
chisq_test##
## Pearson's Chi-squared test with Yates' continuity correction
##
## data: tab_fx_sex
## X-squared = 48.363, df = 1, p-value = 3.542e-12