求适用于双向(2×2)累积链接模型有序ANOVA的效应量方法

双因素设计的秩相关效应量替代方案 #

你提到的rank_eta_squared()和rank_epsilon_squared()确实仅支持单因素设计。针对双因素有序累积链接模型，目前没有直接对应的内置秩效应量函数，但可以通过两种思路手动推导类似指标：

• 基于偏差变化的效应量：对比包含目标效应的模型与不含该效应的嵌套模型的偏差差值，结合总偏差或基线偏差计算，逻辑类似R²，能反映该效应对模型解释力的贡献。
• 伪R²拆分法：利用有序模型的伪R²指标（如McFadden's R²、Nagelkerke's R²），通过嵌套模型对比，拆分出主效应和交互效应各自的解释力占比。

关于Cramers V计算的可行性 #

你用sqrt(chi²/(n*df))计算Cramers V的方法，在**2×2设计（每个效应的df=1）**下是完全可行的：

• 当df=1时，Cramers V的公式会简化为sqrt(χ²/n)，和你写的sqrt(anova(mod)[2]/nrow(dat))逻辑一致（anova输出的第二列是卡方统计量）。
• 需要注意：这个指标基于Wald卡方，仅衡量分类变量间的关联强度，没有利用因变量的有序属性，属于粗略的效应量参考，但在2×2设计下可以正常使用。

优化后的示例代码 #

library(ordinal)
set.seed(123) # 添加随机种子保证结果可复现
n <- 100
dat <- data.frame(
  IV1 = sample(LETTERS[1:2], n, replace=TRUE),
  IV2 = sample(LETTERS[3:4], n, replace=TRUE),
  DV = factor(sample(1:6, n, replace=TRUE), ordered = TRUE)
)

# 拟合不同嵌套模型用于效应量计算
mod_full <- clm(DV ~ IV1 + IV2 + IV1:IV2, data = dat)
mod_no_interaction <- clm(DV ~ IV1 + IV2, data = dat)
mod_no_iv2 <- clm(DV ~ IV1, data = dat)

# 输出Wald卡方检验结果
anova(mod_full)

# 计算各效应的Cramers V
cramers_v_iv1 <- sqrt(anova(mod_full)[2, "Chisq"] / n)
cramers_v_iv2 <- sqrt(anova(mod_full)[3, "Chisq"] / n)
cramers_v_interaction <- sqrt(anova(mod_full)[4, "Chisq"] / n)

cat("IV1的Cramers V:", cramers_v_iv1, "\n")
cat("IV2的Cramers V:", cramers_v_iv2, "\n")
cat("交互效应的Cramers V:", cramers_v_interaction, "\n")

# 用偏差变化计算伪R²拆分效应贡献
# IV2的解释力贡献：无IV2模型 vs 全主效应模型
dev_iv2 <- deviance(mod_no_iv2) - deviance(mod_no_interaction)
pseudo_r2_iv2 <- dev_iv2 / deviance(mod_no_iv2)
# 交互效应的解释力贡献：全主效应模型 vs 全模型
dev_interaction <- deviance(mod_no_interaction) - deviance(mod_full)
pseudo_r2_interaction <- dev_interaction / deviance(mod_no_interaction)

cat("IV2的伪R²贡献:", pseudo_r2_iv2, "\n")
cat("交互效应的伪R²贡献:", pseudo_r2_interaction, "\n")

内容的提问来源于stack exchange，提问作者Adrian