混合模型中的计划对比

Question

我的实验由两个因素组成，一是受试者之间的因素，一是受试者内部的因素。Age是受试者间因素，有两个水平：（低和高）。Time是受试者内因素，包含三个级别：1、4、5。dv是我的因变量，id是每个参与者的标识符。我附上了 6 位第一批参与者的数据。

我使用 R 进行了方差分析，对这两个因素都产生了显着的结果。我有两个计划对比：

1. 在时间 1 中，我想比较两个年龄组（低与高）。即对象之间的比较。

2. 在年龄较低的情况下，我想比较时间 1 和时间 5。即对象内比较。

当然，我可以执行 t 检验，但这似乎不合适，因为我可以将我的标准误差估计基于此处的更多单元格。我的问题是如何进行上述对比，以及适当的自由度是多少？

structure(list(id = structure(c(1L, 1L, 1L, 1L, 1L, 1L, 2L, 2L, 
2L, 2L, 2L, 2L, 3L, 3L, 3L, 3L, 3L, 3L, 4L, 4L, 4L, 4L, 4L, 4L, 
5L, 5L, 5L, 5L, 5L, 5L), .Label = c("1", "2", "3", "5", "6", 
"7", "8", "11", "12", "13", "15", "17", "18", "19", "20", "21", 
"22", "23", "24", "25", "27", "28", "29", "31", "32", "34", "35", 
"36", "37", "38", "39", "40", "42", "43", "44", "45", "46", "47", 
"48", "49", "52", "53", "54", "55", "56", "58", "59", "60", "62", 
"63", "64", "66", "67", "68", "69", "70", "71", "72", "73", "74", 
"75", "77", "79", "80", "81", "83", "84", "85", "86", "87", "88", 
"89", "90", "91", "92", "93", "94", "96", "97", "98", "99", "100", 
"101", "102", "103", "104", "105", "106", "107", "108", "109", 
"110", "111", "112", "113", "114", "115", "116", "117", "118", 
"119", "120", "121", "122", "123", "124", "125", "126", "127", 
"128", "129", "130", "132", "133", "134", "135", "136", "137", 
"138", "139", "140", "142", "143", "144", "145", "146", "147", 
"148", "149", "150", "151", "152", "153", "154", "156", "157", 
"158", "159", "160", "161", "162", "163", "165", "166", "167", 
"168", "169", "171", "172", "174", "175", "176", "177", "178", 
"179", "180", "181", "182", "183", "184", "185", "186", "187", 
"188", "189", "190", "191", "192", "193", "194", "195", "196", 
"200", "201", "202", "203", "204", "205", "206", "208", "209", 
"210", "212", "213", "214", "215", "216", "217", "218", "219", 
"220", "222", "223", "224", "226", "228", "230", "231", "232", 
"233", "234", "236", "237", "238", "239", "240", "241", "242", 
"243", "244", "246", "247", "248", "249", "250", "251", "252", 
"253", "254", "255", "256", "257", "258", "260", "261", "262", 
"263", "266", "267", "269", "270", "271", "272", "273", "274", 
"275", "276", "277", "278", "279", "280", "281", "282", "283", 
"284", "285", "286", "287", "288", "289", "290", "291", "292", 
"293", "294", "295", "296", "298", "299", "300"), class = "factor"), 
    age = structure(c(2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 
    2L, 2L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 2L, 
    2L, 2L, 2L, 2L, 2L), .Label = c("high", "low"), class = "factor"), 
    time = structure(c(3L, 3L, 2L, 2L, 1L, 1L, 3L, 3L, 2L, 2L, 
    1L, 1L, 3L, 3L, 2L, 2L, 1L, 1L, 3L, 3L, 2L, 2L, 1L, 1L, 3L, 
    3L, 2L, 2L, 1L, 1L), .Label = c("1", "4", "5"), class = "factor"), 
    dv = c(104, 102, 104, 103, 104, 104, 102, 102, 102, 102, 
    106, 106, 106, 106, 107, 107, 106, 106, 106, 107, 105, 104, 
    106, 107, 104, 101, 104, 101, 104, 106)), row.names = c(NA, 
-30L), class = c("tbl_df", "tbl", "data.frame"))

Answer 1

加载lme4适合多级模型（重复测量方差分析的替代方法）的包，并且emmeans可以进行对比。

\n\n

library(lme4)\nlibrary(emmeans)\n

\n\n

这适合一个模型，其中 是dv通过交互作用预测的（R 自动填充主效应），加上随机截距和time\xe2\x80\x94 的随机效应，两者都嵌套在 中id。dat是我从dput你的帖子中保存的内容。

\n\n

mod <- lmer(dv ~ age * time + (1 + time | id), dat)\n

\n\n

对比令人困惑，我总是担心自己会弄错。所以我们可以用来emmeans找到它们。我们可以拟合一个对象，它获取和emmeans的每个组合的值：timeage

\n\n

emm_mod <- emmeans(mod, ~ time + age)\n

\n\n

我们想要的对比是结果中的第三个和第十四个pairs()（您自己运行它看看它是什么样子）。coef()您可以通过放置物体周围来获得您想要的特定对比度pairs()。您只需要两列\xe2\x80\x943rd 和 14th：

\n\n

(contr_mat <- coef(pairs(emm_mod))[, c("c.3", "c.14")])\n

\n\n

返回：

\n\n

       c.3 c.14\n1,high   1    0\n4,high   0    0\n5,high   0    0\n1,low   -1    1\n4,low    0    0\n5,low    0   -1\n

\n\n

您可以通过在 中指定这两个对比来重点关注它们contr。您还可以使用您选择的 p 值调整来调整以adjust任何 \xe2\x80\x94I\ 开头"holm"：

\n\n

emmeans(mod, ~ time + age, contr = contr_mat, adjust = "holm")\n

\n\n

该contrasts位将为您提供您感兴趣的 p 值：

\n\n

$emmeans\n time age    emmean        SE df lower.CL upper.CL\n 1    high 106.2500 0.6518719  3 104.1755 108.3245\n 4    high 105.7500 0.8544406  3 103.0308 108.4692\n 5    high 106.2500 0.4759431  3 104.7353 107.7647\n 1    low  105.0000 0.5322511  3 103.3061 106.6939\n 4    low  102.6667 0.6976478  3 100.4464 104.8869\n 5    low  102.5000 0.3886059  3 101.2633 103.7367\n\nDegrees-of-freedom method: kenward-roger \nConfidence level used: 0.95 \n\n$contrasts\n contrast estimate        SE df t.ratio p.value\n c.3          1.25 0.8415630  3   1.485  0.2341\n c.14         2.50 0.7104068  3   3.519  0.0779\n\nP value adjustment: holm method for 2 tests \n

\n\n

您还可以尝试adjust = "none"：

\n\n

> emmeans(mod, ~ time + age, contr = contr_mat, adjust = "none")$contrasts\n contrast estimate        SE df t.ratio p.value\n c.3          1.25 0.8415630  3   1.485  0.2341\n c.14         2.50 0.7104068  3   3.519  0.0389\n

\n