Я пытаюсь перейти от использования ez
пакета lme
для повторных измерений ANOVA (как я надеюсь, я смогу использовать пользовательские контрасты с lme
).
Следуя совету из этого поста в блоге, я смог настроить одну и ту же модель, используя обе функции aov
(как это делается ez
при запросе) и lme
. Однако, в то время как в примере , приведенном в этой должности в F -значения делать совершенно согласен между aov
и lme
(я проверил это, и они делают), это не относится к моим данным. Хотя значения F похожи, они не одинаковы.
aov
возвращает значение f 1,3399, lme
возвращает 1,3664. Я готов принять aov
результат как «правильный», поскольку это также то, что возвращает SPSS (и это то, что имеет значение для моего поля / руководителя).
Вопросов:
Было бы здорово, если бы кто-то мог объяснить, почему существует такая разница и как я могу использовать ее
lme
для получения достоверных результатов. (Я также хотел бы использоватьlmer
вместоlme
этого типа материала, если он дает «правильный» результат. Однако я до сих пор не использовал его.)После решения этой проблемы я хотел бы провести анализ контраста. Особенно меня будет интересовать контраст объединения первых двух уровней фактора (т. Е.
c("MP", "MT")
) И сравнение его с третьим уровнем фактора (т. Е."AC"
). Кроме того, тестирование третьего по сравнению с четвертым уровнем фактора (то есть, по"AC"
сравнению с"DA"
).
Данные:
tau.base <- structure(list(id = structure(c(1L, 2L, 3L, 4L, 5L, 6L, 7L, 8L,
9L, 10L, 11L, 12L, 13L, 14L, 15L, 16L, 17L, 18L, 19L, 20L, 21L,
22L, 1L, 2L, 3L, 4L, 5L, 6L, 7L, 8L, 9L, 10L, 11L, 12L, 13L,
14L, 15L, 16L, 17L, 18L, 19L, 20L, 21L, 22L, 1L, 2L, 3L, 4L,
5L, 6L, 7L, 8L, 9L, 10L, 11L, 12L, 13L, 14L, 15L, 16L, 17L, 18L,
19L, 20L, 21L, 22L, 1L, 2L, 3L, 4L, 5L, 6L, 7L, 8L, 9L, 10L,
11L, 12L, 13L, 14L, 15L, 16L, 17L, 18L, 19L, 20L, 21L, 22L), .Label = c("A18K",
"D21C", "F25E", "G25D", "H05M", "H07A", "H08H", "H25C", "H28E",
"H30D", "J10G", "J22J", "K20U", "M09M", "P20E", "P26G", "P28G",
"R03C", "U21S", "W08A", "W15V", "W18R"), class = "factor"), factor = structure(c(1L,
1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L,
1L, 1L, 1L, 1L, 1L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L,
2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 3L, 3L, 3L, 3L, 3L,
3L, 3L, 3L, 3L, 3L, 3L, 3L, 3L, 3L, 3L, 3L, 3L, 3L, 3L, 3L, 3L,
3L, 4L, 4L, 4L, 4L, 4L, 4L, 4L, 4L, 4L, 4L, 4L, 4L, 4L, 4L, 4L,
4L, 4L, 4L, 4L, 4L, 4L, 4L), .Label = c("MP", "MT", "AC", "DA"
), class = "factor"), value = c(0.9648092876, 0.2128662077, 1,
0.0607615485, 0.9912814024, 3.22e-08, 0.8073856412, 0.1465590332,
0.9981672618, 1, 1, 1, 0.9794401938, 0.6102546108, 0.428651501,
1, 0.1710644881, 1, 0.7639763913, 1, 0.5298989196, 1, 1, 0.7162733447,
0.7871177434, 1, 1, 1, 0.8560509327, 0.3096989662, 1, 8.51e-08,
0.3278862311, 0.0953598576, 1, 1.38e-08, 1.07e-08, 0.545290432,
0.1305621416, 2.61e-08, 1, 0.9834051136, 0.8044114935, 0.7938839461,
0.9910112678, 2.58e-08, 0.5762677121, 0.4750002288, 1e-08, 0.8584252623,
1, 1, 0.6020385797, 8.51e-08, 0.7964935271, 0.2238374288, 0.263377904,
1, 1.07e-08, 0.3160751898, 5.8e-08, 0.3460325565, 0.6842217296,
1.01e-08, 0.9438301877, 0.5578367224, 2.18e-08, 1, 0.9161424562,
0.2924856039, 1e-08, 0.8672987992, 0.9266688748, 0.8356425464,
0.9988463913, 0.2960361777, 0.0285680426, 0.0969063841, 0.6947998266,
0.0138254805, 1, 0.3494775301, 1, 2.61e-08, 1.52e-08, 0.5393467752,
1, 0.9069223275)), .Names = c("id", "factor", "value"), class = "data.frame", row.names = c(1L,
6L, 10L, 13L, 16L, 17L, 18L, 22L, 23L, 24L, 27L, 29L, 31L, 33L,
42L, 43L, 44L, 45L, 54L, 56L, 58L, 61L, 64L, 69L, 73L, 76L, 79L,
80L, 81L, 85L, 86L, 87L, 90L, 92L, 94L, 96L, 105L, 106L, 107L,
108L, 117L, 119L, 121L, 124L, 127L, 132L, 136L, 139L, 142L, 143L,
144L, 148L, 149L, 150L, 153L, 155L, 157L, 159L, 168L, 169L, 170L,
171L, 180L, 182L, 184L, 187L, 190L, 195L, 199L, 202L, 205L, 206L,
207L, 211L, 212L, 213L, 216L, 218L, 220L, 222L, 231L, 232L, 233L,
234L, 243L, 245L, 247L, 250L))
И код:
require(nlme)
summary(aov(value ~ factor+Error(id/factor), data = tau.base))
anova(lme(value ~ factor, data = tau.base, random = ~1|id))
lme
результатах из стандартного учебника ANOVA (предоставленногоaov
и который мне нужен), это не вариант для меня. В своей статье я хочу сообщить о ANOVA, а не о чем-то вроде ANOVA. Интересно, что Venables & Ripley (2002, стр. 285) показывают, что оба подхода приводят к одинаковым оценкам. Но различия в значениях F оставляют у меня плохое предчувствие. Кроме того,Anova()
(fromcar
) возвращает только значения Chi² дляlme
объектов. Поэтому для меня мой первый вопрос еще не ответил.lme
; но для контрастов,glht
работает наlm
подгонки тоже, не толькоlme
подходит. (Кроме того,lme
результаты являются стандартными результатами учебника тоже.)lm
для анализа повторных измерений. Толькоaov
может иметь дело с повторными мерами, но вернет объект класса,aovlist
который, к сожалению, не обрабатываетсяglht
.lm
использует остаточную ошибку в качестве члена ошибки для всех эффектов; когда есть эффекты, которые должны использовать другой термин ошибки,aov
необходимо (или вместо этого, используя результатыlm
для вычисления F-статистики вручную). В вашем примере термин ошибки дляfactor
- этоid:factor
взаимодействие, которое является термином остаточной ошибки в аддитивной модели. Сравните свои результаты сanova(lm(value~factor+id))
.Ответы:
Они отличаются, потому что модель Ime заставляет компонент дисперсии
id
быть больше нуля. Глядя на необработанную таблицу anova для всех терминов, мы видим, что среднеквадратичная ошибка для id меньше, чем для остаточных значений.Когда мы вычисляем компоненты дисперсии, это означает, что дисперсия из-за id будет отрицательной. Моя память об ожидаемых средних квадратах шаткая, но вычисление что-то вроде
Это звучит странно, но может случиться. Это означает, что средние значения для каждого идентификатора ближе друг к другу, чем вы ожидаете, учитывая количество остаточных изменений в модели.
Напротив, использование Ime заставляет эту дисперсию быть больше нуля.
Это сообщает о стандартных отклонениях, возводя в квадрат, чтобы получить выход дисперсии
9.553e-10
для идентификатора дисперсии и0.1586164
для остаточной дисперсии.Теперь вы должны знать, что использование
aov
для повторных измерений уместно, только если вы считаете, что корреляция между всеми парами повторных измерений идентична; это называется сложной симметрией. (Технически, сферичности требуется , но этого достаточно для теперь.) Одной из причин для использованияlme
болееaov
, что он может обрабатывать различные виды корреляционных структур.В этом конкретном наборе данных оценка для этой корреляции является отрицательной; это помогает объяснить, как среднеквадратичная ошибка для идентификатора была меньше, чем остаточная квадратичная ошибка. Отрицательная корреляция означает, что если бы первое измерение индивидуума было ниже среднего, в среднем его второе было бы выше среднего, что делало бы средние значения для индивидов менее изменчивыми, чем мы ожидали, если бы была нулевая корреляция или положительная корреляция.
Использование
lme
со случайным эффектом эквивалентно подгонке модели составной симметрии, где эта корреляция должна быть неотрицательной; мы можем подобрать модель, в которой корреляция может быть отрицательной, используяgls
:Этот стол ANOVA согласуется со столом по
aov
форме иlm
форме.Ок, и что? Что ж, если вы считаете, что отклонение
id
и корреляция между наблюдениями должны быть неотрицательными, тоlme
подгонка на самом деле является более подходящей, чем использование подгонки,aov
илиlm
ее оценка остаточной дисперсии немного лучше. Тем не менее, если вы считаете , что корреляция между наблюдениями может быть отрицательным,aov
илиlm
илиgls
лучше.Вы также можете быть заинтересованы в дальнейшем изучении корреляционной структуры; чтобы посмотреть на общую структуру корреляции, вы бы сделали что-то вроде
Здесь я ограничиваю вывод только структурой корреляции. Значения от 1 до 4 представляют четыре уровня
factor
; мы видим, что фактор 1 и фактор 4 имеют довольно сильную отрицательную корреляцию:Одним из способов выбора между этими моделями является тест отношения правдоподобия; это показывает, что модель случайных эффектов и модель общей структуры корреляции статистически значимо не различаются; когда это происходит, более простая модель обычно предпочтительнее.
источник
lme
чтобы получить те же результаты, что и сaov
(и, следовательно, включитьlme
все ANOVA), а именно, использовать аргумент корреляции в вызовеlme
:anova(lme(value ~ factor, data = tau.base, random = ~1|id, correlation = corCompSymm(form = ~1|id)))
aov
,lme
без соединения симметрии, иlme
с соединением симметрии) имеют точно такое же количество ДФСА.anova.lme()
. Из вашего ответа я понял, что связь между ANOVA и смешанными моделями заключается в их корреляционной структуре. Затем я прочитал, что введение симметричной структуры корреляции компоновки приводит к равенству между двумя подходами. Поэтому я навязал это. Я понятия не имею, если это съедает другой дф. Вывод однако не согласуется с этой интерпретацией.aov()
соответствует модели сlm()
использованием наименьших квадратов,lme
соответствует максимальной вероятности. Эта разница в том, как оцениваются параметры линейной модели, вероятно, объясняет (очень маленькую) разницу в ваших f-значениях.На практике (например, для проверки гипотез) эти оценки одинаковы, поэтому я не понимаю, как одно можно считать «более достоверным», чем другое. Они приходят из разных модельных парадигм.
Для контрастов вам нужно установить контрастную матрицу для ваших факторов. Venebles и Ripley показывают, как это сделать, на стр. 143, стр. 146 и стр. 293-294 4-го издания.
источник
lme
илиlmer
для расчета ANOVA (строго говоря), поскольку он использует метод, который похож, но не идентичен. Так нет ли способа вычисления контрастов для ANOVA с повторным измерением в R?