use "http://biostat.mc.vanderbilt.edu/wiki/pub/Main/CourseBios312/students.dta"

gen total = daysatt + daysabs

summ

* Histograms
hist daysabs
hist daysabs, by(male)

* Unadjusted model, classical standard errors
poisson daysabs male, exposure(total)

* Unadjusted model, robust standard errors
poisson daysabs male, exposure(total) robust

* Male adjusted for math scores
poisson daysabs male math, expsoure(total) robust

* Male adjusted for language arts scores
poisson daysabs male langarts, exposure(total) robust

* Multivariable model with male, math scores, and language arts scores
poisson daysabs male math langarts, exposure(total) robust

* For predictions in Stata, the following commands will display the
*   regression coefficients from the previosly fit model
di _b[male]
di _b[math]
di _b[langarts]
di _b[_cons]

* I will calculate the predictions by hand for females (yhat0) and males (yhat1)
*  The prediction is just exp('linear predictor')
gen yhat0 = exp(_b[_cons] + _b[male]*0 + _b[math]*48 + _b[langarts]*langarts)*100
gen yhat1 = exp(_b[_cons] + _b[male]*1 + _b[math]*48 + _b[langarts]*langarts)*100

* Scatter plot of the prediction
twoway (line yhat0 langarts, sort) (line yhat1 langarts, sort), ytitle(Expected Number of Days Absent) xtitle(Language Arts Score) legend(on order(1 "Female" 2 "Male"))