Tap huan Phan tich du lieu
Duy Nguyen
4/29/2020
Đọc dữ liệu
t = "C:\\Users\\Admin\\Desktop\\Data thuc hanh.csv"
fx = read.csv(t)
head(fx)## id gioitinh tuoi cannang chieucao gayxuong matdoxuong hutthuoc
## 1 1 Nam 73 98 175 Co 1.08 Khong
## 2 2 Nam 68 72 166 Khong 0.97 Co
## 3 3 Nam 68 87 184 Khong 1.01 Co
## 4 4 Nu 62 72 173 Khong 0.84 Khong
## 5 5 Nam 61 72 173 Khong 0.81 Khong
## 6 6 Nu 76 57 156 Khong 0.74 CoTạo biến định lượng mới
fx$bmi = fx$cannang/fx$chieucao/fx$chieucao*100*100
head(fx)## id gioitinh tuoi cannang chieucao gayxuong matdoxuong hutthuoc bmi
## 1 1 Nam 73 98 175 Co 1.08 Khong 32.00000
## 2 2 Nam 68 72 166 Khong 0.97 Co 26.12861
## 3 3 Nam 68 87 184 Khong 1.01 Co 25.69707
## 4 4 Nu 62 72 173 Khong 0.84 Khong 24.05693
## 5 5 Nam 61 72 173 Khong 0.81 Khong 24.05693
## 6 6 Nu 76 57 156 Khong 0.74 Co 23.42209Tạo biến định tính từ biến định lượng
fx$bmi_pl[fx$bmi <= 18] <- "Gay"
fx$bmi_pl[fx$bmi > 18 & fx$bmi < 25] <- "Binh thuong"
fx$bmi_pl[fx$bmi >= 25] <- "Beo phi"
head(fx)## id gioitinh tuoi cannang chieucao gayxuong matdoxuong hutthuoc bmi
## 1 1 Nam 73 98 175 Co 1.08 Khong 32.00000
## 2 2 Nam 68 72 166 Khong 0.97 Co 26.12861
## 3 3 Nam 68 87 184 Khong 1.01 Co 25.69707
## 4 4 Nu 62 72 173 Khong 0.84 Khong 24.05693
## 5 5 Nam 61 72 173 Khong 0.81 Khong 24.05693
## 6 6 Nu 76 57 156 Khong 0.74 Co 23.42209
## bmi_pl
## 1 Beo phi
## 2 Beo phi
## 3 Beo phi
## 4 Binh thuong
## 5 Binh thuong
## 6 Binh thuongTạo biến định tính từ biến định lượng: ifelse
fx$beophi = ifelse(fx$bmi_pl == "Beo phi", "Co", "Khong")
head(fx)## id gioitinh tuoi cannang chieucao gayxuong matdoxuong hutthuoc bmi
## 1 1 Nam 73 98 175 Co 1.08 Khong 32.00000
## 2 2 Nam 68 72 166 Khong 0.97 Co 26.12861
## 3 3 Nam 68 87 184 Khong 1.01 Co 25.69707
## 4 4 Nu 62 72 173 Khong 0.84 Khong 24.05693
## 5 5 Nam 61 72 173 Khong 0.81 Khong 24.05693
## 6 6 Nu 76 57 156 Khong 0.74 Co 23.42209
## bmi_pl beophi
## 1 Beo phi Co
## 2 Beo phi Co
## 3 Beo phi Co
## 4 Binh thuong Khong
## 5 Binh thuong Khong
## 6 Binh thuong KhongTạo tập con
nam = subset(fx, gioitinh == "Nam")
nu = subset(fx, gioitinh == "Nu")
head(nam)## id gioitinh tuoi cannang chieucao gayxuong matdoxuong hutthuoc bmi
## 1 1 Nam 73 98 175 Co 1.08 Khong 32.00000
## 2 2 Nam 68 72 166 Khong 0.97 Co 26.12861
## 3 3 Nam 68 87 184 Khong 1.01 Co 25.69707
## 5 5 Nam 61 72 173 Khong 0.81 Khong 24.05693
## 7 7 Nam 63 97 173 Khong 1.01 Khong 32.41004
## 14 14 Nam 62 97 171 Khong 1.16 Co 33.17260
## bmi_pl beophi
## 1 Beo phi Co
## 2 Beo phi Co
## 3 Beo phi Co
## 5 Binh thuong Khong
## 7 Beo phi Co
## 14 Beo phi CoPhân tích mô tả 1 biến định tính
#Sau khi hoàn thiện file xử lý, có thể dùng lệnh attach() để hạn chế việc dùng "$"
attach(fx)
#Gọi package "DescTools)
library(DescTools)
Desc(bmi_pl)## ------------------------------------------------------------------------------
## bmi_pl (character)
##
## length n NAs unique levels dupes
## 2'200 2'145 55 3 3 y
## 97.5% 2.5%
##
## level freq perc cumfreq cumperc
## 1 Beo phi 1'139 53.1% 1'139 53.1%
## 2 Binh thuong 964 44.9% 2'103 98.0%
## 3 Gay 42 2.0% 2'145 100.0%
Phân tích mô tả 2 biến định tính
#Mô tả biến “gayxuong” theo “hutthuoc”
Desc(gayxuong ~ hutthuoc)## ------------------------------------------------------------------------------
## gayxuong ~ hutthuoc
##
##
## Summary:
## n: 2'199, rows: 2, columns: 2
##
## Pearson's Chi-squared test (cont. adj):
## X-squared = 3.0376, df = 1, p-value = 0.08135
## Fisher's exact test p-value = 0.07564
## McNemar's chi-squared = 436.26, df = 1, p-value < 2.2e-16
##
## estimate lwr.ci upr.ci'
##
## odds ratio 1.195 0.983 1.453
## rel. risk (col1) 1.076 0.995 1.163
## rel. risk (col2) 0.900 0.800 1.012
##
##
## Phi-Coefficient 0.038
## Contingency Coeff. 0.038
## Cramer's V 0.038
##
##
## hutthuoc Co Khong Sum
## gayxuong
##
## Co freq 348 220 568
## perc 15.8% 10.0% 25.8%
## p.row 61.3% 38.7% .
## p.col 27.3% 23.9% .
##
## Khong freq 929 702 1'631
## perc 42.2% 31.9% 74.2%
## p.row 57.0% 43.0% .
## p.col 72.7% 76.1% .
##
## Sum freq 1'277 922 2'199
## perc 58.1% 41.9% 100.0%
## p.row . . .
## p.col . . .
##
##
## ----------
## ' 95% conf. level
Phân tích mô tả 1 biến định lượng
#Kiểm tra phân bố chuẩn bằng kiểm định
shapiro.test(matdoxuong)##
## Shapiro-Wilk normality test
##
## data: matdoxuong
## W = 0.99372, p-value = 7.891e-08#Chuyển p-value về giá trị thập phân
options(scipen = 999)shapiro.test(matdoxuong)##
## Shapiro-Wilk normality test
##
## data: matdoxuong
## W = 0.99372, p-value = 0.00000007891#Kiểm tra phân bố chuẩn bằng biểu đồ
Desc(matdoxuong)## ------------------------------------------------------------------------------
## matdoxuong (numeric)
##
## length n NAs unique 0s mean meanCI
## 2'200 2'111 89 99 0 0.83 0.82
## 96.0% 4.0% 0.0% 0.84
##
## .05 .10 .25 median .75 .90 .95
## 0.59 0.64 0.73 0.82 0.93 1.03 1.08
##
## range sd vcoef mad IQR skew kurt
## 1.23 0.15 0.19 0.15 0.20 0.28 0.58
##
## lowest : 0.28, 0.34, 0.37 (2), 0.38, 0.39 (2)
## highest: 1.38, 1.39 (2), 1.41, 1.47, 1.51
#Kiểm tra phân bố chuẩn bằng biểu đồ
qqnorm(matdoxuong)
qqline(matdoxuong)
#Mô tả biến “matdoxuong”
Desc(matdoxuong)## ------------------------------------------------------------------------------
## matdoxuong (numeric)
##
## length n NAs unique 0s mean meanCI
## 2'200 2'111 89 99 0 0.83 0.82
## 96.0% 4.0% 0.0% 0.84
##
## .05 .10 .25 median .75 .90 .95
## 0.59 0.64 0.73 0.82 0.93 1.03 1.08
##
## range sd vcoef mad IQR skew kurt
## 1.23 0.15 0.19 0.15 0.20 0.28 0.58
##
## lowest : 0.28, 0.34, 0.37 (2), 0.38, 0.39 (2)
## highest: 1.38, 1.39 (2), 1.41, 1.47, 1.51
Phân tích mô tả biến định lượng - định tính
Desc(matdoxuong ~ gioitinh)## ------------------------------------------------------------------------------
## matdoxuong ~ gioitinh
##
## Summary:
## n pairs: 2'200, valid: 2'111 (96.0%), missings: 89 (4.0%), groups: 2
##
##
## Nam Nu
## mean 0.909 0.778
## median 0.900 0.775
## sd 0.153 0.132
## IQR 0.185 0.170
## n 823 1'288
## np 38.986% 61.014%
## NAs 32 57
## 0s 0 0
##
## Kruskal-Wallis rank sum test:
## Kruskal-Wallis chi-squared = 366.99, df = 1, p-value < 0.00000000000000022
Sử dụng package table1
library(table1)## Warning: package 'table1' was built under R version 3.6.3##
## Attaching package: 'table1'## The following objects are masked from 'package:base':
##
## units, units<-table1(~tuoi + matdoxuong + beophi |gioitinh, data = fx)| Nam (N=855) |
Nu (N=1345) |
Overall (N=2200) |
|
|---|---|---|---|
| tuoi | |||
| Mean (SD) | 70.4 (6.44) | 71.1 (7.58) | 70.8 (7.17) |
| Median [Min, Max] | 69.0 [59.0, 92.0] | 70.0 [57.0, 96.0] | 70.0 [57.0, 96.0] |
| matdoxuong | |||
| Mean (SD) | 0.909 (0.153) | 0.778 (0.132) | 0.829 (0.155) |
| Median [Min, Max] | 0.900 [0.340, 1.51] | 0.775 [0.280, 1.31] | 0.820 [0.280, 1.51] |
| Missing | 32 (3.7%) | 57 (4.2%) | 89 (4.0%) |
| beophi | |||
| Co | 506 (59.2%) | 633 (47.1%) | 1139 (51.8%) |
| Khong | 336 (39.3%) | 670 (49.8%) | 1006 (45.7%) |
| Missing | 13 (1.5%) | 42 (3.1%) | 55 (2.5%) |
So sánh tỷ lệ 2 nhóm độc lập
#So sánh tỷ lệ gãy xương ở nhóm hút thuốc/không hút thuốc
Desc(gayxuong ~ hutthuoc)## ------------------------------------------------------------------------------
## gayxuong ~ hutthuoc
##
##
## Summary:
## n: 2'199, rows: 2, columns: 2
##
## Pearson's Chi-squared test (cont. adj):
## X-squared = 3.0376, df = 1, p-value = 0.08135
## Fisher's exact test p-value = 0.07564
## McNemar's chi-squared = 436.26, df = 1, p-value < 0.00000000000000022
##
## estimate lwr.ci upr.ci'
##
## odds ratio 1.195 0.983 1.453
## rel. risk (col1) 1.076 0.995 1.163
## rel. risk (col2) 0.900 0.800 1.012
##
##
## Phi-Coefficient 0.038
## Contingency Coeff. 0.038
## Cramer's V 0.038
##
##
## hutthuoc Co Khong Sum
## gayxuong
##
## Co freq 348 220 568
## perc 15.8% 10.0% 25.8%
## p.row 61.3% 38.7% .
## p.col 27.3% 23.9% .
##
## Khong freq 929 702 1'631
## perc 42.2% 31.9% 74.2%
## p.row 57.0% 43.0% .
## p.col 72.7% 76.1% .
##
## Sum freq 1'277 922 2'199
## perc 58.1% 41.9% 100.0%
## p.row . . .
## p.col . . .
##
##
## ----------
## ' 95% conf. level
So sánh tỷ lệ > 2 nhóm độc lập
#So sánh tỷ lệ gãy xương ở các phân nhóm BMI
Desc(gayxuong ~ bmi_pl)## ------------------------------------------------------------------------------
## gayxuong ~ bmi_pl
##
##
## Summary:
## n: 2'145, rows: 2, columns: 3
##
## Pearson's Chi-squared test:
## X-squared = 16.623, df = 2, p-value = 0.0002457
## Likelihood Ratio:
## X-squared = 16.159, df = 2, p-value = 0.0003098
## Mantel-Haenszel Chi-squared:
## X-squared = 15.971, df = 1, p-value = 0.00006431
##
## Phi-Coefficient 0.088
## Contingency Coeff. 0.088
## Cramer's V 0.088
##
##
## bmi_pl Beo phi Binh thuong Gay Sum
## gayxuong
##
## Co freq 249 272 17 538
## perc 11.6% 12.7% 0.8% 25.1%
## p.row 46.3% 50.6% 3.2% .
## p.col 21.9% 28.2% 40.5% .
##
## Khong freq 890 692 25 1'607
## perc 41.5% 32.3% 1.2% 74.9%
## p.row 55.4% 43.1% 1.6% .
## p.col 78.1% 71.8% 59.5% .
##
## Sum freq 1'139 964 42 2'145
## perc 53.1% 44.9% 2.0% 100.0%
## p.row . . . .
## p.col . . . .
## 
Phân tích hậu kiểm (post-hoc)
#Tạo bảng 2x3
bang = table(bmi_pl, gayxuong)
bang## gayxuong
## bmi_pl Co Khong
## Beo phi 249 890
## Binh thuong 272 692
## Gay 17 25library(rcompanion)## Warning: package 'rcompanion' was built under R version 3.6.3pairwiseNominalIndependence(bang)## Comparison p.Fisher p.adj.Fisher p.Gtest p.adj.Gtest p.Chisq
## 1 Beo phi : Binh thuong 0.000819 0.00246 0.000782 0.00235 0.000924
## 2 Beo phi : Gay 0.007680 0.01150 0.007980 0.01200 0.008100
## 3 Binh thuong : Gay 0.115000 0.11500 0.095800 0.09580 0.122000
## p.adj.Chisq
## 1 0.00277
## 2 0.01220
## 3 0.12200So sánh trung bình 2 nhóm độc lập
#Kiểm tra giả định phương sai đồng nhất
var.test(matdoxuong ~ beophi)##
## F test to compare two variances
##
## data: matdoxuong by beophi
## F = 1.0258, num df = 1113, denom df = 990, p-value = 0.6806
## alternative hypothesis: true ratio of variances is not equal to 1
## 95 percent confidence interval:
## 0.9085945 1.1577297
## sample estimates:
## ratio of variances
## 1.025848#Kiểm tra giả định phương sai đồng nhất
Desc(matdoxuong ~ beophi)## ------------------------------------------------------------------------------
## matdoxuong ~ beophi
##
## Summary:
## n pairs: 2'200, valid: 2'105 (95.7%), missings: 95 (4.3%), groups: 2
##
##
## Co Khong
## mean 0.875 0.778
## median 0.860 0.770
## sd 0.148 0.146
## IQR 0.200 0.190
## n 1'114 991
## np 52.922% 47.078%
## NAs 25 15
## 0s 0 0
##
## Kruskal-Wallis rank sum test:
## Kruskal-Wallis chi-squared = 193.76, df = 1, p-value < 0.00000000000000022
## Warning:
## Grouping variable contains 55 NAs (2.5%).
So sánh trung bình mật độ xương của nhóm nam và nữ
t.test(matdoxuong ~ gioitinh)##
## Welch Two Sample t-test
##
## data: matdoxuong by gioitinh
## t = 20.264, df = 1568.7, p-value < 0.00000000000000022
## alternative hypothesis: true difference in means is not equal to 0
## 95 percent confidence interval:
## 0.1185807 0.1439978
## sample estimates:
## mean in group Nam mean in group Nu
## 0.9093560 0.7780668So sánh trung vị 2 nhóm độc lập
wilcox.test(matdoxuong ~ gioitinh)##
## Wilcoxon rank sum test with continuity correction
##
## data: matdoxuong by gioitinh
## W = 791619, p-value < 0.00000000000000022
## alternative hypothesis: true location shift is not equal to 0Desc(matdoxuong ~ gioitinh)## ------------------------------------------------------------------------------
## matdoxuong ~ gioitinh
##
## Summary:
## n pairs: 2'200, valid: 2'111 (96.0%), missings: 89 (4.0%), groups: 2
##
##
## Nam Nu
## mean 0.909 0.778
## median 0.900 0.775
## sd 0.153 0.132
## IQR 0.185 0.170
## n 823 1'288
## np 38.986% 61.014%
## NAs 32 57
## 0s 0 0
##
## Kruskal-Wallis rank sum test:
## Kruskal-Wallis chi-squared = 366.99, df = 1, p-value < 0.00000000000000022
So sánh trung bình > 2 nhóm độc lập
#Kiểm tra giả định phương sai đồng nhất
bartlett.test(matdoxuong ~ bmi_pl)##
## Bartlett test of homogeneity of variances
##
## data: matdoxuong by bmi_pl
## Bartlett's K-squared = 5.0994, df = 2, p-value = 0.07811#Kiểm tra giả định phương sai đồng nhất
sosanh = aov(matdoxuong ~ bmi_pl)
plot(sosanh)
So sánh trung bình mật độ x ư ơng giữa các phân nhóm BMI
sosanh = aov(matdoxuong ~ bmi_pl)
summary(sosanh)## Df Sum Sq Mean Sq F value Pr(>F)
## bmi_pl 2 5.29 2.6449 123.5 <0.0000000000000002 ***
## Residuals 2102 45.00 0.0214
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 95 observations deleted due to missingnessPhân tích hậu kiểm (post hoc)
TukeyHSD(sosanh)## Tukey multiple comparisons of means
## 95% family-wise confidence level
##
## Fit: aov(formula = matdoxuong ~ bmi_pl)
##
## $bmi_pl
## diff lwr upr p adj
## Binh thuong-Beo phi -0.09164596 -0.1068051 -0.07648682 0.0000000
## Gay-Beo phi -0.19691417 -0.2508543 -0.14297399 0.0000000
## Gay-Binh thuong -0.10526820 -0.1593785 -0.05115793 0.0000159So sánh trung vị > 2 nhóm độc lập
Desc(matdoxuong ~ bmi_pl)## ------------------------------------------------------------------------------
## matdoxuong ~ bmi_pl
##
## Summary:
## n pairs: 2'200, valid: 2'105 (95.7%), missings: 95 (4.3%), groups: 3
##
##
## Beo phi Binh thuong Gay
## mean 0.875 0.783 0.678
## median 0.860 0.780 0.650
## sd 0.148 0.143 0.179
## IQR 0.200 0.190 0.250
## n 1'114 949 42
## np 52.922% 45.083% 1.995%
## NAs 25 15 0
## 0s 0 0 0
##
## Kruskal-Wallis rank sum test:
## Kruskal-Wallis chi-squared = 205.68, df = 2, p-value < 0.00000000000000022
## Warning:
## Grouping variable contains 55 NAs (2.5%).