PRE

library(Keng)
library(effectsize)
library(car)
#> Loading required package: carData
data("depress")

PRE is called partial R-squared in regression, and partial Eta-squared in ANOVA. This vignette will examine their equivalence using the internal data depress.

depress collected depression, gender, and class at Time 1. Traditionally, we examine the effect of gender and class using anova. We firstly let R know gender and class are factors (i.e., categorical variables). Then we conduct anova using car::Anova() and compute partial Eta-squared using effectsize::eta_squared().

# factor gender and class
depress_factor <- depress
depress_factor$class <- factor(depress_factor$class, labels = 1:4)
depress_factor$gender <- factor(depress_factor$gender, labels = c(0,1))

anova.fit <- lm(depr1 ~ gender + class, depress_factor)
Anova(anova.fit, type = 3)
#> Anova Table (Type III tests)
#> 
#> Response: depr1
#>              Sum Sq  Df   F value  Pr(>F)    
#> (Intercept) 142.363   1 1040.8995 < 2e-16 ***
#> gender        0.134   1    0.9811 0.32335    
#> class         1.461   3    3.5603 0.01554 *  
#> Residuals    23.114 169                      
#> ---
#> Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
cat("\n\n")
print(eta_squared(Anova(anova.fit, type = 3), partial = TRUE), digits = 6)
#> # Effect Size for ANOVA (Type III)
#> 
#> Parameter | Eta2 (partial) |               95% CI
#> -------------------------------------------------
#> gender    |       0.005772 | [0.000000, 1.000000]
#> class     |       0.059444 | [0.006571, 1.000000]
#> 
#> - One-sided CIs: upper bound fixed at [1.000000].

Then we conduct regression analysis and compute PRE. For class with four levels: 3, 5, 9, and 12, we dummy-code it using ifelse() with the class12 as the reference group.

# class3 indicates whether the class is class3 
depress$class3 <- ifelse(depress$class == 1, 1, 0)
# class5 indicates whether the class is class5 
depress$class5 <- ifelse(depress$class == 2, 1, 0)
# class9 indicates whether the class is class9 
depress$class9 <- ifelse(depress$class == 3, 1, 0)

We compute the PRE of gender though comparing Model A with gender against Model C without gender.

fitC <- lm(depr1 ~ class3 + class5 + class9, depress)
fitA <- lm(depr1 ~ class3 + class5 + class9 + gender, depress)
print(compare_lm(fitC, fitA), digits = 3)
#>                      Baseline        C        A   A vs. C
#> SSE                      24.7  23.2482  23.1140  1.34e-01
#> n                       174.0 174.0000 174.0000  1.74e+02
#> Number of parameters      1.0   4.0000   5.0000  1.00e+00
#> df                      173.0 170.0000 169.0000  1.00e+00
#> R_squared                  NA   0.0581   0.0635  5.44e-03
#> f_squared                  NA   0.0616   0.0678  5.81e-03
#> R_squared_adj              NA   0.0414   0.0413        NA
#> PRE                        NA   0.0581   0.0635  5.77e-03
#> F(PA-PC,n-PA)              NA   3.4929   2.8646  9.81e-01
#> p                          NA   0.0170   0.0249  3.23e-01
#> PRE_adj                    NA   0.0414   0.0413 -1.11e-04
#> power_post                 NA   0.7720   0.7683  1.66e-01

Compare gender’s PRE and partial Eta-squared. They should be equal.

We compute the PRE of class. Note that in regression, the PRE of class is the PRE of all class’s dummy codes: class3, class5, and class9.

fitC <- lm(depr1 ~ gender, depress)
fitA <- lm(depr1 ~ class3 + class5 + class9 + gender, depress)
print(compare_lm(fitC, fitA), digits = 3)
#>                      Baseline         C        A  A vs. C
#> SSE                      24.7  24.57485  23.1140   1.4608
#> n                       174.0 174.00000 174.0000 174.0000
#> Number of parameters      1.0   2.00000   5.0000   3.0000
#> df                      173.0 172.00000 169.0000   3.0000
#> R_squared                  NA   0.00431   0.0635   0.0592
#> f_squared                  NA   0.00433   0.0678   0.0632
#> R_squared_adj              NA  -0.00148   0.0413       NA
#> PRE                        NA   0.00431   0.0635   0.0594
#> F(PA-PC,n-PA)              NA   0.74436   2.8646   3.5603
#> p                          NA   0.38947   0.0249   0.0155
#> PRE_adj                    NA  -0.00148   0.0413   0.0427
#> power_post                 NA   0.13765   0.7683   0.7806

Compare class’s PRE and partial Eta-squared. They should be equal.

We compute the PRE of the full model(Model A). The PRE (partial R-squared or partial Eta-squared) of the full model is commonly known as the R-squared or Eta-squared of the full model.

fitC <- lm(depr1 ~ 1, depress)
fitA <- lm(depr1 ~ class3 + class5 + class9 + gender, depress)
print(compare_lm(fitC, fitA), digits = 3)
#>                      Baseline     C        A  A vs. C
#> SSE                      24.7  24.7  23.1140   1.5672
#> n                       174.0 174.0 174.0000 174.0000
#> Number of parameters      1.0   1.0   5.0000   4.0000
#> df                      173.0 173.0 169.0000   4.0000
#> R_squared                  NA   0.0   0.0635   0.0635
#> f_squared                  NA   0.0   0.0678   0.0678
#> R_squared_adj              NA   0.0   0.0413       NA
#> PRE                        NA   0.0   0.0635   0.0635
#> F(PA-PC,n-PA)              NA    NA   2.8646   2.8646
#> p                          NA    NA   0.0249   0.0249
#> PRE_adj                    NA   0.0   0.0413   0.0413
#> power_post                 NA    NA   0.7683   0.7683

As shown, the PRE of Model A against Model C is equal to Model A’s R_squared. Taken the loss of precision into consideration, Model C’s R_squared is zero.