Assessing homogeneity of variances assumption for numeric variables

Question

I'm fitting a model with a 2-level categorical and a numerical predictor. As far as I'm concerned, I cannot use car::LeveneTest() to assess homoscedasticity since I have a numerical predictor, right? What should I use instead to assess homoscedasticity on both numerical and categorical predictors?

• expected error:

car::leveneTest(lmer(SCORE ~ EXAM_GRADE * TASK_TYPE + (1|ID), data = dfModels1, REML = F))

Error in leveneTest.default(lmer(SCORE ~ EXAM_GRADE  * TASK_TYPE + (1 |  : 
  lmer(SCORE ~ EXAM_GRADE  * TASK_TYPE + (1 | ID), data = dfModels1, REML = F) is not a numeric variable

• head(dfModels1):

Each ID has a SCORE on two different TASKs. They also have a score on another Exam Grade. The basic idea is that Exam Grade predicts SCORE independenly of the task type.

ID      SCORE    EXAM_GRADE    TASK_TYPE    
  <chr>  <dbl>    <dbl>            <fct>         
1 1_101   7.1      100             A          
2 1_101   4.92     100             N   
1 2_101   7.4      36              A         
2 2_101   7        36              N   

• Maybe the Breusch-Pagan test (such as performance::check_heteroscedasticity() does?)

• Obs: I'm aware that I should prefer visual inspection instead, I'm also doing that, but I'd like to find a way to assess the assumption numerically as well just in case.

Answer 1

"Detecting heteroscedasticity" is a kind of magical thinking.

As you might expect, there's only one way that variance can be homoscedastic in a sample, but there are infinite ways it may be heteroscedastic. The assumption that heteroscedasticity would be apparent in measured covariates is tenuous enough. More realistically, there are exogenous trends that don't show up in the covariates but do impact analyses.

Add the issue of calibrating an actual hypothesis test and identifying a suitable number of subgroups, it quickly slips out of hand. Consider the large sample size dilemma - at large enough $N$ virtually every assumption-based test fails, but this barely precludes the value of the analysis. Essentially, you punish yourself for having good data.

I often find I can simulate data with a known failure of assumptions yet with standard parametric models, I obtain reasonable precise estimates, well conserved interval estimates, and powerful and well calibrated tests.

Does this mean one should not assess their assumptions? Clearly not. But somehow in the statistics pedagogy, "assessing assumptions" has been equated to conducting tests, which apropos of nothing we can't rely on those $p$-values at all.

Infinitely more valuable are the residual plots - residual versus covariate, and residual versus fitted, residual versus leverage, and so on. Correlation between sites or over time can be compared using the ICC or variograms (note: correlation functionally is similar to heteroscedasticity). On the chance that variance does depend on any such feature, a visual or quantitative summary should suffice. I've never had to report the actual results in the form of a statistical test for any published research manuscript.

If the goal of model fitting is inference rather than prediction, the robust variance estimate, also called Huber White or Sandwich variance errors, provide accurate 95% CIs and p-values even in the presence of heteroscedasticity, at barely any impact to power. I would say that this issue has been almost completely resolved as a consequence of these error estimators. In stata, you enter , robust=TRUE. In R, using the sandwich and lmtest package, you can enter waldtest with vcov=VCOVHC to obtain robust errors.