> print("R Descriptive Statistics")
[1] "R Descriptive Statistics"
> describe(x = Lecture2[, c("langpost", "homework", "mixgrd")])
Lecture2[, c("langpost", "homework", "mixgrd")] 

 3  Variables      3566  Observations
--------------------------------------------------------------------------------------------------------------------------------------------
langpost 
       n  missing distinct     Info     Mean  pMedian      Gmd      .05      .10      .25      .50      .75      .90      .95 
    3566        0       49    0.999    41.46       42    9.936       25       29       36       42       48       52       54 

lowest :  9 11 12 13 14, highest: 54 55 56 57 58
--------------------------------------------------------------------------------------------------------------------------------------------
homework 
       n  missing distinct     Info     Mean  pMedian      Gmd 
    3566        0        9    0.964    2.443      2.5   0.5933 
                                                                
Value      1.333 1.667 2.000 2.333 2.667 3.000 3.333 3.500 3.667
Frequency    190   306   542   790   887   610   163    21    57
Proportion 0.053 0.086 0.152 0.222 0.249 0.171 0.046 0.006 0.016
--------------------------------------------------------------------------------------------------------------------------------------------
mixgrd 
       n  missing distinct     Info      Sum     Mean 
    3566        0        2    0.549      860   0.2412 

--------------------------------------------------------------------------------------------------------------------------------------------
> print("Single-Level Empty Means, No Random Intercept")
[1] "Single-Level Empty Means, No Random Intercept"
> Model1 = lm(data = Lecture2, formula = langpost ~ 1)
> summary(Model1)

Call:
lm(formula = langpost ~ 1, data = Lecture2)

Residuals:
   Min     1Q Median     3Q    Max 
-32.46  -5.46   0.54   6.54  16.54 

Coefficients:
            Estimate Std. Error t value Pr(>|t|)
(Intercept)   41.464      0.148     281   <2e-16

Residual standard error: 8.81 on 3565 degrees of freedom

> print("Get total variance")
[1] "Get total variance"
> summary(Model1)$sigma^2
[1] 77.691
> print("Two-Level Empty Means, Random Intercept")
[1] "Two-Level Empty Means, Random Intercept"
> Model2 = lmer(data = Lecture2, REML = TRUE, formula = langpost ~ 
+ 1 + (1 | schoolID))
> summary(Model2)
Linear mixed model fit by REML. t-tests use Satterthwaite's method ['lmerModLmerTest']
Formula: langpost ~ 1 + (1 | schoolID)
   Data: Lecture2

REML criterion at convergence: 25191

Scaled residuals: 
   Min     1Q Median     3Q    Max 
-4.208 -0.646  0.096  0.739  2.491 

Random effects:
 Groups   Name        Variance Std.Dev.
 schoolID (Intercept) 17.8     4.22    
 Residual             62.2     7.89    
Number of obs: 3566, groups:  schoolID, 191

Fixed effects:
            Estimate Std. Error      df t value Pr(>|t|)
(Intercept)   41.079      0.337 174.914     122   <2e-16
> print("Show ICC1 and its LRT")
[1] "Show ICC1 and its LRT"
> performance::icc(Model2)
# Intraclass Correlation Coefficient

    Adjusted ICC: 0.223
  Unadjusted ICC: 0.223
> lmerTest::ranova(Model2)
ANOVA-like table for random-effects: Single term deletions

Model:
langpost ~ 1 + (1 | schoolID)
               npar logLik   AIC LRT Df Pr(>Chisq)
<none>            3 -12595 25197                  
(1 | schoolID)    2 -12821 25647 452  1     <2e-16
> print("Show ICC2 using psychometric package")
[1] "Show ICC2 using psychometric package"
> psychometric::ICC2.lme(data = Lecture2, dv = langpost, grp = schoolID, 
+ weighted = TRUE)
[1] 0.84647
> print("Add 2 School-Level Predictors")
[1] "Add 2 School-Level Predictors"
> Model3 = lmer(data = Lecture2, REML = TRUE, formula = langpost ~ 
+ 1 + hw2 + mixgrd + (1 | schoolID))
> summary(Model3)
Linear mixed model fit by REML. t-tests use Satterthwaite's method ['lmerModLmerTest']
Formula: langpost ~ 1 + hw2 + mixgrd + (1 | schoolID)
   Data: Lecture2

REML criterion at convergence: 25180

Scaled residuals: 
   Min     1Q Median     3Q    Max 
-4.170 -0.648  0.094  0.736  2.532 

Random effects:
 Groups   Name        Variance Std.Dev.
 schoolID (Intercept) 17.2     4.14    
 Residual             62.2     7.89    
Number of obs: 3566, groups:  schoolID, 191

Fixed effects:
            Estimate Std. Error      df t value Pr(>|t|)
(Intercept)   41.548      0.505 167.500   82.26   <2e-16
hw2            0.507      0.635 172.279    0.80   0.4261
mixgrd        -1.993      0.708 188.544   -2.81   0.0054

Correlation of Fixed Effects:
       (Intr) hw2   
hw2    -0.589       
mixgrd -0.548  0.139
> Lecture2$Pred1c = predict(Model3, re.form = NA)
> rModel3 = cor.test(Lecture2$Pred1c, Lecture2$langpost, method = "pearson")
> print("Total R2")
[1] "Total R2"
> rModel3$estimate^2
      cor 
0.0062052 
> print("Show conditional ICC1 and its LRT")
[1] "Show conditional ICC1 and its LRT"
> performance::icc(Model3)
# Intraclass Correlation Coefficient

    Adjusted ICC: 0.216
  Unadjusted ICC: 0.214
> lmerTest::ranova(Model3)
ANOVA-like table for random-effects: Single term deletions

Model:
langpost ~ 1 + hw2 + mixgrd + (1 | schoolID)
               npar logLik   AIC LRT Df Pr(>Chisq)
<none>            5 -12590 25190                  
(1 | schoolID)    4 -12811 25630 442  1     <2e-16
> print("Add 2 School-Level Predictors \n
+   WITHOUT Random Intercept")
[1] "Add 2 School-Level Predictors \n      WITHOUT Random Intercept"
> Model4 = lm(data = Lecture2, formula = langpost ~ 1 + hw2 + mixgrd)
> summary(Model4)

Call:
lm(formula = langpost ~ 1 + hw2 + mixgrd, data = Lecture2)

Residuals:
    Min      1Q  Median      3Q     Max 
-31.437  -5.911   0.952   6.837  17.426 

Coefficients:
            Estimate Std. Error t value Pr(>|t|)
(Intercept)   41.637      0.216  192.50  < 2e-16
hw2            0.411      0.280    1.47     0.14
mixgrd        -1.474      0.347   -4.25 0.000022

Residual standard error: 8.79 on 3563 degrees of freedom
Multiple R-squared:  0.00621,	Adjusted R-squared:  0.00565 
F-statistic: 11.1 on 2 and 3563 DF,  p-value: 0.0000152

> print("Residual variance")
[1] "Residual variance"
> summary(Model4)$sigma^2
[1] 77.251
