BÀI TẬP 1

library(tidyverse)

## Warning: package 'purrr' was built under R version 4.3.1

## ── Attaching core tidyverse packages ──────────────────────── tidyverse 2.0.0 ──
## ✔ dplyr     1.1.2     ✔ readr     2.1.4
## ✔ forcats   1.0.0     ✔ stringr   1.5.0
## ✔ ggplot2   3.4.2     ✔ tibble    3.2.1
## ✔ lubridate 1.9.2     ✔ tidyr     1.3.0
## ✔ purrr     1.0.1     
## ── Conflicts ────────────────────────────────────────── tidyverse_conflicts() ──
## ✖ 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

library (readxl)

## Warning: package 'readxl' was built under R version 4.3.1

library (ggplot2)
library (DescTools)

## Warning: package 'DescTools' was built under R version 4.3.3

library(epitools)

## Warning: package 'epitools' was built under R version 4.3.1

library(caret)

## Warning: package 'caret' was built under R version 4.3.3

## Loading required package: lattice
## 
## Attaching package: 'caret'
## 
## The following objects are masked from 'package:DescTools':
## 
##     MAE, RMSE
## 
## The following object is masked from 'package:purrr':
## 
##     lift

library(DT)

## Warning: package 'DT' was built under R version 4.3.1

1 Dữ liệu nghiên cứu

Dữ liệu được thu thập bằng cách khảo sát dân số hiện tại tháng 3 năm 1988 của Cục điều tra dân số Hoa Kỳ. Tổng cộng có 28155 quan sát trên 7 biến.

data <- read_excel("C:/Users/Admin/Downloads/CPS1988.xlsx")
str(data)

## tibble [28,155 × 8] (S3: tbl_df/tbl/data.frame)
##  $ rownames  : num [1:28155] 1 2 3 4 5 6 7 8 9 10 ...
##  $ wage      : num [1:28155] 355 123 370 755 594 ...
##  $ education : num [1:28155] 7 12 9 11 12 16 8 12 12 14 ...
##  $ experience: num [1:28155] 45 1 9 46 36 22 51 34 0 18 ...
##  $ ethnicity : chr [1:28155] "cauc" "cauc" "cauc" "cauc" ...
##  $ smsa      : chr [1:28155] "yes" "yes" "yes" "yes" ...
##  $ region    : chr [1:28155] "northeast" "northeast" "northeast" "northeast" ...
##  $ parttime  : chr [1:28155] "no" "yes" "no" "no" ...

datatable(data)

## Warning in instance$preRenderHook(instance): It seems your data is too big for
## client-side DataTables. You may consider server-side processing:
## https://rstudio.github.io/DT/server.html

Bộ dữ liệu gồm có 7 biến, trong đó có 4 biến định tính và 3 biến định lượng

• Biến Wage: Tiền lương của người quan sát (tính bằng đô la mỗi tuần).

• Biến education: Số năm đi học.

• Biến experience: Số năm kinh nghiệm làm việc.

• Biến ethnicity: Chủng tộc, gồm 2 giá trị “cauc”và “afam”.

• Biến smsa: Cá nhân đó có cư trú trong Khu vực Thống kê Đô thị Tiêu chuẩn (SMSA) không?, gồm 2 giá trị: “Yes” và “No”.

• Biến region: Vùng đất, gồm 4 giá trị “northeast”, “midwest”, “south”, “west”.

• Biến parttime: Cá nhân có làm việc bán thời gian không?, gồm 2 giá trị “Yes” và “No”.

2 Thống kê mô tả

2.1 Thống kê mô tả cho biến định tính

d <- data[,c("smsa","parttime","ethnicity","region")]

Biến smsa - Cá nhân có nơi cư trú trong khu vực đô thị tiêu chuẩn không

table(d$smsa)

## 
##    no   yes 
##  7223 20932

prop.table(table(d$smsa))

## 
##        no       yes 
## 0.2565441 0.7434559

