Which power test should I use?

[Nate_L]

I'm trying to do a simple calculation of power for future data, but recently learned that both the stats and pwrss package have a test with the same name. They produce different power calculations though, so I need help deciding which one fits my problem best.

I want to calculate the power to detect a 30% (hypothetical) decline in density (count/sample area), but don't fully understand either function. I'm also probably getting the "one.sample"/"less than" jargon mixed up...(any help there is also appreciated).

stats::power.t.test(n = 47, 
                      sd = 0.4352092, 
                      delta = 0.1139639, 
                      alternative = "one.sided", 
                      type = "one.sample")

pwrss::power.t.test(ncp = (0.2659157*sqrt(47)), # <- tells the direction! (positive must = "greater"!)
                    #  + ...
                    #  + alternative = "less",
                    # Warning message:
                    #   alternative = 'less' but non-centrality parameter is positive 
             df = 46,
             alpha = 0.05,
             alternative = "greater",
             plot = FALSE)

Source for ncp definition: statistical power - Noncentrality Parameter - what is it, what does it do, what would be a suggested value? - Cross Validated

Data:
N=47 (Same 47 sites each time = paired test)
mean_log_den = mean (log transformed) density of 47 sites
sd = standard deviation
ln_down = hypothetical 30% decline of mean (log transformed) density
delta = change in delta

structure(list(Season = "DRY", assem = "Far", n = 47L, mean_log_den = 0.379879532043208, 
    sd = 0.435209207239454, ln_down = 0.265915672430246, delta = 0.113963859612963), class = c("grouped_df", 
"tbl_df", "tbl", "data.frame"), row.names = c(NA, -1L), groups = structure(list(
    Season = "DRY", .rows = structure(list(1L), ptype = integer(0), class = c("vctrs_list_of", 
    "vctrs_vctr", "list"))), class = c("tbl_df", "tbl", "data.frame"
), row.names = c(NA, -1L), .drop = TRUE))

[MetinBulus]

Hi @Nate_L,

They should produce the same result. However, the non-centrality parameter should take into account the correlation between the dependent measures. By default, R codes assume a correlation of 0.50 (in stats::power.t.test() and pwrss::pwrss.t.2means() functions).

pwrss::power.t.test() function is a generic function. It can be used to calculate power for any type of t test (e.g., independent and dependent t tests, as well as to test regression coefficient) as long as the non-centrality parameter is provided.

Using all three approaches, please find corrected calculations below:

# make sure to use paired test
stats::power.t.test(n = 47, 
                    sd = 0.4352092, 
                    delta = 0.1139639, 
                    alternative = "one.sided", 
                    type = "paired")

# using the specific function tailored for dependent and independent t tests
pwrss::pwrss.t.2means(mu1 = 0.1139639, # mu2 = 0,
                      sd1 = 0.4352092, # sd2 = 0.4352092,
                      paired = TRUE, n = 47, 
                      alpha = 0.05, alternative = "greater")

# using the generic function
n <- 47
d <- 0.1139639
sd <- 0.4352092
cor <- 0.50 # correlation between repeated measures
lambda <- d * sqrt(n) / sqrt(2*sd^2 * (1 - cor))

pwrss::power.t.test(ncp = lambda, 
                    df = n - 1,
                    alpha = 0.05,
                    alternative = "greater",
                    plot = FALSE)

The "less", "greater", "not equal" options relates to how you form your alternative hypothesis in relation to the null hypothesis based on the theory. If the theory suggest that the mean for the first time point is less than the mean for the second time point alternative = 'less' should be used (one-sided). In contrast, if the theory suggests that the mean for the first time point is greater than the mean for the second time point alternative = 'greater' should be used (one-sided). In cases where the theory is not clear about the direction of the difference alternative = 'not equal' should be used (two-sided).

I hope this helps.
All the best.

P.S. Edited on 2023.11.08 based on the comment received from Daniel Lakens

Daniel Lakens (@lakens.bsky.social)

[Nate_L]

Thank you very much!! Just a quick follow-up question, if I may! Thinking on it further, the assumption that the standard deviation and/or correlation between between repeated measures will remain the same might be too strong (this is fish density and it is really all over the place). Would you recommend that I stick with the first approach (stats::power.t.test)?

Also, I have to do many of these power analyses (i.e. what would the power be for species A, B, C, and what would happen is we doubled, or tripled the sample size?), do I need to adjust alpha in that case (aplha/the number of tests I do to increase N)?

[MetinBulus]

If I understood the problem correctly, in your case, I believe it might be better to calculate power for a range of scenarios. For example, you may want to consider a range of standard deviations for various time points. You may also want to consider a range of correlations between pairs of time points. They may or may not be the same for species A, B, and C.

You have all three options available to calculate power for each scenario; (i) calculate the non-centrality parameter for each scenario and use pwrss::power.t.test() function, (ii) calculate the standard deviation of the difference for each scenario and use stats::power.t.test() function, (iii) specify standard deviations for each time point and the correlation between them for each scenario and use the pwrss::pwrss.t.2means() function.

Regarding adjusting alpha, it depends on what you will do with the evidence. If there will be a decision like “works - effective” or “does not work - not effective” based on the findings, you may need to adjust alpha. You also have the option to run more general tests instead of running many small tests. Perhaps check pwrss::pwrss.f.rmanova() function and see if it is relevant.