提建议的解决方法

• 缺乏生物学意义
• 需要后续研究
• 机制不明确，等等

increase sample size:但是必须给出统计学上的理由

• Power estimation

1. Chi-sq
   X是我要计算的数据集，可以是向量也可以是matrix

model <- chisq.test(x)
chisq_test <- model$statistic
N <- sum(X)
R <- row_numbers
C <- col_numbers
V <- sqrt(chisq_stat / (N * min(R-1, C-1)))  # 计算Cramer's V


current <-pwr.chisq.test(w=V,df=(R-1)*(C-1),sig.level = 0.05,N = N)
current
expected <- pwr.chisq.test(w=V, df=(R-1)*(C-1), sig.level=0.05, power=0.8)
expected

情况一：满足了power
Based on the power analysis for the chi-square test, at a significance level of 0.05, the current power is current$power. This power value has reached the commonly used standard of 0.8 in statistical analysis, indicating that with the current sample size, we have sufficient ability to detect an actual effect if it exists.

Since the power has reached the desired level, there is no need to increase the sample size. We can proceed with the actual chi-square test and subsequent statistical analyses based on the current sample size.

情况二：不满足power
Based on the power analysis for the chi-square test, at a significance level of 0.05, the current power is current$power, which is lower than the commonly used standard of 0.8 in statistical analysis. This indicates that with the current sample size, our ability to detect an actual effect, if it exists, is relatively low.

To increase the power to the desired level of 0.8, it is recommended to increase the sample size. Based on the expected power of 0.8, the minimum required sample size is expected\(n. Therefore, I recommend increasing the sample size to at least expected\)n to ensure sufficient statistical power to detect an actual effect.

2. Anova （只能做one-way的，两变量不考虑这个）

library(pwr)

# 计算eta squared效应量
x <- c(10, 12, 14, 16, 18, 20, 22, 24, 8, 10, 12, 14, 16, 18, 20, 22, 6, 8, 10, 12, 14, 16, 18, 20)  # 输入所有数据
group <- c(rep("A", 8), rep("B", 8), rep("C", 8))  #


anova_model <- aov(x ~ group)
anova_summary <- anova(anova_model)
ss_total <- sum(anova_summary$"Sum Sq")
ss_group <- anova_summary$"Sum Sq"[1]
ss_residual <- anova_summary$"Sum Sq"[2]
eta_squared <- ss_group / ss_total


# 计算所需样本量
k <- length(unique(group))  # 组数
current <- pwr.anova.test(f = eta_squared / (1 - eta_squared),
                          k = k, n= sum(x)/k,
                          sig.level = 0.05)
expected <- pwr.anova.test(f = eta_squared / (1 - eta_squared),
                           sig.level = 0.05, k = k,
                           power = 0.8)
current
expected

3. t test

library(pwr)

# 计算Cohen's d效应量
x <- c(12,212,3232,323,31,31,43)  # 输入第一组数据
y <- c(12,33,22,4,42,111,3,222)  # 输入第二组数据
n1 <- length(x)
n2 <- length(y)
d <- (mean(x) - mean(y)) / sqrt(((n1-1)*var(x) + (n2-1)*var(y)) / (n1 + n2 - 2))

# 计算所需样本量
current <- pwr.t.test(n = n1, d = d,
                      sig.level = 0.05,
                      type = "two.sample")
expected <- pwr.t.test(d = d,
                       sig.level = 0.05,
                       power = 0.8,
                       type = "two.sample")

cat("Current power is:", current$power, "\n")
cat("With an expected power of 0.8, the minimum sample size required is:", expected$n, "\n")

Based on your request, here is a possible report for the power analysis of a two-sample t-test, discussing two scenarios:

Scenario 1: Power meets the desired level

Report:
The power analysis for the two-sample t-test was conducted to assess the ability to detect a significant difference in means between two independent groups. The effect size, Cohen's d, was calculated based on the provided data.

At a significance level of 0.05, the current power is current$power. This power value meets the commonly accepted standard of 0.8 in statistical analysis, indicating that with the current sample size, we have sufficient ability to detect an actual effect if it exists.

Recommendation:
Since the power has reached the desired level of 0.8, there is no need to increase the sample size. We can proceed with the actual two-sample t-test and subsequent statistical analyses based on the current sample size of n1 = [current sample size for group 1] and n2 = [current sample size for group 2].

Scenario 2: Power does not meet the desired level

Report:
The power analysis for the two-sample t-test was conducted to assess the ability to detect a significant difference in means between two independent groups. The effect size, Cohen's d, was calculated based on the provided data.

At a significance level of 0.05, the current power is current$power, which is lower than the commonly accepted standard of 0.8 in statistical analysis. This indicates that with the current sample size, our ability to detect an actual effect, if it exists, is relatively low.

