Việc 1: Đánh giá tiến bộ của học sinh
wil = read.csv("C:\\Thach\\VN trips\\VN trip 5 (Apr 2023)\\Datasets\\willett.csv", header = T)
1.1 Vẽ biểu đồ bánh tằm đánh giá tiến bộ của học sinh
library(ggplot2)
ggplot(data = wil, aes(x = time, y = y, group = id, col = id)) + geom_line() + theme(legend.position = "none")

1.2 Thực hiện mô hình tuyến tính
m.12 = lm(y ~ time, data = wil)
summary(m.12)
##
## Call:
## lm(formula = y ~ time, data = wil)
##
## Residuals:
## Min 1Q Median 3Q Max
## -88.374 -25.584 1.186 28.926 64.746
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 164.374 5.035 32.65 <2e-16 ***
## time 26.960 2.691 10.02 <2e-16 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 35.6 on 138 degrees of freedom
## Multiple R-squared: 0.421, Adjusted R-squared: 0.4168
## F-statistic: 100.3 on 1 and 138 DF, p-value: < 2.2e-16
1.3 Thưc hiện mô hình tuyến tính cho từng học sinh
library(tidyverse)
## ── Attaching core tidyverse packages ──────────────────────── tidyverse 2.0.0 ──
## ✔ dplyr 1.1.1 ✔ readr 2.1.4
## ✔ forcats 1.0.0 ✔ stringr 1.5.0
## ✔ lubridate 1.9.2 ✔ tibble 3.2.1
## ✔ purrr 1.0.1 ✔ tidyr 1.3.0
## ── 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
wil %>% group_by(id) %>% mutate(beta = lm(y ~ time)$coefficients[2]) %>% mutate(alpha = lm(y ~ time)$coefficients[1]) %>% summarize(mbeta=mean(beta), malpha=mean(alpha))
## # A tibble: 35 × 3
## id mbeta malpha
## <int> <dbl> <dbl>
## 1 1 34.2 197.
## 2 2 28.5 218
## 3 3 47.9 151.
## 4 4 21.9 206.
## 5 5 9.9 200.
## 6 6 29.6 166.
## 7 7 32.6 160.
## 8 8 37.9 92.9
## 9 9 22.7 139.
## 10 10 26.8 220.
## # ℹ 25 more rows
1.4 Thực hiện mô hình ảnh hưởng hỗn hợp
library(lme4)
## Loading required package: Matrix
##
## Attaching package: 'Matrix'
## The following objects are masked from 'package:tidyr':
##
## expand, pack, unpack
m.14 = lmer(y ~ time + (1 + time | id), data = wil)
summary(m.14)
## Linear mixed model fit by REML ['lmerMod']
## Formula: y ~ time + (1 + time | id)
## Data: wil
##
## REML criterion at convergence: 1266.8
##
## Scaled residuals:
## Min 1Q Median 3Q Max
## -2.27205 -0.60035 0.02632 0.60730 1.58439
##
## Random effects:
## Groups Name Variance Std.Dev. Corr
## id (Intercept) 1198.8 34.62
## time 132.4 11.51 -0.45
## Residual 159.5 12.63
## Number of obs: 140, groups: id, 35
##
## Fixed effects:
## Estimate Std. Error t value
## (Intercept) 164.374 6.119 26.86
## time 26.960 2.167 12.44
##
## Correlation of Fixed Effects:
## (Intr)
## time -0.489
Việc 2: Điểm năm học trước có ảnh hưởng như thế nào đến tiến bộ của
học sinh?
2.1 Thực hiện mô hình tuyến tính
wil$ccovar = wil$covar - mean(wil$covar)
m.21 = lm(y ~ time + ccovar, data = wil)
summary(m.21)
##
## Call:
## lm(formula = y ~ time + ccovar, data = wil)
##
## Residuals:
## Min 1Q Median 3Q Max
## -86.522 -24.591 3.165 28.033 70.018
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 164.3743 4.9646 33.109 <2e-16 ***
## time 26.9600 2.6537 10.159 <2e-16 ***
## ccovar 0.5357 0.2410 2.223 0.0278 *
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 35.11 on 137 degrees of freedom
## Multiple R-squared: 0.4412, Adjusted R-squared: 0.433
## F-statistic: 54.08 on 2 and 137 DF, p-value: < 2.2e-16
2.2 Thực hiện mô hình tuyến tính có yếu tố tương tác
m.22 = lm(y ~ time + ccovar + time:ccovar, data = wil)
summary(m.22)
##
## Call:
## lm(formula = y ~ time + ccovar + time:ccovar, data = wil)
##
## Residuals:
## Min 1Q Median 3Q Max
## -88.767 -26.040 3.676 26.224 73.433
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 164.3743 4.9089 33.485 <2e-16 ***
## time 26.9600 2.6239 10.275 <2e-16 ***
## ccovar -0.1136 0.3987 -0.285 0.7762
## time:ccovar 0.4329 0.2131 2.031 0.0442 *
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 34.71 on 136 degrees of freedom
## Multiple R-squared: 0.4576, Adjusted R-squared: 0.4457
## F-statistic: 38.25 on 3 and 136 DF, p-value: < 2.2e-16
2.3 Thực hiện mô hình ảnh hưởng hỗn hợp
m.23 = lmer(y ~ time + ccovar + time:ccovar + ( 1 + time | id), data = wil)
summary(m.23)
## Linear mixed model fit by REML ['lmerMod']
## Formula: y ~ time + ccovar + time:ccovar + (1 + time | id)
## Data: wil
##
## REML criterion at convergence: 1260.3
##
## Scaled residuals:
## Min 1Q Median 3Q Max
## -2.24817 -0.61873 0.00428 0.61472 1.55688
##
## Random effects:
## Groups Name Variance Std.Dev. Corr
## id (Intercept) 1236.4 35.16
## time 107.2 10.36 -0.49
## Residual 159.5 12.63
## Number of obs: 140, groups: id, 35
##
## Fixed effects:
## Estimate Std. Error t value
## (Intercept) 164.3743 6.2061 26.486
## time 26.9600 1.9939 13.521
## ccovar -0.1136 0.5040 -0.225
## time:ccovar 0.4329 0.1619 2.673
##
## Correlation of Fixed Effects:
## (Intr) time ccovar
## time -0.522
## ccovar 0.000 0.000
## time:ccovar 0.000 0.000 -0.522
Việc 3: Yếu tố ảnh hưởng đến điểm số học sinh
pop = read.csv("C:\\Thach\\VN trips\\VN trip 5 (Apr 2023)\\Datasets\\popular.csv")
3.1 Mô tả bằng biểu đồ
a. Phân bố điểm của học sinh
ggplot(data = pop, aes(popular)) + geom_histogram(fill = "blue", col = "white")
## `stat_bin()` using `bins = 30`. Pick better value with `binwidth`.

