Create a gist now

Instantly share code, notes, and snippets.

What would you like to do?
Basic econometrics of counts
#-----------------------------------
# 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())
Sign up for free to join this conversation on GitHub. Already have an account? Sign in to comment