glmer拟合二项式混合模型时PIRLS步长减半错误及dredge诉求

我尝试运行以下glmer代码拟合二项式广义线性混合模型：

M_MAR <- glmer(CLASS.1_2_3. ~ DISTANCE * TREATMENT * FAxSmean * SRxS + ROUND + 
    (1 | SITE/PLANT_CODE), family = binomial, 
     na.action = "na.fail", data = STR_SRFA_MAR)

运行后出现错误：

Error: (maxstephalfit) PIRLS step-halvings failed to reduce deviance in pwrssUpdate

其中CLASS.1_2_3.是二项式因变量，其余解释变量包含连续型和分类型。我后续计划使用dredge函数进行模型筛选；目前移除交互项中的两个变量后模型可以正常运行，但希望保留更多交互项的同时使用dredge。

数据集STR_SRFA_MAR的汇总信息如下：

APE    :194   Min.   :  0.00   Length:1519        APE_FRA_B_NET_D5_W1 :  12   Min.   : 0.020  
 HOOG   :194   1st Qu.: 25.00   Class :character   HOUT_FRA_A_NET_D5_W2:  11   1st Qu.: 1.170  
 ZUN    :175   Median : 50.00   Mode  :character   HOOG_FRA_B_NET_D2_W2:  10   Median : 7.535  
 MAR    :159   Mean   : 48.66                      ROO_FRA_A_NET_D1_W2 :  10   Mean   : 9.010  
 APP    :150   3rd Qu.: 80.00                      ROO_FRA_B_NO_D4_W1  :  10   3rd Qu.:14.557  
 ROO    :138   Max.   :105.00                      ROO_FRA_B_NO_D5_W1  :  10   Max.   :55.340  
 (Other):509                                       (Other)             :1456   NA's   :1       
      TSS        MALFORMATIONS.A_B_C.  CLASS.1_2_3.        DATE           EATEN.1_2_3_4.  OTHERS 
 Min.   : 1.10   A   :616             Min.   :0.0000   Length:1519        Min.   :1.000   0:927  
 1st Qu.: 5.80   B   :264             1st Qu.:0.0000   Class :character   1st Qu.:1.000   1:592  
 Median : 7.60   C   :153             Median :1.0000   Mode  :character   Median :2.000          
 Mean   : 7.56   NA's:486             Mean   :0.5286                      Mean   :1.944          
 3rd Qu.: 8.90                        3rd Qu.:1.0000                      3rd Qu.:3.000          
 Max.   :26.60                        Max.   :1.0000                      Max.   :4.000          
 NA's   :638                                                              NA's   :1447           
 TREATMENT ROUND     D            AB             DISTANCE_TOG      DISTANCE_F     FAxSmean     
 NET:749   W1:711   D1:356   Length:1519        Min.   :  0.00   30     :200   Min.   : 104.1  
 NO :770   W2:808   D2:350   Class :character   1st Qu.: 25.00   0      :184   1st Qu.: 393.2  
                    D3:261   Mode  :character   Median : 50.00   5      :172   Median : 704.8  
                    D4:275                      Mean   : 46.17   25     :150   Mean   :1630.1  
                    D5:277                      3rd Qu.: 75.00   75     :148   3rd Qu.:1519.1  
                                                Max.   :100.00   50     :140   Max.   :8506.5  
                                                                 (Other):525                    
    FAxDsum             SRxS            AGRI            GRAS             SN              WB        
 Min.   :    0.2   Min.   : 4.00   Min.   : 3.00   Min.   : 7.00   Min.   : 4.00   Min.   : 1.000  
 1st Qu.:   21.0   1st Qu.: 9.00   1st Qu.:14.00   1st Qu.:13.00   1st Qu.: 8.00   1st Qu.: 1.000  
 Median :  128.6   Median :13.00   Median :24.00   Median :21.00   Median : 8.00   Median : 4.000  
 Mean   : 1571.5   Mean   :13.56   Mean   :23.62   Mean   :24.69   Mean   :14.31   Mean   : 5.562  
 3rd Qu.:  727.1   3rd Qu.:18.00   3rd Qu.:30.00   3rd Qu.:40.00   3rd Qu.:15.00   3rd Qu.: 5.000  
 Max.   :23389.5   Max.   :23.00   Max.   :52.00   Max.   :49.00   Max.   :37.00   Max.   :31.000  
                                                                                                   
    BUILT_UP    
 Min.   : 3.00  
 1st Qu.: 5.00  
 Median :11.00  
 Mean   :19.73  
 3rd Qu.:29.00  
 Max.   :56.00  

解决建议 #

1. 处理模型收敛问题 #

你遇到的PIRLS迭代失败，大多是模型过于复杂（比如四阶交互项带来的参数爆炸），或者连续变量尺度差异过大导致数值不稳定。可以试试这些方法：

• 标准化连续变量：把DISTANCE、FAxSmean、SRxS这类连续变量用scale()函数做标准化，缩小尺度差异，降低迭代时的数值波动。
• 简化高阶交互项：四阶交互项在实际研究中极少有统计学意义，还会大幅增加模型复杂度。可以先保留二阶或三阶交互项，验证模型能收敛后，再逐步尝试加入更高阶交互。
• 调整迭代控制参数：换用更适合混合模型的优化器，同时增加最大迭代次数：

M_MAR <- glmer(CLASS.1_2_3. ~ DISTANCE * TREATMENT * FAxSmean * SRxS + ROUND + 
    (1 | SITE/PLANT_CODE), family = binomial, 
     na.action = "na.fail", data = STR_SRFA_MAR,
     control = glmerControl(optimizer = "bobyqa", optCtrl = list(maxfun = 1e5)))

2. 配合dredge的实用策略 #

dredge要求所有候选模型都能正常拟合，全交互模型收敛失败的话直接用会报错，可按以下方式调整：

• 限定候选模型复杂度：先拟合包含所有主效应和二阶交互的基础模型（确保能收敛），再用dredge筛选，之后根据结果逐步加入高阶交互验证。
• 手动筛选高价值交互项：结合领域知识或单变量分析，挑出可能有意义的交互项纳入候选集，减少dredge需要拟合的模型数量，同时降低收敛风险。
• 用subset参数限制模型范围：比如只允许最多包含1个三阶交互，避免拟合过于复杂的模型：

library(MuMIn)
# 先拟合基础模型
M_base <- glmer(CLASS.1_2_3. ~ (DISTANCE + TREATMENT + FAxSmean + SRxS)^2 + ROUND + 
    (1 | SITE/PLANT_CODE), family = binomial, 
     na.action = "na.fail", data = STR_SRFA_MAR)
# 筛选时限制最多一个三阶交互
dredge(M_base, subset = !((term(DISTANCE:TREATMENT:FAxSmean) && term(DISTANCE:TREATMENT:SRxS))),
       rank = AICc)

3. 排查数据问题 #

从数据集汇总看，部分变量有大量NA值（比如MALFORMATIONS.A_B_C.有486个NA，TSS有638个NA），虽然用了na.action = "na.fail"，但可以检查SITE/PLANT_CODE的分组样本量，确保每个分组有足够的观测值，避免因分组样本过少导致模型收敛失败。

内容的提问来源于stack exchange，提问作者Sara Piccini