无法通过 R 中的 emmeans 对比截距

Posted 2023-02-18

【中文标题】无法通过 R 中的 emmeans 对比截距

【英文标题】：Failed to contrast intercepts through emmeans in R

【发布时间】：2021-07-23 20:29:08

【问题描述】：

我想测试观察者对具有不同极性（正极（白色）和负极（黑色））的对比刺激的反应中的对称性。我将反应时间 (RT) 作为因变量，以及四个不同的对比。众所周知，响应时间遵循 Pieron 曲线，其渐近线位于 (1) 观察者阈值 (Inf) 和 (2) 位于介于 250 和 450 毫秒之间的基础 RT 处。 知识允许我们将自变量（有效对比度 EC）转换为 1/EC^2 (tEC) 的关系线性化，因此将 RT 与 EC 联系起来的方程变为：

RT = m * tEC + RT0

为了测试对称性，我建立了标准：两个极性的相同斜率和相同截距意味着对称。 为了获得系数，我制作了一个具有交互作用的线性模型（通过极性的虚拟变量进行编码：正或负）。 lm 的输出对我来说很清楚，但一些同事更喜欢更类似于 ANOVA 输出的东西。所以我决定用emmeans来做对比。斜率是可以的，但是在计算截距时会出现问题。 lm计算的截距和emmeans的输出有很大的不同，得出的结论也不同。在下文中，我重现了这个例子。 问题有两个：可以使用 emmeans 来解决我的问题吗？如果没有，是否可以通过其他包（哪一个）进行对比？

数据

RT1000	EC	tEC	Polarity
596.3564	-25	0.001600	Negative
648.2471	-20	0.002500	Negative
770.7602	-17	0.003460	Negative
831.2971	-15	0.004444	Negative
1311.3331	15	0.004444	Positive
1173.8942	17	0.003460	Positive
1113.7240	20	0.002500	Positive
869.3635	25	0.001600	Positive

代码

# Model
model <- lm(RT1000 ~ tEC * Polarity, data = Data)

# emmeans
library(emmeans)
# Slopes
m.slopes <- lstrends(model, "Polarity", var="tEC")
# Intercepts
m.intercept <- lsmeans(model, "Polarity")

# Contrasts
pairs(m.slopes)
pairs(m.intercept)

输出

模型

term	estimate	std.error	statistic	p.value
(Intercept)	449.948	66.829	6.733	0.003
tEC	87205.179	20992.976	4.154	0.014
PolarityPositive	230.946	94.511	2.444	0.071
tEC:PolarityPositive	58133.172	29688.551	1.958	0.122

斜坡（可以）

Polarity	tEC.trend	SE	df	lower.CL	upper.CL
Negative	87205.18	20992.98	4	28919.33	145491.0
Positive	145338.35	20992.98	4	87052.51	203624.2

contrast	estimate	SE	df	t.ratio	p.value
Negative - Positive	-58133.17	29688.55	4	-1.958101	0.12182

拦截（问题）

Polarity	lsmean	SE	df	lower.CL	upper.CL
Negative	711.6652	22.2867	4	649.7874	773.543
Positive	1117.0787	22.2867	4	1055.2009	1178.957

contrast	estimate	SE	df	t.ratio	p.value
Negative - Positive	-405.4135	31.51816	4	-12.86285	0.000211

通过 emmeans 计算的截距与通过 lm 计算的截距不同。我认为问题在于没有为 EC = 0 定义模型。但我不确定。

【参考方案1】：

你所说的拦截不是;它们是平均值为tEC 的模型预测。如果您想要拦截，请改用：

m.intercept <- lsmeans(model, "Polarity", at = list(tEC = 0))

您可以通过

了解正在使用的参考水平

ref_grid(model)     # or   str(m.intercept)

请注意这里拟合的模型由两条斜率不同的线组成；因此预测之间的差异取决于tEC 的值。因此，我强烈建议反对测试截距的比较；这些是tEC 值的预测，正如您所说，甚至不可能发生。相反，尽量少做一个数学家，在tEC的几个代表值上进行比较，例如，

LSMs <- lsmeans(model, "Polarity", at = list(tEC = c(0.001, 0.003, 0.005)))
pairs(LSMs, by = tEC)

您还可以轻松地可视化拟合线：

emmip(model, Polarity ~ tEC, cov.reduce = range)

【讨论】：

谢谢拉斯。在提交问题之前，我阅读了很多您的帖子。真的信息丰富....但我没有设法解决问题。在本例中，0 处的截距具有精确的含义，因为它是基础 RT。正如您所建议的那样，其他点可能会更有问题。

如果我的回答有帮助，我将不胜感激。