use stepAIC в списке моделей

Question

Я хочу сделать ступенчатую регрессию с помощью AIC в списке линейных моделей. идея состоит в том, чтобы использовать e список линейных моделей, а затем применить stepAIC для каждого элемента списка. Это не удается.

Привет, ребята, я попытался отследить проблему. Я думаю, что нашел проблему. Однако я не понимаю причину. Попробуйте код, чтобы увидеть разницу между тремя случаями.

require(MASS) 
n<-30 
x1<-rnorm(n, mean=0, sd=1) #create rv x1 
x2<-rnorm(n, mean=1, sd=1) 
x3<-rnorm(n, mean=2, sd=1) 
epsilon<-rnorm(n,mean=0,sd=1) # random error variable 
dat<-as.data.frame(cbind(x1,x2,x3,epsilon)) # combine to a data frame 
dat$id<-c(rep(1,10),rep(2,10),rep(3,10)) 
# y is combination from all three x and a random uniform variable 
dat$y<-x1+x2+x3+epsilon 
# apply lm() only resulting in a list of models 
dat.lin.model.lst<-lapply(split(dat,dat$id),function(d) lm(y~x1+x2+x3,data=d)) 
stepAIC(dat.lin.model.lst[[1]]) # FAIL!!! 
# apply function stepAIC(lm())- works 
dat.lin.model.stepAIC.lst<-lapply(split(dat,dat$id),function(d) stepAIC(lm(y~x1+x2+x3,data=d))) 
# create model for particular group with id==1 
k<-which(dat$id==1) # manually select records with id==1 
lin.model.id1<-lm(dat$y[k]~dat$x1[k]+dat$x2[k]+dat$x3[k]) 
stepAIC(lin.model.id1) # check stepAIC - works! 

Я уверен, что stepAIC() нуждается в исходных данных из data.frame "dat". Это то, о чем я думал раньше. (Надеюсь, что я прав) Но в stepAIC() нет параметра, где я могу передать исходный фрейм данных. Очевидно, что для простых моделей, не завернутых в список, достаточно пройти модель. (последние три строки в коде) Итак, мне интересно:
Q1: Как stepAIC знает, где найти исходные данные «dat» (не только данные модели, которые передаются как параметр)?
Q2: Как я могу знать, что есть еще один параметр в stepAIC(), который явно не указан на страницах справки? (возможно, мой английский слишком плохо, чтобы найти)
Q3: Как передать этот параметр stepAIC()?

Это должно быть где-то в среде функции приложения и передачи данных. Либо lm(), либо stepAIC(), и указатель/ссылка на необработанные данные должны где-то потеряться. Я не очень хорошо понимаю, что такое среда в R. Для меня это было как бы изолирование локальных от глобальных переменных. Но может быть, это более сложно. Любой, кто может объяснить это мне в связи с проблемой выше? Честно говоря, я не много читаю из R documentation. Любое лучшее понимание поможет мне. Спасибо.

OLD: У меня есть данные в dataframe df, которые можно разбить на несколько подгрупп. Для этого я создал идентификатор groupID с именем df $ id. lm() возвращает коэффициент, как и ожидалось для первой подгруппы. Я хочу сделать ступенчатую регрессию с использованием AIC в качестве критерия для каждой подгруппы отдельно. Я использую lmList {lme4}, что приводит к модели для каждой подгруппы (id). Но если я использую stepAIC {MASS} для элементов списка, он выдает ошибку. Смотри ниже.

Вопрос: Какая ошибка в моей процедуре/синтаксисе? Я получаю результаты для отдельных моделей, но не для тех, которые созданы с помощью lmList. LmList() хранит различную информацию о модели, чем lm()?
Но в его помощи говорится: класс «lmList»: список объектов класса lm с общей моделью.

>lme4.list.lm<-lmList(formula=Scherkraft.N~Gap.um+Standoff.um+Voidflaeche.px |df$id,data = df) 
>lme4.list.lm[[1]] 
Call: lm(formula = formula, data = data) 
Coefficients: 
(Intercept)   Gap.um  Standoff.um Voidflaeche.px 
    62.306133  -0.009878  0.026317  -0.015048 

>stepAIC(lme4.list.lm[[1]], direction="backward") 
#stepAIC on first element on the list of linear models 
Start: AIC=295.12 
Scherkraft.N ~ Gap.um + Standoff.um + Voidflaeche.px 
       Df Sum of Sq RSS AIC 
- Standoff.um  1  2.81 7187.3 293.14 
- Gap.um   1  29.55 7214.0 293.37 
<none>      7184.4 295.12 
- Voidflaeche.px 1 604.38 7788.8 297.97 

Error in terms.formula(formula, data = data) : 
'data' argument is of the wrong type 

Очевидно, что-то не работает со списком. Но я понятия не имею, что это может быть. Поскольку я пытался сделать то же самое с базовым пакетом, который создает ту же модель (по крайней мере, те же коэффициенты). Результаты приведены ниже:

>lin.model<-lm(Scherkraft.N ~ Gap.um + Standoff.um + Voidflaeche.px,df[which(df$id==1),]) 
# id is in order, so should be the same subgroup as for the first list element in lmList 

Coefficients: 
(Intercept) Gap.um Standoff.um Voidflaeche.px 
    62.306133 -0.009878  0.026317  -0.015048 

Ну, это то, что я получаю с помощью stepAIC на моем линейном мотеле. Насколько я знаю, критерий информации akaike можно использовать для оценки оптимальной балансировки модели между подгонкой и обобщением с учетом некоторых данных.

