Regresión logística: Scikit Learn vs Statsmodels

Question

Estoy tratando de entender por qué el resultado de la regresión logística de estos dos bibliotecas da resultados diferentes.

Estoy utilizando el conjunto de datos de UCLA idre tutorial , prediciendo admit basado en en gre , gpa y rank . rank se trata como variable categórica, por lo que se convierte primero en una variable ficticia con rank_1 se cayó. También se añade una columna de intercepción También se añade una columna de intercepción.

py
from patsy import dmatrices
from sklearn.linear_model import LogisticRegression
import pandas as pd
import statsmodels.api as sm

df = pd.read_csv("https://stats.idre.ucla.edu/stat/data/binary.csv")
y, X = dmatrices('admit ~ gre + gpa + C(rank)', df, return_type = 'dataframe')
X.head()
>  Intercept  C(rank)[T.2]  C(rank)[T.3]  C(rank)[T.4]  gre   gpa
0          1             0             1             0  380  3.61
1          1             0             1             0  660  3.67
2          1             0             0             0  800  4.00
3          1             0             0             1  640  3.19
4          1             0             0             1  520  2.93

# Output from scikit-learn
model = LogisticRegression(fit_intercept = False)
mdl = model.fit(X, y)
model.coef_
> array([[-1.35417783, -0.71628751, -1.26038726, -1.49762706,  0.00169198,
     0.13992661]]) 
# corresponding to predictors [Intercept, rank_2, rank_3, rank_4, gre, gpa]

# Output from statsmodels
logit = sm.Logit(y, X)
logit.fit().params
> Optimization terminated successfully.
     Current function value: 0.573147
     Iterations 6
Intercept      -3.989979
C(rank)[T.2]   -0.675443
C(rank)[T.3]   -1.340204
C(rank)[T.4]   -1.551464
gre             0.002264
gpa             0.804038
dtype: float64

La salida de statsmodels es el mismo que se muestra en el sitio web de idre, pero yo no estoy seguro de por qué scikit-learn produce un conjunto diferente de coeficientes. ¿Acaso minimiza alguna función de pérdida diferente? ¿Hay alguna documentación que que indique la implementación?

Answer 1

Su pista para averiguar esto debería ser que las estimaciones de los parámetros de la estimación de scikit-learn son uniformemente más pequeñas en magnitud que la contraparte de statsmodels. Esto podría llevarle a pensar que scikit-learn aplica algún tipo de regularización de los parámetros. Usted puede confirmar esto mediante la lectura de la documentación de scikit-learn .

No hay forma de desactivar la regularización en scikit-learn, pero se puede hacer que no sea efectiva estableciendo el parámetro de ajuste C a un gran número. Así es como funciona en su caso:

# module imports
from patsy import dmatrices
import pandas as pd
from sklearn.linear_model import LogisticRegression
import statsmodels.discrete.discrete_model as sm

# read in the data & create matrices
df = pd.read_csv("http://www.ats.ucla.edu/stat/data/binary.csv")
y, X = dmatrices('admit ~ gre + gpa + C(rank)', df, return_type = 'dataframe')

# sklearn output
model = LogisticRegression(fit_intercept = False, C = 1e9)
mdl = model.fit(X, y)
model.coef_

# sm
logit = sm.Logit(y, X)
logit.fit().params

ACTUALIZACIÓN : Como se ha señalado correctamente en los comentarios, ahora se puede desactivar la relularización en scikit-learn configurando penalty='none' (véase el docs ).

Answer 2

Lo que me hizo tropezar:

• desactivar la regularización de sklearn LogisticRegression(C=1e9)

• añadir intercepción de modelos de estadísticas sm.Logit(y, sm.add_constant(X)) O desactivar la interceptación de sklearn LogisticRegression(C=1e9, fit_intercept=False)

• sklearn devuelve la probabilidad de cada clase, por lo que model_sklearn.predict_proba(X)[:, 1] == model_statsmodel.predict(X)

• uso de la función de predicción model_sklearn.predict(X) == (model_statsmodel.predict(X) > 0.5).astype(int)

Ahora veo los mismos resultados en ambas bibliotecas.

Answer 3

Otra diferencia es que has puesto fit_intercept=False, lo que efectivamente es un modelo diferente. Puede ver que Statsmodel incluye el intercepto. No tener un intercepto seguramente cambia los pesos esperados en las características. Pruebe lo siguiente y vea cómo se compara:

model = LogisticRegression(C=1e9)