Python中使用rpy2实现的序列多项式回归: 关于共线预测变量的问题

Question

我试图使用rpy2（Python界面的R语言）调用R的mass.polr函数在Python中执行顺序Logistic回归。然而，当我的预测因子中存在一些共线或几乎共线的列时，mass.polr会自动丢弃其中的一些列来进行拟合，这将导致我尝试在训练数据上进行预测时出错。

以下是一个最小示例：

from rpy2.robjects import r, pandas2ri
from rpy2.robjects.packages import importr

pandas2ri.activate()

mass = importr("MASS")

# dataframe with two collinear predictors (x1 and x2)
df = pd.DataFrame(columns = ['target', 'x1', 'x2', 'x3'],
                  data    = [[   0   ,  0  ,  0  ,  1  ],
                             [   1   ,  1  ,  1  ,  0  ],
                             [   2   ,  1  ,  1  ,  1  ]])

model = mass.polr('as.factor(target) ~ .', df, Hess = True) # gives warning below
'''
Warning message:
In polr(as.factor(target) ~ ., data = df, Hess = TRUE) :
  design appears to be rank-deficient, so dropping some coefs

'''

r.predict(model, df, type = "class").__array__() # gives error below
'''
Error in X %*% object$coefficients : non-conformable arguments
'''

同样的错误在R中也会发生，但是在那里我至少可以通过查看 "summary(model)" 来看到哪些列已被丢弃。而在Python中，"r.summary(model).rx2('coefficients')"（应该显示与R中的 "summary(model)" 相同的输出）不显示系数名称，只显示裸值。

array([[4.57292582e+01, 8.25605929e+02, 5.53887231e-02],
       [2.11604944e+01, 2.85721885e+02, 7.40597606e-02],
       [3.19476895e+01, 3.60605165e+02, 8.85946531e-02],
       [5.66312792e+01, 8.93862000e+02, 6.33557296e-02]])

有没有办法在Python中检索系数名称？还是有其他的解决方法吗？

- Elena Fortina

Answer 1

即使没有使用 pandas2ri.activate()，从 r.summary(model).rx2('coefficients') 返回的 FloatMatrix 中也不包括变量名。然而，我们可以使用 R 的 dimnames 函数提取这些名称。完整示例如下：

import pandas as pd
import rpy2.robjects as ro
from rpy2.robjects import r, pandas2ri
from rpy2.robjects.packages import importr
from rpy2.robjects.conversion import localconverter
mass = importr("MASS")

df = pd.DataFrame(columns = ['target', 'x1', 'x2', 'x3'],
                  data    = [[   0   ,  0  ,  0  ,  1  ],
                             [   1   ,  1  ,  1  ,  0  ],
                             [   2   ,  1  ,  1  ,  1  ]])

with localconverter(ro.default_converter + pandas2ri.converter):
    df = ro.conversion.py2rpy(df)

model = mass.polr('as.factor(target) ~ .', df, Hess = True)

coefs = r.summary(model).rx2('coefficients')

[x for x in r('dimnames')(coefs)[0]]

返回结果为['x1', 'x3', '0|1', '1|2']，表示x2已被删除。

或者，您可以使用r.print(r.summary(model))打印完整的模型输出。

- muchomaas

Answer 2

r.summary(model).rx2('coefficients') 返回一个没有名称的对象，因为您在脚本中较早的位置（行 pandas2ri.activate()）请求将 R 对象转换为 pandas（和隐式的 numpy）对象。Numpy 数组没有命名元素。

不再建议使用 activate。考虑在上下文中使用本地转换器代替（例如，在文档中使用 pandas 的示例：https://rpy2.github.io/doc/v3.3.x/html/generated_rst/pandas.html）。

- lgautier