Подтвердить что ты не робот

Большая биномиальная регрессия фиксированных эффектов в R

Мне нужно запустить логическую регрессию на относительно большом кадре данных с 480 000 записей с тремя фиксированными переменными эффекта. Фиксированный эффект var A имеет 3233 уровня, var B имеет 2326 уровней, var C имеет 811 уровней. Итак, у меня есть 6370 фиксированных эффектов. Данные являются поперечными. Если я не могу запустить эту регрессию, используя обычную функцию glm, потому что матрица регрессии кажется слишком большой для моей памяти (я получаю сообщение "Error: cannot allocate vector of size 22.9 Gb" ). Я ищу альтернативные способы запуска этой регрессии на моем Macbook Air (OS X 10.9.5 8GB RAM). У меня также есть доступ к серверу с ОЗУ 16 ГБ.

Я попытался решить проблему несколькими способами, но пока ни один не привел к удовлетворительным результатам:

LFE/felm: Используя функцию регрессии felm lfe package, которая вычитает фиксированные эффекты перед запуском регрессии. Это работает отлично и позволило мне запустить вышеупомянутую регрессию как обычную линейную модель всего за несколько минут. Однако lfe не поддерживает логистические регрессии и glms. Таким образом, felm отлично подходит для получения идеи о пригодности модели для разных моделей, но не работает для окончательных моделей логистической регрессии.

biglm/bigglm: Я думал об использовании bigglm, чтобы разбить мою функцию на более управляемые куски. Однако несколько источников (например, link1, link2, link3) упомянуть, что для того, чтобы это работало, уровни факторов должны быть согласованными по кускам, т.е. каждый кусок должен содержать по крайней мере один из каждого коэффициента каждой факторной переменной. Факторы A и B содержат уровни, которые появляются только один раз, поэтому я не могу разделить наборы на разные куски с согласованными уровнями. Если я удалю 10 факторов фиксированного эффекта A и 8 факторов B (незначительное изменение), у меня будут только факторы с оставшимися 4 + уровнями, а разделение моих данных на 4 куска сделает его намного более управляемым. Тем не менее, мне все равно нужно выяснить, как сортировать мой df таким образом, чтобы гарантировать, что мои 480 000 записей отсортированы на 4 куска, в которых каждый факторный уровень каждого из трех факторов появляется хотя бы один раз.

GlmmGS/glmgs: Функция glmmgs в пакете с тем же именем выполняет вычитание с фиксированными эффектами, как пакет lfe для логистических регрессий с использованием "Гаусса-Зейделя". К сожалению, пакет больше не разрабатывается. Будучи относительно новым для R и не имея глубокого опыта в статистике, я не могу понять выход и не знаю, как его преобразовать таким образом, чтобы дать мне нормальный "размер эффекта", "подходящий образ", Значения интервала ", которые дают сводки регрессии glm.

Я отправил сообщение авторам пакета. Они любезно ответили следующим образом:

Пакет не предоставляет выход в том же формате объекта glm. Однако вы может легко вычислить большую часть статистики соответствия (стандартная ошибка оценки, добротность соответствия) с учетом текущего выпуска (в CRAN версии, я считаю, что текущий вывод является вектором оценки коэффициенты и связанный вектор стандартных ошибок; то же самое для компоненты ковариации, но вам не нужно беспокоиться о них, если вы подходят модели без случайных эффектов). Остерегайтесь только того, что ковариационных матриц, используемых для расчета стандартных ошибок, являются обратная к диагональным блокам матрицы точности, связанная с алгоритм Гаусса-Зейделя, и поэтому они склонны недооценивать стандартные ошибки совместного правдоподобия. Я не поддерживаю пакет, и у меня нет времени, чтобы Детали; семантическая теория, лежащая в основе пакета, может быть найдена в бумага, на которую ссылается в руководстве, все остальное нужно проработать вы с ручкой и бумагой:).

Если кто-нибудь может объяснить, как "легко вычислить большую часть статистики соответствия" таким образом, чтобы кто-то без какого-либо образования в статистике мог понять это (может быть невозможно) или предоставить R-код, который показывает на примере того, как это можно сделать Я был бы очень обязан!