ggplot(data,aes(smsa))+
  geom_bar(color = "blue", fill = "green")+
   geom_text(aes(label = scales :: percent(after_stat(count/sum(count)))), stat=  'count', color = 'black', vjust = -.5)+
  ylab("count")+ xlab("smsa") + labs(title = 'Đồ thị về nơi cư trú tiêu chuẩn của 28155 người được khảo sát')

Trong tổng số 28155 người tham gia khảo sát thì:

• Có 7223 số người không có cư trú trong khu vực thống kê đô thị đạt tiêu chuẩn, chiếm 26%.

• Có 20932 số người có cư trú trong khu vực thống kê đô thị đạt tiêu chuẩn, chiếm 74%.

Biến parttime - Cá nhân có làm việc bán thời gian không

table(d$parttime)

## 
##    no   yes 
## 25631  2524

prop.table(table(d$parttime))

## 
##        no       yes 
## 0.9103534 0.0896466

ggplot(data,aes(parttime))+
  geom_bar(color = "blue", fill = "green")+
   geom_text(aes(label = scales :: percent(after_stat(count/sum(count)))), stat=  'count', color = 'black', vjust = -.5)+
  ylab("count")+ xlab("parttime") + labs(title = 'Đồ thị về việc làm bán thời gian của 28155 người được khảo sát')

Trong tổng số 28155 người tham gia khảo sát thì có 25631 số người không có công việc bán thời gian, chiếm 91%. Có 2524 số người có công việc bán thời gian, chiếm 9%.

Biến ethnicity - Chủng tộc

table(d$ethnicity)

## 
##  afam  cauc 
##  2232 25923

prop.table(table(d$ethnicity))

## 
##       afam       cauc 
## 0.07927544 0.92072456

ggplot(data,aes(ethnicity))+
  geom_bar(color = "blue", fill = "green")+
   geom_text(aes(label = scales :: percent(after_stat(count/sum(count)))), stat=  'count', color = 'black', vjust = -.5)+
  ylab("count")+ xlab("ethnicity") + labs(title = 'Đồ thị về chủng tộc của 28155 người được khảo sát')

Trong tổng số 28155 người tham gia khảo sát thì có 2232 số người thuộc chủng tộc afam, chiếm 8%. Có 25923 số người thuộc chủng tộc cauc, chiếm 92%.

Biến region - Vùng đất

table(d$region)

## 
##   midwest northeast     south      west 
##      6863      6441      8760      6091

prop.table(table(d$region))

## 
##   midwest northeast     south      west 
## 0.2437578 0.2287693 0.3111348 0.2163381

ggplot(data,aes(region))+
  geom_bar(color = "blue", fill = "green")+
   geom_text(aes(label = scales :: percent(after_stat(count/sum(count)))), stat=  'count', color = 'black', vjust = -.5)+
  ylab("count")+ xlab("region") + labs(title = 'Đồ thị về vùng đất của 28155 người được khảo sát')

Trong tổng số 28155 người tham gia khảo sát thì:

• Có 6863 số người sống ở phía Trung Tây, chiếm 24.4%,

• Có 6441 số người sống ở phía Đông Bắc, chiếm 22.9%,

• Có 8760 số người sống ở phía Nam, chiếm 31.1%,

• Có 6091 số người sống ở phía Tây, chiếm 21.6%.

2.2 Thống kê mô tả cho biến định lượng

d <- data[,c("wage","education","experience")]

Biến Wage - Tiền lương

summary(data$wage)

##     Min.  1st Qu.   Median     Mean  3rd Qu.     Max. 
##    50.05   308.64   522.32   603.73   783.48 18777.20

Trong 28155 người được khảo sát cho biết tiền lương:

• Thấp nhất là 50.05

• Cao nhất là 18777.20

• Trung bình là 603.73

• 50% người trong nhóm khỏa sát có tiền lương thấp hơn 522.32

• 25% người trong nhóm khảo sát có tiền lương thấp hơn 308.64

• 75% người trong nhóm khảo sát có tiền lương thấp hơn 783.48

sd(data$wage)

## [1] 453.5474

Độ lệch chuẩn của tiền lương là 453.5474, cho biết mức độ phân tán của tiền lương của từng người được khảo sát so với giá trị trung bình là 453.5474.

