Оптимизируйте время для запуска фиксированных эффектов в модели R lm ()

Question

Я пытаюсь запустить регрессионную модель, которая включает фиксированные эффекты для городов в Соединенных Штатах. У меня более 10 000 000 миллионов строк и 600 городов. Код ниже работает, но он действительно медленный. При включении множителя для переменной с большим количеством уровней, есть ли способ запустить модель быстрее. часть as.factor (). @StupidWolf предоставил решение в комментариях. Он предложил функцию slm в library( SparseM ) Это гораздо быстрее, см .:

> system.time(a1 <- lm( a~cityfips+d  , x ) )
   user  system elapsed 
   0.22    0.00    0.22 
> system.time(a2 <- lm( a~as.factor(cityfips) + d  , x ) )
   user  system elapsed 
  95.65    0.97   96.62 
> system.time(a3 <- slm( a~as.factor(cityfips) + d  , x ) )
   user  system elapsed 
   4.58    2.06    6.65 

Answer 1

Когда у вас столько факторов, построение model.matrix в lm() займет большую часть времени, один из способов - использовать sparseMatrix, как в glmnet, и есть два пакета: sparseM, MatrixModels который позволяет lm на sparseMatrix:

set.seed(111)
x <- data.frame(
    a = sample( 1:1000, 1000000 , replace=T),
    cityfips = sample( 1:250, 1000000 , replace=T),
    d = sample( 1:4, 1000000 , replace=T)
)

library(SparseM)
library(MatrixModels)
library(Matrix)

system.time(f_lm <- lm( a~as.factor(cityfips) + d  , x )  )
   user  system elapsed 
 75.720   2.494  79.365  
system.time(f_sparseM <- slm(a~as.factor(cityfips) + d  , x ))
   user  system elapsed 
  5.373   3.952  10.646
system.time(f_modelMatrix <- glm4(a~as.factor(cityfips) + d  ,data=x,sparse=TRUE))
   user  system elapsed 
  1.878   0.335   2.219

Наиболее близким, который я могу найти, является glm4 в MatrixModels, но вы можете увидеть, что ниже коэффициенты такие же, как и при использовании lm:

all.equal(as.numeric(f_sparseM$coefficients),as.numeric(f_lm$coefficients))
[1] TRUE
all.equal(as.numeric(f_lm$coefficients),as.numeric(coefficients(f_modelMatrix)))
[1] TRUE

Еще одна опция, кроме glm4 в MatrixModels, заключается в использовании lm.fit (как указано @BenBolker:

lm.fit(x=Matrix::sparse.model.matrix(~as.factor(cityfips) + d,data=x),y=x$a)

. Это дает вам список, как обычно lm.fit(), и вы не может применять такие функции, как summary() et c.

Авторы обоих пакетов предупреждают о том, что он экспериментальный, поэтому могут существовать некоторые различия по сравнению с stats::lm, позаботьтесь о проверке.