lmerTest::lmer：手动创建交互项

Question

精简版

factorC <- paste(factorA, factorB)

factorA*factorB

相当于线性混合效应模型公式中的

factorA+factorB+factorC

吗？

背景

我使用包

ancombc2()

中的函数

ANCOMBC

。

ANCOMBC

是一个用于使用物种频率数据进行差异丰度测试的软件包。有关详细信息，请参阅此处：https://www.nature.com/articles/s41592-023-02092-7但是，我希望我们不必针对这个问题处理此包的特殊性，因为

ANCOMBC

正在使用

lmer()

/

lmerTest

中的

lme4

。

我的问题是关于“手动”包含交互术语。 ANCOMBC不支持交互术语。然而，作者写道：

问：ancombc 或 ancombc2 函数可以处理分析中的交互项吗？

A：不幸的是，在 ancombc 或 ancombc2 的 fix_formula 参数中包含交互项可能会导致多组比较的复杂性和潜在的混乱。为了解决这个问题，建议在公式之外手动创建感兴趣的交互项并进行相应的分析。

通过手动创建交互项，可以确保分析准确捕捉变量之间的交互作用。创建交互项后，您可以将其包含在 fix_formula 参数或分析的任何其他相关部分中，具体取决于您的具体研究问题和设计。

来自：
https://github.com/FrederickHuangLin/ANCOMBC

我的解释和可重现的例子

lmer()

我使用

lmer()

中的

lmeTest

进行测试。我将使用以下示例数据展示我如何解释作者的解释：

在一项实验中，对来自两个治疗组和两个年龄组的样本进行了观察。此外，该实验在两年内重复进行。因此，实验包括以下因素：

“年份”（随机），级别为“1”和“2”

“治疗”（固定），级别为“对照”和“治疗”
“年龄”（固定），级别为“旧”和“新”

# create example data set.seed(123) dat_year1 <- data.frame(year = factor(1), treatment = c(rep("control", 20), rep("treated", 20)), age = c(rep("old", 10), rep("new", 10), rep("old", 10), rep("new", 10)), value = c(rnorm(n = 10, mean = 10, sd = 2), rnorm(n = 10, mean = 10, sd = 2), rnorm(n = 10, mean = 15, sd = 2), rnorm(n = 10, mean = 16, sd = 2))) set.seed(321) dat_year2 <- data.frame(year = factor(2), treatment = c(rep("control", 20), rep("treated", 20)), age = c(rep("old", 10), rep("new", 10), rep("old", 10), rep("new", 10)), value = c(rnorm(n = 10, mean = 10, sd = 2), rnorm(n = 10, mean = 10, sd = 2), rnorm(n = 10, mean = 15, sd = 2), rnorm(n = 10, mean = 18, sd = 2))) dat <- rbind(dat_year1, dat_year2)

及其可视化：

我有兴趣测试“治疗”和“年龄”之间的相互作用。对于

lmer()

我会这样做：

# run lmer()
test <- lmerTest::lmer(value~treatment*age+(1|year), data = dat)
summary(test)

返回：

Linear mixed model fit by REML. t-tests use Satterthwaite's method [lmerModLmerTest] Formula: value ~ treatment * age + (1 | year) Data: dat REML criterion at convergence: 322.2 Scaled residuals: Min 1Q Median 3Q Max -2.41222 -0.69075 0.07493 0.44019 2.22631 Random effects: Groups Name Variance Std.Dev. year (Intercept) 0.05503 0.2346 Residual 3.44757 1.8568 Number of obs: 80, groups: year, 2 Fixed effects: Estimate Std. Error df t value (Intercept) 10.6492 0.4471 7.6713 23.819 treatmenttreated 6.6717 0.5872 75.0000 11.363 ageold -0.4259 0.5872 75.0000 -0.725 treatmenttreated:ageold -2.6641 0.8304 75.0000 -3.208 Pr(>|t|) (Intercept) 1.81e-08 *** treatmenttreated < 2e-16 *** ageold 0.47044 treatmenttreated:ageold 0.00196 ** --- Signif. codes: 0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1 Correlation of Fixed Effects: (Intr) trtmnt ageold tretmnttrtd -0.657 ageold -0.657 0.500 trtmnttrtd: 0.464 -0.707 -0.707

我将创建交互项的建议解释为应手动创建一个新因素，该因素结合了要测试交互的因素的级别：

dat$treatment_age <- paste(dat$treatment, "_", dat$age, sep = "")

现在我将其作为另一个固定因子包含在模型公式中：

test_manual <- lmerTest::lmer(value~treatment+age+treatment_age+(1|year), data = dat) summary(test_manual)

返回：

Linear mixed model fit by REML. t-tests use Satterthwaite's method [lmerModLmerTest] Formula: value ~ treatment + age + treatment_age + (1 | year) Data: dat REML criterion at convergence: 322.2 Scaled residuals: Min 1Q Median 3Q Max -2.41222 -0.69075 0.07493 0.44019 2.22631 Random effects: Groups Name Variance Std.Dev. year (Intercept) 0.05503 0.2346 Residual 3.44757 1.8568 Number of obs: 80, groups: year, 2 Fixed effects: Estimate Std. Error df (Intercept) 10.6492 0.4471 7.6713 treatmenttreated 6.6717 0.5872 75.0000 ageold -3.0900 0.5872 75.0000 treatment_agecontrol_old 2.6641 0.8304 75.0000 t value Pr(>|t|) (Intercept) 23.819 1.81e-08 *** treatmenttreated 11.363 < 2e-16 *** ageold -5.263 1.30e-06 *** treatment_agecontrol_old 3.208 0.00196 ** --- Signif. codes: 0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1 Correlation of Fixed Effects: (Intr) trtmnt ageold tretmnttrtd -0.657 ageold 0.000 -0.500 trtmnt_gcn_ -0.464 0.707 -0.707 fit warnings: fixed-effect model matrix is rank deficient so dropping 2 columns / coefficients

对于“治疗”，所有统计数据都是相同的，对于“年龄”，它们有很大不同，对于“治疗：年龄”，估计值和 t 值现在为负值。此外，还有一个拟合警告“固定效应模型矩阵存在秩不足，因此删除 2 列/系数”。

问题

如何改进以手动正确包含交互项？

与

ancombc2

一起使用

当我尝试将

ancombc2()

与手动包含的交互项一起使用时，它不会执行但返回错误：

library(ANCOMBC)
library(phyloseq)

ancombc2(data = phyloseq_object,
         tax_level = NULL,
         fix_formula = "treatment+age+treatment_age",
         verbose = TRUE)

返回：

Obtaining initial estimates ... Error: Estimation failed for the following covariates: treated_old, treated_new Please ensure that these covariates do not have missing values and check for multicollinearity before re-estimating the model

Answer 1

factorC

本身相当于

factorA*factorB

- 您不需要包含

factorA+factorB

。这就是为什么您在较长版本中收到排名不足警告的原因，

treatment_age

中的两个级别已被

treatment

和

age

主要效果覆盖。

lmerTest::lmer：手动创建交互项

问题描述投票：0回答：1

精简版

背景

我使用

如何改进以手动正确包含交互项？

一起使用
当我尝试将

1个回答

最新问题

lmerTest::lmer：手动创建交互项

问题描述 投票：0回答：1

精简版

背景

我使用

如何改进以手动正确包含交互项？

一起使用 当我尝试将

1个回答

最新问题

问题描述投票：0回答：1

一起使用
当我尝试将