wage1 <- cut(data$wage, breaks =  c(0,200,400,600,800,1000,1200,1400,1600,18777), labels = c('0 - 200', '201 - 400','401 - 600','601 - 800','801 - 1000','1001 - 1200','1201 - 1400','1401 - 1600','>1600'))

d <- mutate(data, wage1)

d %>% na.omit() %>% ggplot(aes(x = wage1, y = after_stat(count))) + geom_bar(fill = 'skyblue') + geom_text(aes(label= scales::percent(after_stat(count/sum(count)), accuracy= 0.01)), stat = 'count', color = 'darkblue', vjust = -0.5)

Từ đồ thị trên ta thấy, những người được khảo sát có tiền lương chủ yếu từ 201 đến 600, trong đó nhóm từ 1401 đến 1600, chiếm tỷ lệ thấp nhất 1.66%.

Biến education - Trình độ học vấn

summary(data$education)

##    Min. 1st Qu.  Median    Mean 3rd Qu.    Max. 
##    0.00   12.00   12.00   13.07   15.00   18.00

Trong 28155 người được khảo sát cho biết trình độ học vấn:

• Thấp nhất là 0

• Cao nhất là 18

• Trung bình là 13.07

• 50% người trong nhóm khỏa sát có trình độ học vấn thấp hơn 12

• 25% người trong nhóm khảo sát có trình độ học vấn hơn 12

• 75% người trong nhóm khảo sát có trình độ học vấn hơn 15

sd(data$education)

## [1] 2.899682

Độ lệch chuẩn của trình độ học vấn là 2.899682, cho biết mức độ phân tán của trình độ học vấn của từng người được khảo sát so với giá trị trung bình là 2.899682.

table(data$education)

## 
##     0     1     2     3     4     5     6     7     8     9    10    11    12 
##    79    22    65   102    90   104   330   204   663   650  1022  1083 10549 
##    13    14    15    16    17    18 
##  2052  2967  1154  3873   815  2331

prop.table(table(data$education))

## 
##            0            1            2            3            4            5 
## 0.0028058959 0.0007813887 0.0023086486 0.0036228023 0.0031965903 0.0036938377 
##            6            7            8            9           10           11 
## 0.0117208311 0.0072456047 0.0235482152 0.0230864855 0.0362990588 0.0384656367 
##           12           13           14           15           16           17 
## 0.3746759013 0.0728822589 0.1053809270 0.0409873912 0.1375599361 0.0289469011 
##           18 
## 0.0827916889

education1 <- cut(data$education, breaks =  c(0,5,9,12,17), labels = c('0 - 5', '6 - 9','10 - 12','>12'))

d <- mutate(data, education1)

d %>% na.omit() %>% ggplot(aes(x = education1, y = after_stat(count))) + geom_bar(fill = 'skyblue') + geom_text(aes(label= scales::percent(after_stat(count/sum(count)), accuracy= 0.01)), stat = 'count', color = 'darkblue', vjust = -0.5)

Từ đồ thị trên ta thấy, những người được khảo sát có trình độ học vấn chủ yếu từ 10 đến 12 năm chiếm 49.15%, số người có trình độ học vấn từ 0 đến 5 năm chiếm tỷ lệ thấp nhất 1.49%

Biến experience - Kinh nghiệm

summary(data$experience)

##    Min. 1st Qu.  Median    Mean 3rd Qu.    Max. 
##    -4.0     8.0    16.0    18.2    27.0    63.0

Trong 28155 người được khảo sát cho biết số năm kinh nghiệm:

• Thấp nhất là -4

• Cao nhất là 63

• Trung bình là 18.2

• 50% người trong nhóm khỏa sát có số năm kinh nghiệm thấp hơn 16

• 25% người trong nhóm khảo sát có số năm kinh nghiệm hơn 8

• 75% người trong nhóm khảo sát có số năm kinh nghiệm hơn 27

sd(data$experience)

## [1] 13.07923

Độ lệch chuẩn của số năm kinh nghiệm là 13.07923, cho biết mức độ phân tán của số năm kinh nghiệm của từng người được khảo sát so với giá trị trung bình là 13.07923.

