Подходит для многих моделей GLM: улучшить скорость - PullRequest
Новая
       4

Подходит для многих моделей GLM: улучшить скорость

1 голос
/ 18 октября 2019

Я пишу функцию, подходящую для многих glm моделей. Чтобы просто дать вам некоторое представление о функции, я включил небольшой раздел своего кода. С помощью нескольких пользователей SO функция теперь работает для моего анализа. Однако иногда, особенно когда размер образца относительно невелик, может потребоваться довольно много времени для завершения всего процесса. Чтобы сократить время, я рассматриваю возможность изменения некоторых деталей итеративной максимизации, таких как максимальное количество итераций. Я не нашел способа сделать это, возможно, потому что я все еще не знаком с R терминологией. Будем благодарны за любые предложения сделать это или другие способы сократить время. 

all_glm <- function(crude, xlist, data, family = "binomial", ...) {
  # md_lst include formula for many models to be fitted  
  comb_lst <- unlist(lapply(1:n, function(x) combn(xlist, x, simplify=F)), recursive=F)
  md_lst   <- lapply(comb_lst,function(x) paste(crude, "+", paste(x, collapse = "+")))
  models  <- lapply(md_lst, function(x) glm(as.formula(x), family = family, data = data))
  OR      <- unlist(lapply(models, function(x) broom::tidy(x, exponentiate = TRUE)$estimate[2]))

}

EDIT Благодаря @BenBolker, который направил меня к пакету fastglm, я получил несколько r пакетов, которые могли бы предоставить более быстрые альтернативы glm. Я пробовал fastglm и speedglm. Похоже, что оба быстрее, чем glm на моей машине. 

library(fastglm)
library(speedglm)
# from 
set.seed(1)
n <- 25000
k <- 500
y <- rbinom(n, size = 1, prob = 0.5)
x <- round( matrix(rnorm(n*k),n,k),digits=3)
colnames(x) <-paste("s",1:k,sep = "")
df <- data.frame(y,x)
fo <- as.formula(paste("y~",paste(paste("s",1:k,sep=""),collapse="+")))   

# Fit three models: 
system.time(m_glm <- glm(fo, data=df, family = binomial))
system.time(m_speedglm <- speedglm(fo, data= df, family = binomial()))
system.time(m_fastglm <- fastglm(x, y, family = binomial()))

> system.time(m_glm <- glm(fo, data=df, family = binomial))
   user  system elapsed 
  56.51    0.22   58.73 
> system.time(m_speedglm <- speedglm(fo, data= df, family = binomial()))
   user  system elapsed 
  17.28    0.04   17.55 
> system.time(m_fastglm <- fastglm(x, y, family = binomial()))
   user  system elapsed 
  23.87    0.09   24.12

1 Ответ

1 голос
/ 02 ноября 2019

Алгоритм IRLS, обычно используемый для подгонки glms, требует инверсии / разложения матрицы на каждой итерации. fastglm предлагает несколько различных опций для декомпозиции, и выбор по умолчанию - более медленный, но более стабильный вариант (QR с поворотом столбца). Если вас интересует только скорость, то одна из двух доступных декомпозиций типа Холецкого значительно улучшит скорость, что было бы более целесообразно, чем просто изменение количества итераций IRLS. Другим заметным отличием между fastglm и стандартными реализациями IRLS является осторожное использование полушагов для предотвращения расхождения (IRLS может на практике расходиться во многих случаях).

Аргумент method fastglm позволяет изменить декомпозицию. вариант 2 дает ванильное разложение Холецкого, а вариант 3 дает немного более устойчивую версию этого. На моем компьютере время для вашего предоставленного примера: 

> system.time(m_glm <- glm(fo, data=df, family = binomial))
   user  system elapsed 
 23.206   0.429  23.689 

> system.time(m_speedglm <- speedglm(fo, data= df, family = binomial()))
   user  system elapsed 
 15.448   0.283  15.756 

> system.time(m_fastglm <- fastglm(x, y, family = binomial(), method = 2))
   user  system elapsed 
  2.159   0.055   2.218 

> system.time(m_fastglm <- fastglm(x, y, family = binomial(), method = 3))
   user  system elapsed 
  2.247   0.065   2.337

Что касается использования метлы с объектами fastglm, я могу рассмотреть это.

Еще одно замечание о разложении: Когда fastglm использует QR-декомпозицию, он работает с матрицей дизайна напрямую. Хотя speedglm технически предлагает QR-декомпозицию, он сначала вычисляет $ X ^ TX $ и декомпозирует его, что является более численно нестабильным, чем QR на X.

...