b. Điểm học sinh nam có khác học sinh nữ không?
ggplot(data = pop, aes(x = sex, y = popular, col = sex)) + geom_boxplot() + geom_jitter(alpha = 0.1)

c.1 Ngoại cảnh ảnh hưởng thế nào đến điểm học sinh?
ggplot(data = pop, aes(x = factor(extrov), y = popular, col = factor(extrov))) + geom_boxplot() + geom_jitter(alpha = 0.1)

c.2 Ảnh hưởng của ngoại cảnh có khác biệt giữa nam và nữ không?
ggplot(data = pop, aes(x=factor(extrov), y = popular, col = sex)) + geom_boxplot() + geom_jitter(alpha = 0.1)

3.2 Dữ liệu này có gì đặc biệt?
3.3 Điểm của học sinh có khác nhau giữa các lớp không?
m.33 = lmer(popular ~ 1 + (1 | class), data = pop)
summary(m.33)
## Linear mixed model fit by REML ['lmerMod']
## Formula: popular ~ 1 + (1 | class)
## Data: pop
##
## REML criterion at convergence: 6330.5
##
## Scaled residuals:
## Min 1Q Median 3Q Max
## -3.5655 -0.6975 0.0020 0.6758 3.3175
##
## Random effects:
## Groups Name Variance Std.Dev.
## class (Intercept) 0.7021 0.8379
## Residual 1.2218 1.1053
## Number of obs: 2000, groups: class, 100
##
## Fixed effects:
## Estimate Std. Error t value
## (Intercept) 5.07786 0.08739 58.1
3.4 Điểm học sinh nam có khác nữ không?
a. Không khác biệt giữa các lớp
m.34a = lmer(popular ~ 1 + Csex + (1 | class), data = pop)
summary(m.34a)
## Linear mixed model fit by REML ['lmerMod']
## Formula: popular ~ 1 + Csex + (1 | class)
## Data: pop
##
## REML criterion at convergence: 5564.1
##
## Scaled residuals:
## Min 1Q Median 3Q Max
## -3.6890 -0.6415 0.0040 0.6788 3.1586
##
## Random effects:
## Groups Name Variance Std.Dev.
## class (Intercept) 0.4839 0.6957
## Residual 0.8308 0.9115
## Number of obs: 2000, groups: class, 100
##
## Fixed effects:
## Estimate Std. Error t value
## (Intercept) 5.06951 0.07252 69.91
## Csex 1.35010 0.04403 30.66
##
## Correlation of Fixed Effects:
## (Intr)
## Csex -0.004
b. Có khác biệt giữa các lớp
m.34b = lmer(popular ~ 1 + Csex + (1 + Csex | class), data = pop)
summary(m.34b)
## Linear mixed model fit by REML ['lmerMod']
## Formula: popular ~ 1 + Csex + (1 + Csex | class)
## Data: pop
##
## REML criterion at convergence: 5559.1
##
## Scaled residuals:
## Min 1Q Median 3Q Max
## -3.5848 -0.6405 0.0062 0.6664 3.1420
##
## Random effects:
## Groups Name Variance Std.Dev. Corr
## class (Intercept) 0.48299 0.6950
## Csex 0.05461 0.2337 -0.30
## Residual 0.81896 0.9050
## Number of obs: 2000, groups: class, 100
##
## Fixed effects:
## Estimate Std. Error t value
## (Intercept) 5.07291 0.07254 69.93
## Csex 1.35218 0.05030 26.88
##
## Correlation of Fixed Effects:
## (Intr)
## Csex -0.144
3.5 Ngoại cảnh ảnh hưởng như thế nào đến điểm số?
a. Không có khác biệt giữa các lớp
m.35a = lmer(popular ~ 1 + Cextrav + (1 | class), data = pop)
summary(m.35a)
## Linear mixed model fit by REML ['lmerMod']
## Formula: popular ~ 1 + Cextrav + (1 | class)
## Data: pop
##
## REML criterion at convergence: 5832.6
##
## Scaled residuals:
## Min 1Q Median 3Q Max
## -3.0644 -0.7267 0.0165 0.7088 3.3587
##
## Random effects:
## Groups Name Variance Std.Dev.
## class (Intercept) 0.8406 0.9168
## Residual 0.9304 0.9646
## Number of obs: 2000, groups: class, 100
##
## Fixed effects:
## Estimate Std. Error t value
## (Intercept) 5.07824 0.09421 53.90
## Cextrav 0.48631 0.02015 24.13
##
## Correlation of Fixed Effects:
## (Intr)
## Cextrav 0.000
b. Có khác biệt giữa các lớp
m.35b = lmer(popular ~ 1 + Cextrav + (1 + Cextrav | class), data = pop)
summary(m.35b)
## Linear mixed model fit by REML ['lmerMod']
## Formula: popular ~ 1 + Cextrav + (1 + Cextrav | class)
## Data: pop
##
## REML criterion at convergence: 5779.4
##
## Scaled residuals:
## Min 1Q Median 3Q Max
## -3.1961 -0.7291 0.0146 0.6816 3.2217
##
## Random effects:
## Groups Name Variance Std.Dev. Corr
## class (Intercept) 0.89178 0.9443
## Cextrav 0.02599 0.1612 -0.88
## Residual 0.89492 0.9460
## Number of obs: 2000, groups: class, 100
##
## Fixed effects:
## Estimate Std. Error t value
## (Intercept) 5.03127 0.09702 51.86
## Cextrav 0.49286 0.02546 19.36
##
## Correlation of Fixed Effects:
## (Intr)
## Cextrav -0.552
3.6 Ngoại cảnh và giới tính của học sinh có ảnh hưởng như thế nào
đến điểm số
m.36 = lmer(popular ~ 1 + Cextrav + Csex + (1 + Cextrav + Csex | class), data = pop)
## boundary (singular) fit: see help('isSingular')
summary(m.36)
## Linear mixed model fit by REML ['lmerMod']
## Formula: popular ~ 1 + Cextrav + Csex + (1 + Cextrav + Csex | class)
## Data: pop
##
## REML criterion at convergence: 4870.5
##
## Scaled residuals:
## Min 1Q Median 3Q Max
## -3.01902 -0.64955 -0.01055 0.67099 3.11761
##
## Random effects:
## Groups Name Variance Std.Dev. Corr
## class (Intercept) 0.673941 0.82094
## Cextrav 0.029847 0.17276 -0.73
## Csex 0.005365 0.07325 -0.65 -0.03
## Residual 0.552883 0.74356
## Number of obs: 2000, groups: class, 100
##
## Fixed effects:
## Estimate Std. Error t value
## (Intercept) 5.02051 0.08411 59.69
## Cextrav 0.44300 0.02343 18.91
## Csex 1.24483 0.03728 33.39
##
## Correlation of Fixed Effects:
## (Intr) Cextrv
## Cextrav -0.533
## Csex -0.127 -0.065
## optimizer (nloptwrap) convergence code: 0 (OK)
## boundary (singular) fit: see help('isSingular')