Я хотел бы знать, как правильно обратно преобразовать выходные данные одномерной линейной модели смешанных эффектов, чтобы интерпретировать ее. Я не публиковал данные, чтобы ответить на мой вопрос, потому что на мой вопрос можно было бы ответить без данных.
Моя модель (упрощенная для целей этого вопроса):
library(lme4)
m1<-lmer(activity ~ sex + BirthDate+ (1|id), data=merge.data)
> m1
Linear mixed model fit by REML ['lmerMod']
Formula: activity ~ sex + BirthDate + (1 | id)
Data: merge.data
REML criterion at convergence: 572.0483
Random effects:
Groups Name Std.Dev.
id (Intercept) 0.7194
Residual 1.4651
Number of obs: 150, groups: id, 89
Fixed Effects:
(Intercept) sexM BirthDate
-0.08661 0.20718 0.43022
Где:
activity
— непрерывная переменная откликаsex
— категориальная переменная с 2 уровнями (женский и мужской)BirthDate
— непрерывная переменная; BirthDate
— это количество дней, прошедших с 1 января, а затем оно центрировано по среднему значению и стандартизировано до одного стандартного отклонения.id
— это случайный эффект индивидуальной идентичностиmerge.data
— это имя моего набора данныхДо того, как BirthDate
будет среднецентрированным и стандартизированным до одного стандартного отклонения:
> summary(merge.data$BirthDate)
Min. 1st Qu. Median Mean 3rd Qu. Max.
94.96 115.96 121.96 122.67 127.96 138.96
После того, как BirthDate
центрирован по среднему значению и стандартизирован до одного стандартного отклонения:
merge.data<-merge.data %>%
mutate(BirthDate = ((BirthDate-mean(BirthDate))/(1*(sd(BirthDate)))))
> summary(merge.data$BirthDate)
Min. 1st Qu. Median Mean 3rd Qu. Max.
-3.09082 -0.74816 -0.07883 0.00000 0.59050 1.81761
Я хотел бы знать, каково среднее значение для sex
и BirthDate
. Основываясь на чтении Книга R Кроули, я могу получить среднее значение из моей модели m1
с помощью следующего кода:
tapply(predict(m1,type = "response"), merge.data$sex,mean) #gives you the back-transformed mean for sex from the model "m1"
F M
-0.08334649 0.11199685
Это говорит о том, что средний показатель активности для женщин составляет -0,083, а для мужчин - 0,11.
Когда я пробую это для BirthDate
, вот так:
tapply(predict(m1,type = "response"), merge.data$BirthDate,mean)
-3.09082367412411 -1.6406056364576 -1.52905040279094 #mean centered birth date
-0.79030344 -0.87012920 -0.44792213 #activity score
and so on...
В итоге я получаю 1 среднее значение для каждой даты рождения (BirthDate
среднее значение центрировано и стандартизировано до одного стандартного отклонения). В отличие от sex
, я ничего не могу сделать с этой информацией... Я пытаюсь представить эффект (размер эффекта) увеличения даты рождения на активность.
Что я хотел бы в конечном итоге сделать, так это сказать, что на каждый 1 день увеличения даты рождения происходит [число из модели] увеличение показателя активности.
BirthDate на самом деле не дата, а дни из какой-то эпохи?
@ R5W Ничего не произойдет, если я наберу ml
. @ Р.С. Я уточню, что BirthDate
в вопросе.
Я думаю, @G5W означает ввести название модели m1
с «1», а не с «l».
@ R5W Я тоже добавил этот вывод.
Это нормально, если вы называете меня Р, но мои друзья зовут меня Г.
@ G5W Извините за это!
Когда вы распечатываете модель, набрав m1
, эта часть:
Fixed Effects:
(Intercept) sexM BirthDate
-0.08661 0.20718 0.43022
сообщает вам наклоны, то есть насколько результат изменится в зависимости от изменения входных данных. В частности, если вы увеличите BirthDate на единицу (и оставите все остальное без изменений), прогнозируемый показатель активности увеличится на 0,43022.
Вы не предоставляете никаких данных, поэтому я не могу работать напрямую с вашими данными и вашей моделью. Вместо этого я проиллюстрирую некоторые данные, встроенные в R, данные радужной оболочки.
## Build a linear model
Mod1 = lm(Petal.Length ~ ., data=iris[,1:4])
Теперь мы могли бы просто набрать Mod1
, но это дает больше, чем я хочу видеть. Мы можем ограничить наше внимание интересной частью, используя
Mod1$coefficients
(Intercept) Sepal.Length Sepal.Width Petal.Width
-0.2627112 0.7291384 -0.6460124 1.4467934
Это дает наклон для каждой из переменных-предикторов (и точки пересечения).
Я хочу проиллюстрировать, как ответ Petal.Length
зависит от входных данных.
Я просто возьму какую-то точку, изменю один предиктор и посмотрю на результат.
NewPoint = iris[30,1:4]
NewPoint[,1] = NewPoint[,1]+1
iris[30, 1:4]
Sepal.Length Sepal.Width Petal.Length Petal.Width
30 4.7 3.2 1.6 0.2
NewPoint
Sepal.Length Sepal.Width Petal.Length Petal.Width
30 5.7 3.2 1.6 0.2
Вы можете видеть, что NewPoint
совпадает с исходной точкой iris[30,1:4]
за исключением того, что Sepal.Length был увеличен на 1. Как это влияет на прогноз?
predict(Mod1, newdata=iris[30,1:4])
30
1.386358
predict(Mod1, newdata=NewPoint)
30
2.115497
predict(Mod1, newdata=NewPoint) - predict(Mod1, newdata=iris[30,1:4])
30
0.7291384
Разница в предсказанных значениях составляет 0,7291384, что является коэффициентом для Sepal.Length, показанным выше.
У меня есть два дополнительных вопроса: 1) влияет ли фактор случайного эффекта на интерпретацию коэффициента модели? 2) влияет ли на интерпретацию модельного коэффициента тот факт, что BirthDate
стандартизован?
1. Нет. Модель обеспечивает прогноз на основе переменных-предикторов. Случайная часть влияет на ошибку в этих прогнозах, но не меняет прогноз. 2. Да. Наклон (скорость изменения), заданный вашей моделью, относится к переменной, на которой вы ее построили, - стандартизованным датам. Таким образом, наклон модели для BirthDate 0.43022
относится к изменению стандартизированный BirthDates. Чтобы найти изменение показателя активности при изменении исходных дат рождения на единицу, необходимо масштабировать результат. Смена будет (1 / sd(BirthDate)) * 0.43022
Ok. В этом есть смысл. Я не понимаю, почему тогда я получаю несколько выходов для строки predict(m1, data=merge.data)
(я думал, что это из-за случайного эффекта)?
(Я должен уточнить, это строка tapply(predict(m1,type = "response"), merge.data$BirthDate,mean)
в моем ОП)
Что вы получите, если наберете
ml
?