GLM logistic regression in Python

From: Bayesian Models for Astrophysical Data, Cambridge Univ. Press

you are kindly asked to include the complete citation if you used this material in a publication

Code 5.2 GLM logistic regression in Python

==================================================

import numpy as np
import statsmodels.api as sm

# Data
x = np.array([13, 10, 15, 9, 18, 22, 29, 13, 17, 11, 27, 21, 16, 14, 18, 8])
y = np.array([1, 1, 1, 0, 0, 1, 0, 1, 1, 0, 0, 0, 1, 1, 0, 0])

X = np.transpose(x)
X = sm.add_constant(X) # add intercept

# Fit
mu = (y + 0.5) / 2 # initialize mu
eta = np.log(mu/(1-mu)) # initialize eta with the Bernoulli link

for i in range(8):
w = mu * (1 - mu); # variance function
z = eta + (y - mu)/(mu * (1 - mu))    # working response
mod = sm.WLS(z, X, weights=w).fit()   # weigthed regression
eta = mod.predict()      # linear predictor
mu = 1/(1 + np.exp(-eta))      # fitted value
print(mod.params)      # print iteration log

# Output
print(mod.summary())

# Write data as dictionary
mydata = {}
mydata['x'] = x
mydata['y'] = y

# fit using glm package

import statsmodels.formula.api as smf

mylogit = smf.glm(formula='y ~ x', data=mydata, family=sm.families.Binomial())
res = mylogit.fit()
print(res.summary())

==================================================

GET SOURCE

Output on screen:

WLS Regression Results
==============================================================================
Dep. Variable: y   R-squared: 0.045
Model: WLS Adj. R-squared: -0.023
Method: Least Squares F-statistic: 0.6658
Date:    Sun, 18 Dec 2016 Prob (F-statistic): 0.428
Time:    20:51:45 Log-Likelihood: -34.192
No. Observations: 16 AIC:    72.38
Df Residuals: 14 BIC: 73.93
Df Model:    1
Covariance Type: nonrobust
==============================================================================
coef    std err t P>|t|    [95.0% Conf. Int.]
------------------------------------------------------------------------------------------------------------------------------------
const 1.2861 1.657    0.776 0.451 -2.269 4.841
x1 -0.0792 0.097 -0.816 0.428 -0.287 0.129
==============================================================================
Omnibus: 15.189 Durbin-Watson: 1.774
Prob(Omnibus): 0.001 Jarque-Bera (JB): 2.191
Skew: -0.089 Prob(JB): 0.334
Kurtosis: 1.196 Cond. No. 51.9
==============================================================================

Warnings:
[1] Standard Errors assume that the covariance matrix of the errors is correctly specified.

Generalized Linear Model Regression Results
==============================================================================
Dep. Variable:   y No. Observations: 16
Model: GLM Df Residuals: 14
Model Family: Binomial Df Model: 1
Link Function: logit Scale: 1.0
Method: IRLS Log-Likelihood: -10.684
Date: Sun, 18 Dec 2016 Deviance: 21.369
Time: 20:51:45 Pearson chi2: 15.9
No. Iterations: 6
==============================================================================
coef std err z P>|z|    [95.0% Conf. Int.]
------------------------------------------------------------------------------------------------------------------------------------
Intercept 1.2861 1.553 0.828 0.408 -1.758 4.330
x -0.0792 0.091 -0.871    0.384 -0.257 0.099
==============================================================================

HSI

HSI