Diferencia en los modelos de estadísticas de Python OLS y R''s lm

ols in python statsmodels

Question

No estoy seguro de por qué obtengo resultados ligeramente diferentes para un OLS simple, dependiendo de si utilizo la interfaz rpy experimental de panda para realizar la regresión en R o si uso statsmodels en Python.

import pandas
from rpy2.robjects import r
from functools import partial
loadcsv = partial(pandas.DataFrame.from_csv,
index_col="seqn", parse_dates=False)
demoq = loadcsv("csv/DEMO.csv")
rxq = loadcsv("csv/quest/RXQ_RX.csv")
num_rx = {}
for seqn, num in rxq.rxd295.iteritems():
try:
val = int(num)
except ValueError:
val = 0
num_rx[seqn] = val
series = pandas.Series(num_rx, name="num_rx")
demoq = demoq.join(series)
import pandas.rpy.common as com
df = com.convert_to_r_dataframe(demoq)
r.assign("demoq", df)
r(''lmout <- lm(demoq$num_rx ~ demoq$ridageyr)'') # run the regression
r(''print(summary(lmout))'') # print from R

De R , obtengo el siguiente resumen:

Call:
lm(formula = demoq$num_rx ~ demoq$ridageyr)
Residuals:
Min 1Q Median 3Q Max
-2.9086 -0.6908 -0.2940 0.1358 15.7003
Coefficients:
Estimate Std. Error t value Pr(>|t|)
(Intercept) -0.1358216 0.0241399 -5.626 1.89e-08 ***
demoq$ridageyr 0.0358161 0.0006232 57.469 < 2e-16 ***
---
Signif. codes: 0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1
Residual standard error: 1.545 on 9963 degrees of freedom
Multiple R-squared: 0.249, Adjusted R-squared: 0.2489
F-statistic: 3303 on 1 and 9963 DF, p-value: < 2.2e-16

Usando statsmodels.api para hacer el OLS:

import statsmodels.api as sm
results = sm.OLS(demoq.num_rx, demoq.ridageyr).fit()
results.summary()

Los resultados son similares a la salida de R pero no son los mismos:

OLS Regression Results
Adj. R-squared: 0.247
Log-Likelihood: -18488.
No. Observations: 9965 AIC: 3.698e+04
Df Residuals: 9964 BIC: 3.698e+04
coef std err t P>|t| [95.0% Conf. Int.]
ridageyr 0.0331 0.000 82.787 0.000 0.032 0.034

El proceso de instalación es un poco engorroso. Pero, hay un cuaderno de ipython here , que puede reproducir la inconsistencia.

Answer 1

Parece que Python no agrega una intercepción de forma predeterminada a tu expresión, mientras que R lo hace cuando usas la interfaz de fórmula ...

Esto significa que encajaste en dos modelos diferentes. Tratar

lm( y ~ x - 1, data)

en R para excluir la intercepción, o en su caso y con una notación algo más estándar

lm(num_rx ~ ridageyr - 1, data=demoq)

Answer 2

Tenga en cuenta que todavía puede usar statsmodels.formula.api desde statsmodels.formula.api :

from statsmodels.formula.api import ols
results = ols(''num_rx ~ ridageyr'', demoq).fit()
results.summary()

Creo que usa patsy en el backend para traducir la expresión de fórmula, y la intercepción se agrega automáticamente.