lme4::lmer 需要很长时间(超过 48 小时)并且意外中止而不报告任何警告或错误

lme4::lmer takes forever (48 hrs +) and aborts unexpectedly without reporting any warnings or errors

数据集dt由16.6万行组成,其列为: log_price:依赖数值回归量 sku:具有 381 个级别的独立分类回归量 年:具有 15 个级别的独立分类回归量 transaction_type:具有 2 个级别的独立分类回归量 购买者:具有 1001 个级别的独立分类回归量 regressor_01, ..., regressor_04:四个独立的数值回归变量

我们只考虑在数据集中出现超过 200 行的 skus。

model_1 回归花费了几分钟并给出了不错的结果:

model_1<- lmer (formula = log_price ~    
                           0 +
                           transaction_type +
                           regressor_01 +
                           regressor_02 +
                           regressor_03 +
                           regressor_04 +
                           year +
                           (1 | sku) +
                           (1 + year | purchaser)
                         , data = dt)

model _2类似于model_1,不同的是我认为sku是固定效应而不是随机效应:

model_2<- lmer (formula = log_price ~    
                           0 +
                           transaction_type +
                           regressor_01 +
                           regressor_02 +
                           regressor_03 +
                           regressor_04 +
                           year +
                           sku +
                           (1 + year | purchaser)
                         , data = dt)

但是,model_2 a) 花费了超过 48 小时 运行,b) 突然结束而没有发送任何警告或错误,并且 c) 继续优化第九次重要的藻类活动(见下面的输出) ). 在其他情况下,我试图通过以下方式加快速度: control = lmerControl(optimizer = "optimx",calc.derivs = FALSE, optCtrl = list(method = "nlminb", starttests = FALSE, kkt = FALSE)。我放弃了优化,因为它突然结束了,所以我把它拿掉以确保它不是由于控制指令。

我不明白为什么它在 sku 是随机效应时收敛得很好,而在 sku 成为固定效应时不收敛。

我可能做错了什么?有什么建议吗?

最后一行输出:

iteration: 19880
    f(x) = 272182.459680
iteration: 19881
    f(x) = 272182.459677
iteration: 19882
    f(x) = 272182.459672
iteration: 19883
    f(x) = 272182.459669
iteration: 19884
    f(x) = 272182.459669
iteration: 19885
    f(x) = 272182.459672
iteration: 19886
    f(x) = 272182.459665
iteration: 19887
    f(x) = 272182.459665

Session 信息:

R version 4.0.2 (2020-06-22)
Platform: x86_64-w64-mingw32/x64 (64-bit)
Running under: Windows 10 x64 (build 19042)

Matrix products: default

locale:
[1] LC_COLLATE=Portuguese_Brazil.1252  LC_CTYPE=Portuguese_Brazil.1252    LC_MONETARY=Portuguese_Brazil.1252
[4] LC_NUMERIC=C                       LC_TIME=Portuguese_Brazil.1252    

attached base packages:
[1] stats     graphics  grDevices utils     datasets  methods   base     

other attached packages:
 [1] stargazer_5.2.2    ivreg_0.5-0        ggplot2_3.3.2      lme4_1.1-25        Matrix_1.2-18      plm_2.2-5         
 [7] future.apply_1.6.0 future_1.20.1      magrittr_2.0.1     data.table_1.13.2 

loaded via a namespace (and not attached):
 [1] Rcpp_1.0.5          bdsmatrix_1.3-4     lattice_0.20-41     listenv_0.8.0       zoo_1.8-8           digest_0.6.27      
 [7] lmtest_0.9-38       parallelly_1.21.0   R6_2.5.0            cellranger_1.1.0    pillar_1.4.7        Rdpack_2.1         
[13] miscTools_0.6-26    rlang_0.4.8         curl_4.3            readxl_1.3.1        rstudioapi_0.13     minqa_1.2.4        
[19] car_3.0-10          nloptr_1.2.2.2      splines_4.0.2       statmod_1.4.35      foreign_0.8-80      munsell_0.5.0      
[25] tinytex_0.27        numDeriv_2016.8-1.1 compiler_4.0.2      xfun_0.19           pkgconfig_2.0.3     globals_0.13.1     
[31] maxLik_1.4-4        tidyselect_1.1.0    tibble_3.0.4        rio_0.5.16          codetools_0.2-16    crayon_1.3.4       
[37] dplyr_1.0.2         withr_2.3.0         MASS_7.3-51.6       rbibutils_2.0       grid_4.0.2          nlme_3.1-148       
[43] gtable_0.3.0        lifecycle_0.2.0     scales_1.1.1        zip_2.1.1           stringi_1.5.3       carData_3.0-4      
[49] ellipsis_0.3.1      generics_0.1.0      vctrs_0.3.5         optimx_2020-4.2     boot_1.3-25         sandwich_3.0-0     
[55] openxlsx_4.2.3      Formula_1.2-4       tools_4.0.2         forcats_0.5.0       glue_1.4.2          purrr_0.3.4        
[61] hms_0.5.3           abind_1.4-5         parallel_4.0.2      yaml_2.2.1          colorspace_2.0-0    gbRd_0.4-11        
[67] haven_2.3.1        

主要问题是固定 sku 会激增固定效应模型矩阵的大小(X,如果您正在阅读 vignette("lmer"))。随机效应模型矩阵(Z)编码为稀疏指标矩阵;固定效应模型矩阵是密集的。添加 sku 将使固定效应模型矩阵增加 381 列,或 (8*381*166e3)/2^20 = 483 Mb。如果你有可用的内存,那么不一定会杀了你,但毫不奇怪,在一个巨大的密集矩阵上所需的矩阵运算比在它的稀疏等效矩阵上要慢得多.

不清楚“突然结束”是指 R 命令退出(出现错误?)并且 returns 您看到提示,还是整个 R 会话停止。整个 R 会话意外退出通常是 运行 内存不足的症状(至少在 Unix 操作系统上,操作系统将终止进程而不是允许它通过抓取来冻结整个 OS更多内存)。

对此你能做些什么?如果能够指定固定效应设计矩阵应该是稀疏的就好了...

  • 有一个 open issue on GitHub since 2012 关于这个(它解释了一些技术问题......)
  • glmmTMB 的最新版本允许稀疏固定效应模型矩阵
  • 有一些讨论,在 this issue 底部附近有一个答案,关于如何使 sku 稀疏但迫使 sku 之间的方差非常大,这正在有效地使其恢复为固定效果

你也可以试试control=lmerControl(calc.derivs=FALSE);您最后看到的缓慢是蛮力 Hessian 计算。由于年份是 分类 预测变量,因此与 (1+year|purchaser) 关联的协方差矩阵为 15x15,具有 (15*16/2)=120 个参数。如果你能减少它,你的计算会加快很多。例如,您可以使 year 数值化,拟合一个非常复杂的样条模型,并且仍然保存大量参数:(1+splines::ns(year,df=4)|purchaser) 只需要 5*6/2=15 个参数。进一步添加 (1|purchaser:year) 会给你 uncorrelated 购买者(在样条曲线附近)的年份之间的变化,并且只会让你多花一个参数 - 你仍然会减少数量顶级参数(问题最重要的维度)的 7.5 倍。