experience1 <- cut(data$experience, breaks =  c(-4,10,25,45,62), labels = c('-4 - 10', '11 - 25','26 - 45','46 - 63'))

d <- mutate(data, experience1)

d %>% na.omit() %>% ggplot(aes(x = experience1, y = after_stat(count))) + geom_bar(fill = 'skyblue') + geom_text(aes(label= scales::percent(after_stat(count/sum(count)), accuracy= 0.01)), stat = 'count', color = 'darkblue', vjust = -0.5)

Từ đồ thị trên ta thấy, những người được khảo sát có số năm kinh nghiệm chủ yếu từ 11 đến 25 năm chiếm 39.24%, số người có trình độ học vấn từ 46 đến 63 năm chiếm tỷ lệ thấp nhất 3.25%

3 Relative Risk/Risk Ratio

3.1 Rủi ro tương đối giữa biến parttime (có công việc bán thời gian) và biến smsa (khu vực cư trú tiêu chuẩn)

risk1 <- table(data$smsa,data$parttime) 
addmargins(risk1)

##      
##          no   yes   Sum
##   no   6591   632  7223
##   yes 19040  1892 20932
##   Sum 25631  2524 28155

RelRisk(risk1)

## [1] 1.003177

Tỷ lệ người có công việc bán thời gian không có nơi cư trú tiêu chuẩn bằng 100,32% tỷ lệ người có công việc bán thời gian nhưng có nơi cư trú tiêu chuẩn.

riskratio(risk1)

## $data
##        
##            no  yes Total
##   no     6591  632  7223
##   yes   19040 1892 20932
##   Total 25631 2524 28155
## 
## $measure
##      risk ratio with 95% C.I.
##       estimate     lower   upper
##   no  1.000000        NA      NA
##   yes 1.033025 0.9479131 1.12578
## 
## $p.value
##      two-sided
##       midp.exact fisher.exact chi.square
##   no          NA           NA         NA
##   yes  0.4595492    0.4736688  0.4585465
## 
## $correction
## [1] FALSE
## 
## attr(,"method")
## [1] "Unconditional MLE & normal approximation (Wald) CI"

Tỷ lệ người có công việc bán thời gian có nơi cư trú tiêu chuẩn bằng 1.033 lần tỷ lệ người có công việc bán thời gian nhưng không có nơi cư trú tiêu chuẩn.

3.2 Rủi ro tương đối giữa biến parttime (có công việc bán thời gian) và biến ethnicity (chủng tộc)

risk2 <- table(data$ethnicity,data$parttime) 
addmargins(risk2)

##       
##           no   yes   Sum
##   afam  1988   244  2232
##   cauc 23643  2280 25923
##   Sum  25631  2524 28155

RelRisk(risk2)

## [1] 0.9765733

Tỷ lệ số người có công việc bán thời gian thuộc chủng tộc người Mỹ gốc Phi bằng 97.66% tỷ lệ số người có công việc bán thời gian thuộc chủng tộc người gốc Âu.

riskratio(risk2)

## $data
##        
##            no  yes Total
##   afam   1988  244  2232
##   cauc  23643 2280 25923
##   Total 25631 2524 28155
## 
## $measure
##       risk ratio with 95% C.I.
##         estimate     lower     upper
##   afam 1.0000000        NA        NA
##   cauc 0.8045517 0.7102013 0.9114366
## 
## $p.value
##       two-sided
##          midp.exact fisher.exact   chi.square
##   afam           NA           NA           NA
##   cauc 0.0009550463   0.00101711 0.0006975761
## 
## $correction
## [1] FALSE
## 
## attr(,"method")
## [1] "Unconditional MLE & normal approximation (Wald) CI"

Tỷ lệ số người có công việc bán thời gian thuộc chủng tộc người gốc Âu bằng 0.8 lần tỷ lệ số người có công việc bán thời gian thuộc chủng tộc người Mỹ gốc Phi.

4 Tỷ lệ chênh (Odd Ratio)

