Все возможные регрессии в R: сохранение коэффициентов в матрице

Я запускаю код для всех возможных моделей филогенетической обобщенной линейной модели. У меня возникла проблема с извлечением и сохранением бета-коэффициентов для каждой модели.

Я хочу сохранить коэффициенты в матрицу, где столбцы соответствуют определенной переменной, а строки соответствуют формуле. Проблема возникает из-за того, что переменные разные для каждой модели. Таким образом, нельзя просто привязать коэффициенты к матрице.

В приведенном ниже примере показано, где я нахожусь в проблеме:

y = rnorm(10)
inpdv = matrix(c(rnorm(10), runif(10), rpois(10, 1)), ncol = 3)
colnames(inpdv) = c("A", "B", "C")
data = cbind(y, inpdv)

model.mat = expand.grid(c(TRUE,FALSE), c(TRUE,FALSE), c(TRUE,FALSE))
names(model.mat) = colnames(inpdv)

formula = apply(model.mat, 1, function(x)
                       paste(colnames(model.mat)[x], collapse=" + "))
formula = paste("y", formula, sep = " ~ ")
formula[8] = paste(formula[8], 1, sep = "")

beta = matrix(NA, nrow = length(formula), ncol = 3)

for(i in 1:length(formula)){
   fit = lm(formula(formula), data)
   ## Here I want to extract the beta coeffecients here into a n * k matrix
   ## However, I cannot find a way to assign the value to the right cell in the matrix

}

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

Настоящий анализ будет проведен примерно 30 000 раз, поэтому любые советы по эффективности также будут оценены.

Редактировать: например, вывод для модели y ~ a + c должен быть в форме

a NA b 

Где буквы обозначают коэффициент для этой модели. Следующей моделью может быть y ~ b + c, которая затем будет добавлена ​​внизу. Таким образом, результат теперь будет выглядеть так

a  NA b
NA b c

person SamPassmore    schedule 11.09.2013    source источник


Ответы (2)


Как насчет использования names и %in% для подмножества правильных столбцов. Извлеките значения коэффициентов, используя coef. Как это:

beta = matrix(NA, nrow = length(formula), ncol = 3)
colnames(beta) <- colnames(inpdv)

for(i in 1:length(formula)){
   fit = lm(formula(formula[i]), data)
    coefs <- coef(fit)
    beta[ i , colnames(beta) %in% names( coefs ) ] <- coefs[ names( coefs ) %in% colnames( beta ) ]
}
#              A          B         C
#[1,] -0.4229837 -0.0519900 0.3787666
#[2,]         NA  0.7015679 0.0555350
#[3,] -0.4165834         NA 0.3692974
#[4,]         NA         NA 0.1346726
#[5,] -0.2035173  0.7049951        NA
#[6,]         NA  0.7978726        NA
#[7,] -0.2229959         NA        NA
#[8,]         NA         NA        NA

Я думаю, что вполне допустимо использовать здесь цикл for. Для начала, используя что-то вроде lapply, иногда увеличивайте использование памяти по мере запуска все большего количества симуляций. R иногда не помечает объекты из более ранних моделей как мусор, пока цикл lapply не завершится, поэтому иногда может возникать ошибка выделения памяти. Используя цикл for, я обнаружил, что R будет повторно использовать память, выделенную для предыдущей итерации цикла, если это необходимо, поэтому, если вы можете запустить цикл один раз, вы можете запускать его много раз.

Другая причина не использовать цикл for — это скорость, но я бы предположил, что время итерации незначительно по сравнению со временем, необходимым для подбора модели, поэтому я бы использовал его.

person Simon O'Hanlon    schedule 11.09.2013

for(i in 1:length(formula)){
    fit = lm(formula(formula), data)
     beta[i, 1:length(fit$coefficients)] <- fit$coefficients
}

Обновить

Идея: называть ваши столбцы после коэффициентов и присваивать значения столбцам по имени.

Это всего лишь фиктивный пример, но он должен вам помочь: Создайте выходную матрицу:

beta <- matrix(NA,  nrow=7, ncol=4)
colnames(beta) <- c("(Intercept)", 'A', 'B', 'C')

Создайте некоторые фиктивные данные:

 A <- rnorm(10)
 B <- rpois(10, 1)
 C <- rnorm(10, 2)
 Y <- rnorm(10, -1)

Теперь вы можете сделать что-то вроде этого:

fit <- lm(Y ~ A + B + C)
beta[1, names(fit$coefficients)] <- fit$coefficients

fit <- lm(Y ~ A + B)
beta[2, names(fit$coefficients)] <- fit$coefficients

fit <- lm(Y ~ A + C)
beta[3, names(fit$coefficients)] <- fit$coefficients

fit <- lm(Y ~ B + C)
beta[4, names(fit$coefficients)] <- fit$coefficients

fit <- lm(Y ~ A)
beta[5, names(fit$coefficients)] <- fit$coefficients

fit <- lm(Y ~ B)
beta[6, names(fit$coefficients)] <- fit$coefficients

fit <- lm(Y ~ C)
beta[7, names(fit$coefficients)] <- fit$coefficients
person zero323    schedule 11.09.2013
comment
Привет, спасибо за ваш ответ, но это не сработает. Поскольку каждая модель имеет разное количество коэффициентов, у вас не может быть матрицы, в которой каждая строка имеет разную длину. Также каждый коэффициент должен быть выровнен с соответствующим столбцом, поэтому модель с y ~ a + c необходимо будет поместить в матрицу в виде NA b - person SamPassmore; 11.09.2013
comment
Привет, еще раз спасибо за ваше редактирование, это все равно не сработает, потому что это запутает порядок коэффициентов. Ибо если модель имеет форму a + b, а затем a + c, они будут помещены друг под другом, а затем, если модель будет иметь 3 переменные, она потерпит неудачу. Я был бы признателен за проработанный ответ, а не за догадки. - person SamPassmore; 11.09.2013
comment
Что вы думаете об этом? - person zero323; 11.09.2013
comment
Еще раз спасибо за ответ. Это конкретно не решает мою проблему, потому что, как я уже упоминал, это нужно будет сделать 32 000 раз, поэтому выписывать модель каждый раз невозможно, однако ваши имена (fit $ coefficent) - интересная идея, которую я должен попробовать из. Там, где я нахожусь, уже довольно поздно, поэтому мне придется попробовать завтра, но спасибо за ваши усилия! - person SamPassmore; 11.09.2013