Recommendation:
To increase the power to the desired level of 0.8, it is recommended to increase the sample size. Based on the expected power of 0.8, the minimum required sample size for each group is expected\(n. Therefore, I recommend increasing the sample size to at least n1 = expected\)n and n2 = expected$n to ensure sufficient statistical power to detect an actual effect.

It is important to note that the power analysis results depend on the estimation of the effect size, Cohen's d. If the actual effect size deviates from the estimated value, further adjustment of the sample size may be necessary.

Please note that in this report, you should replace "current\(power" and "expected\)n" with the actual values obtained from your power analysis. Additionally, you should replace "[current sample size for group 1]" and "[current sample size for group 2]" with the respective sample sizes for each group in the current data.

• 旧版本，不删除了怕还是有点用

library(pwr)

# t检验
power_t_test <- function(m1, m2, sd1, sd2, n1, n2, alpha = 0.05) {
  sd_pooled <- sqrt(((n1 - 1) * sd1^2 + (n2 - 1) * sd2^2) / (n1 + n2 - 2))
  d <- (m1 - m2) / sd_pooled
  power <- pwr.t.test(n = n1 + n2, d = d, sig.level = alpha, type = "two.sample", alternative = "two.sided")$power
  return(list(d = d, power = power))
}

# Wilcoxon秩和检验
power_wilcox_test <- function(z, n, alpha = 0.05) {
  r <- z / sqrt(n)
  power <- pwr.r.test(n = n, r = r, sig.level = alpha)$power
  return(list(r = r, power = power))
}

# 单因素方差分析
power_anova_test <- function(ss_between, ss_total, n, k, alpha = 0.05) {
  eta_sq <- ss_between / ss_total
  f <- sqrt(eta_sq / (1 - eta_sq))
  power <- pwr.anova.test(k = k, n = n, f = f, sig.level = alpha)$power
  return(list(eta_sq = eta_sq, f = f, power = power))
}

# Fisher's精确检验和卡方检验
power_fisher_chisq_test <- function(chi_sq, n, k, alpha = 0.05) {
  if (k == 2) {
    effect_size <- sqrt(chi_sq / n)
  } else {
    effect_size <- sqrt(chi_sq / (n * (k - 1)))
  }
  power <- pwr.chisq.test(w = effect_size, N = n, df = (k - 1) * (k - 1), sig.level = alpha)$power
  return(list(effect_size = effect_size, power = power))
}

# 示例用法
t_test_result <- power_t_test(m1 = 10, m2 = 12, sd1 = 2, sd2 = 2.5, n1 = 30, n2 = 30)
print(t_test_result)

wilcox_test_result <- power_wilcox_test(z = 2.5, n = 50)
print(wilcox_test_result)

anova_test_result <- power_anova_test(ss_between = 100, ss_total = 500, n = 10, k = 3)
print(anova_test_result)

fisher_chisq_test_result <- power_fisher_chisq_test(chi_sq = 8, n = 100, k = 2)
print(fisher_chisq_test_result)

好的,我来解释一下这些函数中使用的统计量:

1. power_t_test函数:

   • d: 这是Cohen's d效应量,用于衡量两组平均值之间的差异大小。它是两组平均值的差值除以pooled标准差。
   • power: 这是在给定效应量和样本量的情况下,t检验能够检测到真实效应的概率。

2. power_wilcox_test函数:

   • r: 这是Wilcoxon秩和检验的效应量,计算方式为标准化的Wilcoxon统计量z除以sqrt(n)。
   • power: 这是在给定效应量和样本量的情况下,Wilcoxon秩和检验能够检测到真实效应的概率。

3. power_anova_test函数:

   • eta_sq: 这是η^2效应量,用于衡量在单因素方差分析中,因子对总变异的解释程度。它等于因子平方和除以总平方和。
   • f: 这是f效应量,与η^2相关,计算方式为sqrt(eta_sq / (1 - eta_sq))。
   • power: 这是在给定效应量和样本量的情况下,单因素方差分析能够检测到真实效应的概率。

4. power_fisher_chisq_test函数:

   • effect_size: 这是效应量,对于2x2列联表,计算方式为sqrt(chi_sq / n);对于其他情况,计算方式为sqrt(chi_sq / (n * (k - 1)))。
   • power: 这是在给定效应量和样本量的情况下,Fisher's精确检验或卡方检验能够检测到真实效应的概率。

总的来说,这些函数中使用的统计量包括Cohen's d、Wilcoxon秩和检验的r、η^2、f和卡方检验的效应量,它们都是衡量效应大小的指标。通过这些指标和样本量,可以计算出相应检验的统计力。

怎么阐述这些理由：
The current study has not enough power - it can detect differences between groups with the probability of xx. That may cause high-level of type-II error. If we want to achieve at least 0.8 power, we need to have that many animals in each group, and the groups must be balanced.

然后下面用power = 0.8作为阈值，看看需要的n大小是多少

pwr.anova.test(k = k, f = f, sig.level = alpha, power=0.8)