4.1 Tỷ lệ chênh giữa biến parttime (có công việc bán thời gian) và biến smsa (khu vực cư trú tiêu chuẩn)

odd1 <- table(data$smsa,data$parttime)
odd1

##      
##          no   yes
##   no   6591   632
##   yes 19040  1892

oddsratio(odd1)

## $data
##        
##            no  yes Total
##   no     6591  632  7223
##   yes   19040 1892 20932
##   Total 25631 2524 28155
## 
## $measure
##      odds ratio with 95% C.I.
##       estimate     lower    upper
##   no  1.000000        NA       NA
##   yes 1.036141 0.9434217 1.139305
## 
## $p.value
##      two-sided
##       midp.exact fisher.exact chi.square
##   no          NA           NA         NA
##   yes  0.4595492    0.4736688  0.4585465
## 
## $correction
## [1] FALSE
## 
## attr(,"method")
## [1] "median-unbiased estimate & mid-p exact CI"

Tỷ lệ người có công việc bán thời gian/người không có công việc bán thời gian mà có nơi cư trú tiêu chuẩn bằng 1.36 lần so với Tỷ lệ người có công việc bán thời gian/người không có công việc bán thời gian mà không có nơi cư trú tiêu chuẩn

4.2 Tỷ lệ chênh giữa biến parttime (có công việc bán thời gian) và biến ethnicity (chủng tộc)

odd2 <- table(data$ethnicity,data$parttime)
odd2

##       
##           no   yes
##   afam  1988   244
##   cauc 23643  2280

oddsratio(odd2)

## $data
##        
##            no  yes Total
##   afam   1988  244  2232
##   cauc  23643 2280 25923
##   Total 25631 2524 28155
## 
## $measure
##       odds ratio with 95% C.I.
##         estimate     lower     upper
##   afam 1.0000000        NA        NA
##   cauc 0.7852903 0.6842119 0.9049881
## 
## $p.value
##       two-sided
##          midp.exact fisher.exact   chi.square
##   afam           NA           NA           NA
##   cauc 0.0009550463   0.00101711 0.0006975761
## 
## $correction
## [1] FALSE
## 
## attr(,"method")
## [1] "median-unbiased estimate & mid-p exact CI"

Tỷ lệ người có công việc bán thời gian/người không có công việc bán thời gian thuộc người gốc Âu bằng 0.79 lần so với Tỷ lệ người có công việc bán thời gian/người không có công việc bán thời gian thuộc người Mỹ gốc Phi

5 Thống kê suy diễn cho biến định tính

5.1 Kiểm tra tính độc lập của 2 biến khu vực cư trú tiêu chuẩn và công việc bán thời gian

Giả thuyết:

• H0: khu vực cư trú tiêu chuẩn và công việc bán thời gian độc lập

• H1: khu vực cư trú tiêu chuẩn và công việc bán thời gian không độc lập nhau

chisq.test(odd1)

## 
##  Pearson's Chi-squared test with Yates' continuity correction
## 
## data:  odd1
## X-squared = 0.5146, df = 1, p-value = 0.4732

Qua kết quả kiểm định ta thấy giá trị p_value = 0.4732 > 0.05 (mức ý nghĩa thông thường), nên chấp nhận giả thuyết H0. Vậy với mức ý nghĩa 5%, khu vực cư trú tiêu chuẩn và công việc bán thời gian độc lập nhau.

5.2 Kiểm tra tính độc lập của 2 biến chủng tộc và công việc bán thời gian

Giả thuyết:

• H0: chủng tộc và công việc bán thời gian gian độc lập

• H1: chủng tộc và công việc bán thời gian không độc lập nhau

chisq.test(odd2)

## 
##  Pearson's Chi-squared test with Yates' continuity correction
## 
## data:  odd2
## X-squared = 11.235, df = 1, p-value = 0.0008025

Qua kết quả kiểm định ta thấy giá trị p_value = 0.0008025 < 0.05 (mức ý nghĩa thông thường), nên bác bỏ giả thuyết H0. Vậy với mức ý nghĩa 5%, chủng tộc và công việc bán thời gian gian không độc lập.