>stepAIC(lin.model,direction="backward") 
Start: AIC=295.12 
Scherkraft.N ~ Gap.um + Standoff.um + Voidflaeche.px 
       Df Sum of Sq RSS AIC 
- Standoff.um  1  2.81 7187.3 293.14 
- Gap.um   1  29.55 7214.0 293.37 
<none>      7184.4 295.12 
- Voidflaeche.px 1 604.38 7788.8 297.97 

Step: AIC=293.14 
Scherkraft.N ~ Gap.um + Voidflaeche.px 
       Df Sum of Sq RSS AIC 
- Gap.um   1  28.51 7215.8 291.38 
<none>      7187.3 293.14 
- Voidflaeche.px 1 717.63 7904.9 296.85 

Step: AIC=291.38 
Scherkraft.N ~ Voidflaeche.px 
       Df Sum of Sq RSS AIC 
<none>      7215.8 291.38 
- Voidflaeche.px 1 795.46 8011.2 295.65 
Call: lm(formula = Scherkraft.N ~ Voidflaeche.px, data = df[which(df$id == 1), ]) 

Coefficients: 
(Intercept) Voidflaeche.px 
    71.7183   -0.0151 

Я прочитал с выхода следует использовать модель: Scherkraft.N ~ Voidflaeche.px, потому что это минимальный AIC. Хорошо, было бы хорошо, если бы кто-то мог коротко описать выход.Мое понимание ступенчатой ​​регрессии (при условии обратной отмены) - все регрессоры включены в исходную модель. Затем устраняется наименее важный. Критерием для решения является AIC. и так далее ... Как-то у меня возникают проблемы, чтобы правильно интерпретировать таблицы. Было бы хорошо, если бы кто-то мог подтвердить мою интерпретацию. «-» (минус) означает исключенный регресс. В верхней части находится «стартовая» модель и в таблице таблицы ниже RSS и AIC рассчитываются для возможных исключений. Итак, первая строка в первой таблице говорит о модели Scherkraft.N ~ Gap.um + Standoff.um + Voidflaeche.px- Standoff.um приведет к созданию AIC 293.14. Выберите один без Standoff.um: Scherkraft.N ~ Gap.um + Voidflaeche.px

EDIT:
Я заменил lmList {lme4} с dlply() для создания списка моделей. Все еще stepAIC не справляется со списком. Он бросает другой Ошибка. На самом деле, я считаю, что это проблема с шагом данных, который должен пройти AIC. Мне было интересно, как он вычисляет значение AIC для каждого шага только из данных модели. I возьмет исходные данные, чтобы построить модели, оставляя один регрессор каждый раз. В этой связи я бы вычислил AIC и сравнил. Итак, как stepAIC работает, если у него нет доступа к исходным данным. (Я не могу увидеть параметр, где я передаю исходные данные stepAIC). Тем не менее, я не знаю, почему он работает с простой моделью, но не с моделью, завернутой в список.

>model.list.all <- dlply(df, .id, function(x) 
    {return(lm(Scherkraft.N~Gap.um+Standoff.um+Voidflaeche.px,data=x)) }) 
>stepAIC(model.list.all[[1]]) 
Start: AIC=295.12 
Scherkraft.N ~ Gap.um + Standoff.um + Voidflaeche.px 
       Df Sum of Sq RSS AIC 
- Standoff.um  1  2.81 7187.3 293.14 
- Gap.um   1  29.55 7214.0 293.37 
<none>      7184.4 295.12 
- Voidflaeche.px 1 604.38 7788.8 297.97 
Error in is.data.frame(data) : object 'x' not found 

Answer 1

Я не уверен, что может быть изменен в версионности, чтобы сделать отладку так трудно, но одно решение было бы использовать do.call, который вычисляет выражения в вызове перед его выполнением. Это означает, что вместо того, чтобы хранить только d в вызове, чтобы update и stepAIC нужно было найти d, чтобы выполнить свою работу, он хранит полное представление самого кадра данных.

То есть, сделать

do.call("lm", list(y~x1+x2+x3, data=d)) 

вместо

lm(y~x1+x2+x3, data=d) 

Вы можете увидеть, что он пытается сделать, глядя на call элемент модели, возможно, как это:

dat.lin.model.lst <- lapply(split(dat, dat$id), function(d) 
          do.call("lm", list(y~x1+x2+x3, data=d))) 
dat.lin.model.lst[[1]]$call 

Также возможно составить свой список фреймов данных в глобальной среде, а затем пропустить toke, так что update и stepAIC в свою очередь ищут каждый фрейм данных, поскольку их цепи окружения всегда возвращаются к глобальной среде; как это:

dats <- split(dat, dat$id) 
dat.lin.model.list <- lapply(seq_along(dats), function(d) 
      do.call("lm", list(y~x1+x2+x3, data=call("[[", quote(dats),i)))) 

Чтобы увидеть, что изменилось, запустите dat.lin.model.lst[[1]]$call снова.

Answer 2

Как кажется, что stepAIC выходит из среды контура (то есть в глобальной среде) искать данные ему нужно, я обмануть его с помощью присваивания функции:

results <- do.call(rbind, lapply(response, function (i) { 
    assign("i", response, envir = .GlobalEnv) 
      mdl <- gls(as.formula(paste0(i,"~",paste(expvar, collapse = "+")), data= parevt, correlation = corARMA(p=1,q=1,form= ~as.integer(Year)), weights= varIdent(~1/Linf_var), method="ML") 
      mdl <- stepAIC(mdl, direction ="backward") 
}))