ordinalNet正则化有序回归系数符号异常问题求助

我有一个包含有序因子响应（"low"、"medium"、"high"）的数据集，还有大量预测变量，正在用ordinalNet包构建正则化有序回归模型。这是我第一次用这个包，所以先在模拟数据上测试，结果发现系数符号和预期相反。

在模拟数据里，我预期x的系数为正——从绘图能看出来，x值越大，出现medium和high的概率越高。用非正则化方法MASS::polr得到了符合预期的正系数，但ordinalNet得出的x系数是负的，求帮忙！

模拟数据与代码 #

library(tidyverse)
library(MASS)
#> 
#> Attaching package: 'MASS'
#> The following object is masked from 'package:dplyr':
#> 
#>     select
library(ordinalNet)

# simulated data
dat <- 
  expand_grid(
    repeats = 1:10,
    id = 1:3
  ) %>%
  mutate(
    x = id + runif(length(id), -1, 1),
    x2 = runif(length(id), -1, 1),
    y = factor(c("low", "med", "high"), c("low", "med", "high"), ordered = TRUE)[id],
  )
  
# inspect
dat
#> # A tibble: 30 × 5
#>    repeats    id     x     x2 y    
#>      <int> <int> <dbl>  <dbl> <ord>
#>  1       1     1 0.498  0.539 low  
#>  2       1     2 1.31  -0.406 med  
#>  3       1     3 2.30   0.774 high 
#>  4       2     1 0.163  0.408 low  
#>  5       2     2 1.99  -0.504 med  
#>  6       2     3 3.86   0.931 high 
#>  7       3     1 0.822 -0.799 low  
#>  8       3     2 1.96  -0.267 med  
#>  9       3     3 2.83   0.769 high 
#> 10       4     1 0.127  0.139 low  
#> # … with 20 more rows

# 验证x值越高，medium和high的概率越高
dat %>%
  ggplot() + 
  aes(x, y) + 
  geom_point()

模型对比结果 #

# MASS::polr模型（非正则化）
model <- polr(y ~ x + x2, data = dat)
coef(model)
#>         x        x2 
#> 3.6007872 0.6991214

# ordinalNet正则化模型
form <- y ~ x + x2
x <- model.matrix(form, dat)[, -1]
y <- model.frame(form, dat)[, 1, drop = TRUE]
rmodel <- ordinalNetCV(x, y)
#> Fitting ordinalNet on full training data
#> Fitting ordinalNet on fold 1 of 5 
#> Fitting ordinalNet on fold 2 of 5 
#> Fitting ordinalNet on fold 3 of 5 
#> Fitting ordinalNet on fold 4 of 5 
#> Fitting ordinalNet on fold 5 of 5 
#> Done
coef(rmodel$fit)
#> (Intercept):1 (Intercept):2             x            x2 
#>     4.2554900     8.4649827    -3.4500300    -0.6367011
# 为什么系数符号不对？

原因与解决方法 #

核心原因

两个包对有序回归的累积概率定义方向不同：

• MASS::polr默认建模累积概率P(Y ≤ k)：即响应变量小于等于类别k的概率。正系数意味着x增大时，Y≤k的概率降低，对应Y取更高类别的概率升高，符合你的预期。
• ordinalNet默认建模累积概率P(Y ≥ k)：即响应变量大于等于类别k的概率。此时负系数表示x增大时，Y≥k的概率升高（模型公式为logit(P(Y≥k)/(1-P(Y≥k))) = 截距 - βx，x增大时线性预测值降低，P(Y≥k)反而升高），本质和polr的预测趋势一致，只是概率定义方向导致系数符号相反。

解决方法

可以通过ordinalNet的reverse=TRUE参数反转响应变量的顺序，让系数符号与polr一致：

rmodel <- ordinalNetCV(x, y, reverse = TRUE)
coef(rmodel$fit)
#> (Intercept):1 (Intercept):2             x            x2 
#>    -8.4649827    -4.2554900     3.4500300     0.6367011

此时x的系数变为正，和polr的结果趋势完全匹配。

也可以不用调整参数，只要理解概率定义的差异，两个模型的预测结论是一致的——x增大时，Y取high的概率都会升高。

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