6 Khoảng ước lượng cho tỷ lệ

6.1 Công thức ước lượng:

Để làm bài toán ước lượng tỷ lệ bằng R, ta thực hiện câu lệnh sau:

B1. install.packages(“DescTools”): cài đặt gói DescTools.

B2. library(DecsTools): gọi gói DescTools.

B3. BinomCI(k, n, conf.level = NULL, method = “wald”): câu lệnh ước lượng tỷ lệ một tổng thể; trong đó, k: số phần tử có tính chất A, n: số phần tử mẫu, conf.level = NULL: độ tin cậy mặc định là 95%, method = wald: phương pháp wald.

6.2 Để ước lượng tỷ lệ những người có số năm đi học dưới 12 năm ta làm như sau:

m <- data$education
m2 <- m[m < 12]
k <- length(m2[m2 < 3])
n <- length(m2)
BinomCI(k, n, conf.level = 0.95, method = "wald")

##             est     lwr.ci     upr.ci
## [1,] 0.03760761 0.03199525 0.04321998

Ta được khoảng ước lượng tỷ lệ những người có số năm đi học dưới 12 năm với độ tin cậy 95% là: 0.03199525 <= p <= 0.04321998

6.3 Ước lượng tỷ lệ người sống ở khu vực cư trú tiêu chuẩn có công việc bán thời gian đồng thời kiểm định xem tỷ lệ (%) người sống ở khu vực cư trú tiêu chuẩn có công việc bán thời gian có bằng 40% không?

#Tập hợp dữ liệu chỉ bao gồm những người có công việc bán thời gian
rm<- data[data$parttime == "yes",]
#Tập hợp con những người có công việc bán thời gian sống ở khu vực cư trú tiêu chuẩn
rm1<-rm[rm$smsa == "yes",]
#kiểm tra tỷ lệ
prop.test(length(rm1$smsa), length(rm$smsa), p = 0.4)

## 
##  1-sample proportions test with continuity correction
## 
## data:  length(rm1$smsa) out of length(rm$smsa), null probability 0.4
## X-squared = 1283.9, df = 1, p-value < 2.2e-16
## alternative hypothesis: true p is not equal to 0.4
## 95 percent confidence interval:
##  0.7321290 0.7663111
## sample estimates:
##         p 
## 0.7496038

Với khoảng tin cậy 95% ước lượng tỷ lệ người sống ở khu vực cư trú tiêu chuẩn có công việc bán thời gian nằm trong khoảng từ 0.7321290 đến 0.7663111

p-value = 2.2e-16 < 0.05,nên bác bỏ giả thuyết H0. Do đó tỷ lệ (%) người sống ở khu vực cư trú tiêu chuẩn có công việc bán thời gian không bằng 40% với mức ý nghĩa 5%

6.4 Ước lượng tỷ lệ người sống ở khu vực cư trú tiêu chuẩn không công việc bán thời gian đồng thời kiểm định xem tỷ lệ (%) người sống ở khu vực cư trú tiêu chuẩn không có công việc bán thời gian có bằng 40% không?

rm<- data[data$parttime == "no",]
rm1<-rm[rm$smsa == "yes",]
prop.test(length(rm1$smsa), length(rm$smsa), p = 0.4)

## 
##  1-sample proportions test with continuity correction
## 
## data:  length(rm1$smsa) out of length(rm$smsa), null probability 0.4
## X-squared = 12552, df = 1, p-value < 2.2e-16
## alternative hypothesis: true p is not equal to 0.4
## 95 percent confidence interval:
##  0.7374440 0.7481838
## sample estimates:
##         p 
## 0.7428505

Với khoảng tin cậy 95% ước lượng tỷ lệ người sống ở khu vực cư trú tiêu chuẩn không có công việc bán thời gian nằm trong khoảng từ 0.7374440 đến 0.7481838

p-value = 2.2e-16 < 0.05,nên bác bỏ giả thuyết H0. Do đó tỷ lệ (%) người sống ở khu vực cư trú tiêu chuẩn không có công việc bán thời gian không bằng 40% với mức ý nghĩa 5%.