1
$\begingroup$

Background

In R, the power analysis package pwr can estimate the required sample size for regression with the following formula:

pwr.f2.test(u = , v = , f2 = , sig.level = , power = )

According to the pwr vignette (https://cran.r-project.org/web/packages/pwr/vignettes/pwr-vignette.html):

The numerator degrees of freedom, u, is the number of coefficients you'll have in your model (minus the intercept).

In R, we can therefore obtain u with the following formula:

u = length(coef(mod))-1

Let's see a simple example:

# Make a simple linear model
mod <- lm(mpg ~ cyl*disp, mtcars)

# Determine u
(u = length(coef(mod))-1)
[1] 3

# Which can be confirmed by looking at the model summary
summary(mod)
Coefficients:
             Estimate Std. Error t value Pr(>|t|)    
(Intercept) 49.037212   5.004636   9.798 1.51e-10 ***
cyl         -3.405244   0.840189  -4.053 0.000365 ***
disp        -0.145526   0.040002  -3.638 0.001099 ** 
cyl:disp     0.015854   0.004948   3.204 0.003369 ** 

# determine v
library(pwr)
(v <- pwr.f2.test(u = u, v = , f2 = 0.25, sig.level = .05, power = .80)$v)
[1] 43.70444

According to the package vignette:

The denominator degrees of freedom, v, is the number of error degrees of freedom: v = n − u − 1. This implies n = v + u + 1.

So the required sample size can be computed as v + u + 1.

# Determine required sample size
ceiling(v + u + 1)
[1] 48

And I can confirm R's results in G*Power when I put "Number of tested predictors" as u (number of coefficients -1):

enter image description here

Notice as well that the denominator df matches v, as expected.

Problem

However, my PI and my university statistician argue that in G*Power, the number of tested coefficients should not be the total number of coefficients minus the intercept, but only 1, as typically you will only be interested in the highest-level interaction. Therefore, it should look like this and the required sample size should only be 34:

enter image description here

We can replicate this in R as well:

# Determine u
u = 1

# determine v
library(pwr)
(v <- pwr.f2.test(u = u, v = , f2 = 0.25, sig.level = .05, power = .80)$v)
[1] 31.42944

(Note that our new v surprisingly does not match the G*Power denominator df exactly, although the difference is small.)

# Determine required sample size
ceiling(v + u + 1)
[1] 34

We do get the same results as G*Power again.

Question

Am I misinterpreting the pwr vignette? Should the numerator degrees of freedom be the number of coefficients minus the intercept or simply 1?

Edit: Another area of confusion is that in G*Power, changing the total number of predictors will change the denominator df, therefore (slightly) changing the required sample size. However, I do not see in pwr how to change the denominator df because it gives it to you directly (as v) without asking anywhere for the total number of predictors.

$\endgroup$

1 Answer 1

1
$\begingroup$

Your PI and you are conflating two tests. The F-test compares nested models. When you test for all variables against an "intercept-only" model, you get numerator DF = number of covariables inc. intercept minus 1. Unintentionally, it seems, this is what you have found here. The null hypothesis is that NONE of the variables are associated of the outcome.

If you test the statistical significance of a single covariate in the presence of other confounds, you test against a null model which adjusts for all variables except the covariate of interest. In that case, num DF = 1 (or more for categorical covariates using dummy encoding, depending the number of levels).

Suppose I have some covariable, like smoking, that I want to test as a risk factor for low-birthweight (measured continuously). Suppose further I adjust for maternal age and maternal weight prior to conception as possible confounds. Say I get 4 variables (including smoking status, and intercept) into the final model.

My scientific question is whether smoking has a statistically significant effect after adjustment for the other variables.

My null model in this case adjusts for all the other confounds as somewhat obvious predictors of birth weight that may be causal of smoking status. If I test all variables, I find out (indirectly) if my 4 variables form a kind of predictive model, but nothing specifically about effect of smoking.

Anyway, returning to the issue of the present analysis. Recall an F distribution with one numerator degree of freedom is the Chi-square distribution. In other words, the test simplifies to the Wald test for the statistical significance of a single covariable in a linear regression model. To remedy your approach to power/sample size calculation, you need only consider one factor in addition to resetting the numerator DF. That is the residual error of the response. When considering the other adjustments, the relative decrease in MSE is much smaller than compared to the intercept only model. You may need to "Return to the drawing board" regarding the literature to understand what effect sizes you can be powered to detect. In the pregnancy outcome, it would not suffice merely to know the variability in birthweight. You need the residual variability after adjustment for the confounds as the residual standard error in the "null" model.

$\endgroup$
5
  • $\begingroup$ So in the case I want to test a single interaction term and not the whole model, the second approach with 1 numerator df is correct, right? Then where do I specify the total number of terms like in G-Power? This seems to change the required sample size in G-Power but I do not see where to add this information in pwr.f2.test. $\endgroup$
    – rempsyc
    May 6, 2021 at 20:21
  • $\begingroup$ @RemPsyc the other terms don't matter, but I did neglect to mention that you need to change both the numerator and denominator degrees of freedom. So Num DF is set to 1. Den DF is set to n-p+1 with p the number of parameters in the biggest model. So for my birthweight example, if I had n=100 and 4 terms (including smoking and intercept) in the biggest model, the denominator DF is 97. $\endgroup$
    – AdamO
    May 6, 2021 at 21:15
  • $\begingroup$ Ok... but how do I determine Den DF if I do not know the sample size? Since this is what I am trying to estimate. $\endgroup$
    – rempsyc
    May 6, 2021 at 22:43
  • $\begingroup$ @RemPsyc As you noted G*power and R give you the Den DF as a result of the power and sample size calculation. You add a person for each parameter in the model. $\endgroup$
    – AdamO
    May 6, 2021 at 22:53
  • $\begingroup$ I am not clear that this is the case as G-Power lets you input the total number of predictors and it does not add 1 to the required sample size per added total number of predictors. Furthermore, that would mean that the pwr is inaccurate in this regard? $\endgroup$
    – rempsyc
    May 6, 2021 at 22:57

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service and acknowledge you have read our privacy policy.

Not the answer you're looking for? Browse other questions tagged or ask your own question.