Революционная аналитика: Я установил предприятие революционной аналитики на виртуальную машину, которая имитирует Windows 7 на моем Mac. Программа имеет функцию под названием RxLogit, оптимизированную для больших логистических регрессий. Используя функцию RxLogit, я получаю the error (Failed to allocate 326554568 bytes. Error in rxCall("RxLogit", params) : bad allocation), поэтому функция также кажется слишком запущенной в проблемах с памятью. Однако программное обеспечение позволяет мне выполнять мою регрессию в распределенном вычислительном кластере. Поэтому я мог просто "убить проблему", купив вычислительное время в кластере с большим количеством памяти. Тем не менее, я задаюсь вопросом, предоставляет ли программа аналитики революции какие-либо формулы или методы, которые я не знаю, что позволило бы мне сделать какую-то операцию вычитания с фиксированным эффектом lfe или операцию tunk-tunk, которая принимает факторов.

MatrixModels/glm4: Один человек предложил использовать функцию glm4 пакета MatrixModels с атрибутом sparse = TRUE для ускорения вычисления. Если я запускаю регрессию glm4 со всеми фиксированными эффектами, я получаю ошибку "Error in Cholesky(crossprod(from), LDL = FALSE) : internal_chm_factor: Cholesky factorization failed ". Если я запустил ее только с фиксированными переменными эффекта B OR A и C, вычисление работает и возвращает объект "glpModel". glmmgs У меня есть некоторые проблемы, связанные с преобразованием этого вывода в форму, которая имеет смысл для меня, поскольку стандартный метод summary(), похоже, не работает на нем.

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

4b9b3361

Ответ 1

Отъезд

glmmboot{glmmML}

http://cran.r-project.org/web/packages/glmmML/glmmML.pdf

Существует также хороший документ Брострома и Холмберга (http://cran.r-project.org/web/packages/eha/vignettes/glmmML.pdf)

Вот пример из их документа:

dat <- data.frame(y = rbinom(5000, size = 1, prob = 0.5),
               x = rnorm(5000), group = rep(1:1000, each = 5))
fit1 <- glm(y ~ factor(group) + x, data = dat, family = binomial)

require(glmmML)
fit2 <- glmmboot(y ~ x, cluster = group,data = dat)

Разница во времени вычислений "огромная"!

Ответ 2

Для потомков я также хотел бы рекомендовать пакет speedglm, который я нашел полезным при попытке выполнить логическую регрессию на больших наборах данных. Кажется, он использует примерно половину объема памяти и заканчивается намного быстрее, чем glm.

Ответ 3

Я согласен с тем, кто (@Ben Bolker, я думаю?) предложил вам использовать функцию glm4 из MatrixModels. Во-первых, он решает проблему с памятью, если вы используете аргумент sparse. Плотная матрица дизайна с 480 000 записей и 6370 фиксированных эффектов потребует 6371 * 480.000 * 8 = 24.464.640.000 байт. Тем не менее, ваша матрица дизайна будет очень скудной (многие нули), поэтому вы можете сделать меньшую (в памяти) матрицу дизайна, если используете редкую. Во-вторых, вы можете использовать разреженность, чтобы ускорить оценку.

Что касается опций, быстрый поиск показывает, что speedglm также имеет аргумент sparse, хотя я его не пробовал. Я ключевую вещь с любым методом, с которым вы заканчиваете, заключается в том, что он должен использовать, что ваша матрица дизайна разрежена как для сокращения времени вычисления, так и для уменьшения требований к памяти.

Вы получите ошибку (Error in Cholesky(crossprod(from), LDL = FALSE) : internal_chm_factor: Cholesky factorization failed" error), потому что ваша матрица дизайна является единственной. В этом случае ваша проблема не имеет уникального решения, и некоторым вариантом является объединение некоторых уровней группы, использование модели штрафов или случайных эффектов.

Вы правы, что не похоже, что существует сводный метод для класса glpModel. Хотя, слоты, кажется, имеют очевидное имя, и не должно занять много времени, чтобы получить, например, стандартные ошибки в вашей оценке, вычислить оценку дисперсии и т.д.