Basic econometrics of counts

poisson reg.py

#-----------------------------------
# poisson regression basic example
#-----------------------------------

# see also: http://econometricsense.blogspot.com/2017/03/count-model-regressions.html 
 
# this code is adapted from a more extensive GEE application 
# set up to allow for repeated measures 
http://nbviewer.jupyter.org/urls/umich.box.com/shared/static/ir0bnkup9rywmqd54zvm.ipynb

# generate count data
import pandas as pd
import numpy as np

# you can simulate data like this
a = np.random.poisson(3, 15) # treatment group
b =  np.random.poisson(5, 15) # control group

# but lets make a quick and dirty data set

data = {'COUNT' :[3,4,2,2,1,6,0,0,1,5,3,2,3,3,4,5,4,7,8,3,7,5,4,7,8,5,3,5,4,6],
        'TRT':[1,1,1,1,1,1,1,1,1,1,1,1,1,1,1,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0],
                              'ID':[1,2,3,4,5,6,7,8,9,10,11,12,13,14,15,16,17,18,19,20,21,22,23,24,25,26,27,28,29,30],
}

count_data = pd.DataFrame(data,columns=['COUNT','TRT','ID'])

# descriptives
count_data.mean() # overall
groupby_trt = count_data.groupby('TRT') # group by treatment
groupby_trt.mean() # means by group

# fit poisson model to test differences in counts for treated vs untreated group

from statsmodels.genmod.generalized_estimating_equations import GEE
from statsmodels.genmod.cov_struct import (Exchangeable,
    Independence,Autoregressive)
from statsmodels.genmod.families import Poisson

fam = Poisson()
ind = Independence()
model1 = GEE.from_formula("COUNT ~ TRT", "ID", data, cov_struct=ind, family=fam)
result1 = model1.fit()

print(result1.summary())