Я провел некоторые поиски по этому поводу, но найденные мной сообщения в списке рассылки связаны с человеком, который не указывает случайный эффект в nlme, тогда как я это сделал.У меня также есть книга «Модели смешанных эффектов в S и S-Plus» Пинхейро и Бейтса, но я не могу решить мою проблему из этой книги.
Я все еще работаю над анализом данных о питательных веществах и имеютеперь перешел на реальные данные.Данные поступают из опроса населения и характеризуются планом повторных измерений, поскольку у каждого респондента есть два 24-часовых отзыва потребления питательных веществ.
Я успешно приспособил модель lme4 к своим данным, и теперь я пытаюсьчтобы выяснить, что произойдет, если я вместо этого использую нелинейный метод.Ниже приведен снимок моих данных:
head(Male.Data)
RespondentID Age SampleWeight IntakeDay IntakeAmt AgeFactor BoxCoxXY
2 100020 12 0.4952835 Day1Intake 12145.852 9to13 15.61196
7 100419 14 0.3632839 Day1Intake 9591.953 14to18 15.01444
8 100459 11 0.4952835 Day1Intake 7838.713 9to13 14.51458
12 101138 15 1.3258785 Day1Intake 11113.266 14to18 15.38541
14 101214 6 2.1198688 Day1Intake 7150.133 4to8 14.29022
18 101389 5 2.1198688 Day1Intake 5091.528 4to8 13.47928
И краткая информация о данных:
str(Male.Data)
'data.frame': 4498 obs. of 7 variables:
$ RespondentID: Factor w/ 4487 levels "100013","100020",..: 2 7 8 12 14 18 19 20 21 22 ...
$ Age : int 12 14 11 15 6 5 10 2 2 9 ...
$ SampleWeight: num 0.495 0.363 0.495 1.326 2.12 ...
$ IntakeDay : Factor w/ 2 levels "Day1Intake","Day2Intake": 1 1 1 1 1 1 1 1 1 1 ...
$ IntakeAmt : num 12146 9592 7839 11113 7150 ...
$ AgeFactor : Factor w/ 4 levels "1to3","4to8",..: 3 4 3 4 2 2 3 1 1 3 ...
$ BoxCoxXY : num 15.6 15 14.5 15.4 14.3 ...
Используя пакет lme4, я успешно установил линейные смешанные эффектыиспользование модели (случайный эффект от субъектов, а IntakeDay - это коэффициент повторного измерения, связанный с BoxCoxXY, который является преобразованием IntakeAmt):
Male.lme1 <- lmer(BoxCoxXY ~ AgeFactor + IntakeDay + (1|RespondentID),
data = Male.Data,
weights = SampleWeight)
Я пытался использоватьпакет nlme для проверки соответствия нелинейной модели для сравнения двух, но я не могу заставить свой синтаксис работать.Моя первоначальная проблема заключалась в том, что, похоже, не существует релевантной модели SelfStart для моих данных, поэтому я использовал geeglm для генерации начальных значений (коэффициенты, сохраненные во фрейме данных с именем Male.nlme.start).Но теперь я просто получаю сообщение об ошибке:
Error in getGroups.data.frame(dataMix, eval(parse(text = paste("~1", deparse(groups[[2]]), :
Invalid formula for groups
Я не могу понять, что я делаю неправильно, я использую синтаксис nlme:
Male.nlme1 <- nlme(BoxCoxXY ~ AgeFactor + IntakeDay + RespondentID , data = Male.Data,
fixed = AgeFactor + IntakeDay ~ 1,
random = RespondentID ~ 1,
start=c(Male.nlme.start[,"Estimate"]))
У меня естьпробовал анализ как с включением RespondentID, так и без него в общую спецификацию модели, и это, похоже, не оказало влияния.
Причина, по которой я пытаюсь продолжать использовать нелинейный метод, заключается в том, что исходный анализ в SASиспользовал нелинейный подход.Хотя мои остатки и т. Д. Выглядят приемлемо хорошо из анализа Ime, мне любопытно посмотреть, какое влияние окажет нелинейный подход.
Если это полезно, traceback() является результатом последней попытки анализа, которая включает RespondentID:
5: stop("Invalid formula for groups")
4: getGroups.data.frame(dataMix, eval(parse(text = paste("~1", deparse(groups[[2]]),
sep = "|"))))
3: getGroups(dataMix, eval(parse(text = paste("~1", deparse(groups[[2]]),
sep = "|"))))
2: nlme.formula(BoxCoxXY ~ AgeFactor + IntakeDay, data = Male.Data,
fixed = AgeFactor + IntakeDay ~ 1, random = RespondentID ~
1, start = c(Male.nlme.start[, "Estimate"]))
1: nlme(BoxCoxXY ~ AgeFactor + IntakeDay, data = Male.Data, fixed = AgeFactor +
IntakeDay ~ 1, random = RespondentID ~ 1, start = c(Male.nlme.start[,
"Estimate"]))
Может кто-нибудь подсказать, где я ошибся?Я начинаю задаваться вопросом, если (1) слишком много уровней факторов для RespondentID, чтобы работать в nlme, или (2) метод будет работать, только если я предоставлю параметр запуска для RespondentID, что кажется бессмысленнымс данными, которые у меня есть, так как это мой идентификатор субъекта.
Обновление: чтобы ответить Бену, модель SAS nlmixed - это общая логарифмическая функция вероятности для фиксированных эффектов:
ll1 <- log(1/sqrt(2*pi*Scale))
ll2 <- as.data.frame(-(BoxCoxXY - Intercept + AgeFactor + IntakeDay + u2)^2)/(2*Scale)+(Lambda.Value-1)*log(IntakeAmt)
ll <- ll1 + ll2
model IntakeAmt ~ general(ll)
где:
Scale =значение дисперсии от geeglm и
Lambda.Value = лямбда-значение, связанное с максимальным логарифмическим выходом правдоподобия от более раннего boxcox(), который использовался для преобразования IntakeAmt в BoxCoxXY по формуле Male.Data$BoxCoxXY <- (Male.Data$IntakeAmt^Lambda.Value-1)/Lambda.Value
Оператор random в коде SAS:
random u1 u2 ~ normal([0,0][&vu1,COV_U1U2,&vu2]) subject=RespondentID
, поэтому в модели есть два условия ошибки, и они оба подходят как случайные эффекты.Вторая квадратная скобка представляет нижний треугольник матрицы дисперсии случайных эффектов, перечисленных в порядке строк, и указывается с использованием макропеременных SAS в синтаксисе SAS.
Краткое изложение модели, которую я получил, представляет собойобычный однострочный обзор, который показывает матрицу ковариат (BX) плюс компонент ошибки, так что здесь нет большой помощи.
Второе обновление: я понял, что не удалил уровни RespondentID, связанные с женщинойсубъекты, как я факторизовал RespondentID по всему фрейму данных, прежде чем разделить его на отдельные фреймы по полу для анализа.Я повторил анализ nlme после удаления неиспользуемых уровней факторов для RespondentID, и я получил ту же ошибку.Результаты lmer одинаковы - это хорошо знать.:)