use https://www3.nd.edu/~rwilliam/statafiles/glm-logit.dta, clear
logit grade gpa tuce i.psi, nolog
margins, dydx(*) atmeans
margins psi, atmeans
display .5632555 - .1067571
margins, at(gpa = (3.117188 4.117188)) atmeans
display .8510027 - .2528205
margins, at(gpa = (3.117188 3.118188)) atmeans noatlegend
display .2533547 - .2528205

* Appendix A: Compute AMEs manually
webuse nhanes2f, clear
* For convenience, keep only nonmissing cases
keep if !missing(diabetes, female, age)
clonevar xage = age
sum xage
gen xdelta = r(sd)/1000
logit diabetes i.female xage, nolog
margins, dydx(xage)
predict xage1
replace xage = xage + xdelta
predict xage2
gen xme = (xage2 - xage1) / xdelta
sum xme


* Appendix B: Technical discussion
use https://www3.nd.edu/~rwilliam/statafiles/glm-logit.dta, clear
logit grade gpa tuce psi, nolog
adjust gpa tuce psi, pr
mfx
use https://www3.nd.edu/~rwilliam/statafiles/glm-logit.dta, clear
probit  grade gpa tuce psi, nolog
mfx
* marginal change for GPA.  The invnorm function gives us the Z-score for the stated
* prob of success.  The normalden function gives us the pdf value for that Z-score.
display invnorm(.2658)
display normalden(invnorm(.2658))
display normalden(invnorm(.2658)) * 1.62581
quietly probit  grade gpa tuce i.psi, nolog
margins, dydx(*) atmeans
margins, dydx(*)
use https://www3.nd.edu/~rwilliam/statafiles/glm-reg.dta, clear
reg income  educ jobexp black
margins, dydx(*) atmeans
margins, dydx(*)
reg income  educ jobexp i.black i.black#c.educ
margins, dydx(*)
reg income  educ jobexp i.black c.educ#c.educ
margins, dydx(*)