dta <- read.csv('Chemistry100.csv',header=T)
head(dta)
## schID stuID gender sciscore
## 1 1 1 M 6.875
## 2 1 2 M 5.875
## 3 1 3 M 8.000
## 4 1 4 F 6.625
## 5 1 5 M 6.125
## 6 1 6 M 7.250
dta$schID <- as.factor(dta$schID)
library(pacman)
## Warning: 套件 'pacman' 是用 R 版本 4.2.3 來建造的
p_load(lme4,lmerTest,ggplot2)
rstsep <- lmList(sciscore ~ gender | schID, data=dta)
IntSlope <- as.data.frame(coef(rstsep))
names(IntSlope) <- c('Intercept','genderSlope')
head(IntSlope)
## Intercept genderSlope
## 1 6.180778 0.43907937
## 2 6.656429 -0.63672269
## 3 5.747071 0.05865584
## 4 7.500000 -0.30084000
## 5 6.309667 -0.65900000
## 6 6.739806 -0.37417398
注意一下第一列
rstfix <- lm(sciscore~schID+gender+schID:gender,data=dta)
summary(rstfix)
##
## Call:
## lm(formula = sciscore ~ schID + gender + schID:gender, data = dta)
##
## Residuals:
## Min 1Q Median 3Q Max
## -3.13630 -0.49724 0.03323 0.52237 2.09105
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 6.180778 0.256983 24.051 < 2e-16 ***
## schID2 0.475651 0.388521 1.224 0.220937
## schID3 -0.433706 0.329385 -1.317 0.188020
## schID4 1.319222 0.339956 3.881 0.000106 ***
## schID5 0.128889 0.406325 0.317 0.751106
## schID6 0.559028 0.287315 1.946 0.051771 .
## schID7 0.619008 0.329385 1.879 0.060289 .
## schID8 0.415347 0.374613 1.109 0.267620
## schID9 -0.087248 0.317809 -0.275 0.783693
## schID10 -0.563340 0.321228 -1.754 0.079568 .
## schID11 -0.673778 0.299691 -2.248 0.024623 *
## schID12 -0.817965 0.321228 -2.546 0.010928 *
## schID13 -0.364298 0.299691 -1.216 0.224227
## schID14 -0.028000 0.363428 -0.077 0.938593
## schID15 1.009222 0.354226 2.849 0.004410 **
## schID16 -0.025678 0.354226 -0.072 0.942216
## schID17 0.770949 0.346515 2.225 0.026154 *
## schID18 -0.185473 0.303120 -0.612 0.540656
## schID19 0.058012 0.311965 0.186 0.852490
## schID20 0.381815 0.296738 1.287 0.198281
## schID21 0.430918 0.303120 1.422 0.155229
## schID22 0.576710 0.283790 2.032 0.042211 *
## schID23 0.050197 0.284428 0.176 0.859923
## schID24 -0.025009 0.334305 -0.075 0.940372
## schID25 0.066139 0.339956 0.195 0.845756
## schID26 0.066739 0.269663 0.247 0.804543
## schID27 0.242299 0.334305 0.725 0.468633
## schID28 -0.119421 0.329385 -0.363 0.716958
## schID29 -0.143013 0.317809 -0.450 0.652740
## schID30 0.110131 0.346515 0.318 0.750636
## schID31 0.205915 0.334305 0.616 0.537969
## schID32 0.214589 0.293005 0.732 0.463990
## schID33 0.728131 0.346515 2.101 0.035686 *
## schID34 0.139431 0.301339 0.463 0.643606
## schID35 0.348068 0.334305 1.041 0.297869
## schID36 0.130125 0.283790 0.459 0.646604
## schID37 0.452089 0.325060 1.391 0.164379
## schID38 0.384847 0.374613 1.027 0.304341
## schID39 0.282495 0.346515 0.815 0.414987
## schID40 -0.088153 0.374613 -0.235 0.813977
## schID41 0.007722 0.354226 0.022 0.982609
## schID42 0.355131 0.346515 1.025 0.305498
## schID43 0.025822 0.354226 0.073 0.941892
## schID44 -0.013349 0.329385 -0.041 0.967675
## schID45 0.242222 0.374613 0.647 0.517938
## schID46 0.337556 0.325060 1.038 0.299137
## schID47 0.859422 0.354226 2.426 0.015308 *
## schID48 0.283722 0.374613 0.757 0.448877
## schID49 -0.120278 0.339956 -0.354 0.723507
## schID50 0.231275 0.311965 0.741 0.458532
## schID51 0.388530 0.298161 1.303 0.192631
## schID52 0.348340 0.317809 1.096 0.273125
## schID53 0.320064 0.311965 1.026 0.304980
## schID54 0.090972 0.339956 0.268 0.789023
## schID55 -0.298325 0.277947 -1.073 0.283203
## schID56 -0.456444 0.406325 -1.123 0.261367
## schID57 0.979097 0.301339 3.249 0.001168 **
## schID58 0.752222 0.363428 2.070 0.038545 *
## schID59 0.263508 0.329385 0.800 0.423765
## schID60 0.059639 0.339956 0.175 0.860751
## schID61 0.102937 0.388521 0.265 0.791068
## schID62 -0.006178 0.354226 -0.017 0.986086
## schID63 0.522556 0.339956 1.537 0.124353
## schID64 0.510294 0.329385 1.549 0.121416
## schID65 0.541972 0.374613 1.447 0.148056
## schID66 0.247389 0.339956 0.728 0.466841
## schID67 0.489677 0.346515 1.413 0.157702
## schID68 0.852556 0.363428 2.346 0.019038 *
## schID69 0.336163 0.317809 1.058 0.290241
## schID70 0.444222 0.374613 1.186 0.235776
## schID71 0.594222 0.430014 1.382 0.167101
## schID72 0.045657 0.303120 0.151 0.880281
## schID73 0.077079 0.307153 0.251 0.801869
## schID74 0.652556 0.339956 1.920 0.054999 .
## schID75 0.346915 0.334305 1.038 0.299473
## schID76 -0.005778 0.354226 -0.016 0.986987
## schID77 0.379607 0.334305 1.136 0.256240
## schID78 0.273738 0.291913 0.938 0.348442
## schID79 0.243060 0.286538 0.848 0.396348
## schID80 -0.140159 0.307153 -0.456 0.648191
## schID81 0.649667 0.363428 1.788 0.073926 .
## schID82 -0.236194 0.339956 -0.695 0.487240
## schID83 -0.346232 0.346515 -0.999 0.317775
## schID84 0.129222 0.339956 0.380 0.703883
## schID85 -0.266515 0.311965 -0.854 0.392992
## schID86 0.510660 0.321228 1.590 0.111991
## schID87 -0.067361 0.339956 -0.198 0.842942
## schID88 0.043365 0.307153 0.141 0.887733
## schID89 0.535722 0.321228 1.668 0.095459 .
## schID90 0.289778 0.363428 0.797 0.425305
## schID91 0.237530 0.334305 0.711 0.477430
## schID92 0.706972 0.339956 2.080 0.037635 *
## schID93 0.174741 0.296738 0.589 0.555985
## schID94 0.737472 0.339956 2.169 0.030126 *
## schID95 -0.028992 0.329385 -0.088 0.929867
## schID96 0.460422 0.354226 1.300 0.193756
## schID97 0.869056 0.339956 2.556 0.010619 *
## schID98 0.783508 0.388521 2.017 0.043809 *
## schID99 0.081583 0.287315 0.284 0.776465
## schID100 0.510313 0.305052 1.673 0.094441 .
## genderM 0.439079 0.329385 1.333 0.182609
## schID2:genderM -1.075802 0.477876 -2.251 0.024434 *
## schID3:genderM -0.380424 0.452749 -0.840 0.400824
## schID4:genderM -0.739919 0.426379 -1.735 0.082766 .
## schID5:genderM -1.098079 0.497170 -2.209 0.027263 *
## schID6:genderM -0.813253 0.395331 -2.057 0.039746 *
## schID7:genderM -0.526365 0.447748 -1.176 0.239842
## schID8:genderM -1.426051 0.478027 -2.983 0.002872 **
## schID9:genderM -0.321819 0.418018 -0.770 0.441429
## schID10:genderM -0.728987 0.424976 -1.715 0.086367 .
## schID11:genderM -0.161437 0.417999 -0.386 0.699363
## schID12:genderM -0.869225 0.412804 -2.106 0.035305 *
## schID13:genderM -0.828559 0.399105 -2.076 0.037963 *
## schID14:genderM -1.168274 0.473354 -2.468 0.013632 *
## schID15:genderM -1.382979 0.444577 -3.111 0.001881 **
## schID16:genderM -0.797329 0.444577 -1.793 0.072986 .
## schID17:genderM -0.584070 0.440238 -1.327 0.184690
## schID18:genderM -0.668191 0.396480 -1.685 0.092018 .
## schID19:genderM -0.902226 0.426885 -2.114 0.034628 *
## schID20:genderM -0.753731 0.377044 -1.999 0.045680 *
## schID21:genderM -0.906346 0.420465 -2.156 0.031184 *
## schID22:genderM -1.101276 0.384390 -2.865 0.004195 **
## schID23:genderM -0.897721 0.369542 -2.429 0.015179 *
## schID24:genderM -1.025582 0.440271 -2.329 0.019893 *
## schID25:genderM -1.088150 0.451381 -2.411 0.015973 *
## schID26:genderM -0.629647 0.348104 -1.809 0.070569 .
## schID27:genderM -0.896862 0.434945 -2.062 0.039280 *
## schID28:genderM -0.451068 0.426885 -1.057 0.290744
## schID29:genderM -1.101844 0.444398 -2.479 0.013207 *
## schID30:genderM -1.067488 0.471130 -2.266 0.023524 *
## schID31:genderM -1.587688 0.451381 -3.517 0.000441 ***
## schID32:genderM -1.037495 0.377893 -2.745 0.006073 **
## schID33:genderM -0.868211 0.442207 -1.963 0.049684 *
## schID34:genderM -0.963636 0.398875 -2.416 0.015748 *
## schID35:genderM -0.894396 0.434945 -2.056 0.039823 *
## schID36:genderM -0.893416 0.366340 -2.439 0.014787 *
## schID37:genderM -0.705846 0.455581 -1.549 0.121393
## schID38:genderM -0.881093 0.464553 -1.897 0.057957 .
## schID39:genderM -0.835424 0.452749 -1.845 0.065089 .
## schID40:genderM -0.536441 0.462679 -1.159 0.246362
## schID41:genderM -1.032125 0.471130 -2.191 0.028536 *
## schID42:genderM -0.773930 0.444398 -1.742 0.081680 .
## schID43:genderM -0.398537 0.442982 -0.900 0.368359
## schID44:genderM -0.466294 0.439776 -1.060 0.289082
## schID45:genderM -0.671541 0.478027 -1.405 0.160164
## schID46:genderM -0.608713 0.421707 -1.443 0.148983
## schID47:genderM -1.304126 0.462223 -2.821 0.004808 **
## schID48:genderM -0.468579 0.478027 -0.980 0.327037
## schID49:genderM -0.914496 0.455581 -2.007 0.044793 *
## schID50:genderM -0.571372 0.404414 -1.413 0.157792
## schID51:genderM -0.827923 0.390617 -2.120 0.034116 *
## schID52:genderM -1.346157 0.408939 -3.292 0.001005 **
## schID53:genderM -1.050614 0.430693 -2.439 0.014763 *
## schID54:genderM -0.574529 0.433293 -1.326 0.184939
## schID55:genderM -0.695967 0.364183 -1.911 0.056082 .
## schID56:genderM -0.389886 0.488709 -0.798 0.425047
## schID57:genderM -0.811954 0.389658 -2.084 0.037254 *
## schID58:genderM -0.752913 0.473354 -1.591 0.111792
## schID59:genderM -0.607365 0.443473 -1.370 0.170911
## schID60:genderM -0.856223 0.460497 -1.859 0.063062 .
## schID61:genderM -0.868579 0.485651 -1.788 0.073785 .
## schID62:genderM -0.468367 0.452855 -1.034 0.301089
## schID63:genderM -0.824231 0.460497 -1.790 0.073561 .
## schID64:genderM -0.932951 0.425049 -2.195 0.028234 *
## schID65:genderM -0.879123 0.466639 -1.884 0.059655 .
## schID66:genderM -0.672913 0.473354 -1.422 0.155236
## schID67:genderM -1.373677 0.452749 -3.034 0.002430 **
## schID68:genderM -0.859623 0.453670 -1.895 0.058199 .
## schID69:genderM -0.712754 0.427879 -1.666 0.095846 .
## schID70:genderM -1.410233 0.478027 -2.950 0.003198 **
## schID71:genderM -0.485785 0.512180 -0.948 0.342958
## schID72:genderM -0.755969 0.384505 -1.966 0.049368 *
## schID73:genderM -1.506937 0.414439 -3.636 0.000281 ***
## schID74:genderM -1.127975 0.441783 -2.553 0.010715 *
## schID75:genderM -0.800301 0.434945 -1.840 0.065852 .
## schID76:genderM -0.841352 0.471130 -1.786 0.074215 .
## schID77:genderM -0.918082 0.414362 -2.216 0.026779 *
## schID78:genderM -0.843295 0.373570 -2.257 0.024045 *
## schID79:genderM -0.770765 0.370782 -2.079 0.037712 *
## schID80:genderM -0.546763 0.394931 -1.384 0.166308
## schID81:genderM -1.532857 0.462773 -3.312 0.000935 ***
## schID82:genderM -1.338218 0.473354 -2.827 0.004724 **
## schID83:genderM -1.512425 0.471130 -3.210 0.001338 **
## schID84:genderM -0.962767 0.441783 -2.179 0.029378 *
## schID85:genderM -0.712161 0.408403 -1.744 0.081287 .
## schID86:genderM -1.487384 0.430425 -3.456 0.000556 ***
## schID87:genderM -0.960829 0.455581 -2.109 0.035014 *
## schID88:genderM -0.714572 0.408063 -1.751 0.080011 .
## schID89:genderM -0.871198 0.417067 -2.089 0.036792 *
## schID90:genderM -0.477635 0.462773 -1.032 0.302089
## schID91:genderM -1.546054 0.469313 -3.294 0.000997 ***
## schID92:genderM -0.531488 0.415356 -1.280 0.200773
## schID93:genderM -0.902646 0.398510 -2.265 0.023571 *
## schID94:genderM -1.017329 0.466327 -2.182 0.029207 *
## schID95:genderM -0.546683 0.452749 -1.207 0.227330
## schID96:genderM -1.044221 0.450436 -2.318 0.020493 *
## schID97:genderM -0.549463 0.433293 -1.268 0.204843
## schID98:genderM -1.083626 0.468236 -2.314 0.020710 *
## schID99:genderM -0.355544 0.381444 -0.932 0.351349
## schID100:genderM -1.550863 0.425712 -3.643 0.000273 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 0.7709 on 3500 degrees of freedom
## Multiple R-squared: 0.2405, Adjusted R-squared: 0.1973
## F-statistic: 5.57 on 199 and 3500 DF, p-value: < 2.2e-16
注意一下截距與斜率
rstmix <- lmer(sciscore~1+gender+(1+gender|schID),data=dta)
summary(rstmix)
## Linear mixed model fit by REML. t-tests use Satterthwaite's method [
## lmerModLmerTest]
## Formula: sciscore ~ 1 + gender + (1 + gender | schID)
## Data: dta
##
## REML criterion at convergence: 8810
##
## Scaled residuals:
## Min 1Q Median 3Q Max
## -3.9698 -0.6464 0.0519 0.6983 2.6584
##
## Random effects:
## Groups Name Variance Std.Dev. Corr
## schID (Intercept) 0.10330 0.3214
## genderM 0.01994 0.1412 0.20
## Residual 0.59449 0.7710
## Number of obs: 3700, groups: schID, 100
##
## Fixed effects:
## Estimate Std. Error df t value Pr(>|t|)
## (Intercept) 6.40435 0.03795 97.92954 168.76 <2e-16 ***
## genderM -0.39534 0.03020 76.38063 -13.09 <2e-16 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Correlation of Fixed Effects:
## (Intr)
## genderM -0.262
summary(IntSlope)
## Intercept genderSlope
## Min. :5.363 Min. :-1.1486
## 1st Qu.:6.179 1st Qu.:-0.6001
## Median :6.423 Median :-0.4107
## Mean :6.416 Mean :-0.4211
## 3rd Qu.:6.660 3rd Qu.:-0.2195
## Max. :7.500 Max. : 0.4391
var(IntSlope)
## Intercept genderSlope
## Intercept 0.14012705 -0.02958089
## genderSlope -0.02958089 0.10575264
newdta1 <- data.frame(schID='50',gender='M')
predict(rstfix,newdta1)
## 1
## 6.27976
predict(rstmix,newdta1,re.form=~(1+gender|schID))
## 1
## 6.191668
IntSlope[50,1] + IntSlope[50,2]*1
## [1] 6.27976
newdta2 <- data.frame(schID='101',gender='M')
predict(rstfix,newdta2)
## Error in model.frame.default(Terms, newdata, na.action = na.action, xlev = object$xlevels): 因子 schID 裡出現了新的層次 101
predict(rstmix,newdta2,allow.new.levels = TRUE,re.form=~(1+gender|schID))
## 1
## 6.009013