program nested, rclass
	syntax [varlist(default=none fv)], [stub(string) DESCription(string) ZSHOW]
	* syntax:
	* yvar xvar othervars otheroptionsifdesired
	* yvar should be a 0/1 dependent variable
	* xvar should be a dichotomous categerorical variable or a continuous variable
	* othervars can be categerorical or continuous
	* stub is user to partially specify the outpot file name, i.e.
	*    the output file is named ssr`stub'.rtf
	* description is used in the table titles
	* zshow is used in the KHB analysis. By default, the
	*   Z-residual variables are temporary variables and are not
	*   shown in the output. If zshow is specified then the
	*   Z-residual variables are included in the output
	* NOTE: Any variables and estimation results created by the program 
	*   are not saved. You can modify the program if you want them.
	* NOTE: rtf Output files will likely need some additional editing.
	
	* Only cases with nonmissing data are kept
	marksample touse
	preserve
	keep if `touse'

	gettoken y z: varlist
	gettoken x z: z
	if "`description'" == "" local description TBA
	if "`stub'" == "" local stub TBA

	* This is a check to see if x is a factor variable and if so,
	* use the 2nd column for its coefficient value
	quietly reg `y' `x'
	if e(b)[1,1] != 0 {
		local startcol = 1
	}
	else local startcol = 2

	********************************************
	* 1. Run Logistic regressions
	display " "
	display as result "1. Logistic regressions"
	display " "

	logit `y' `x', nolog
	AddInfo `startcol'
	est store Xonly
	local xbefore = e(xcoef)

	logit `y' `z', nolog
	AddInfo `startcol'
	est store Zonly

	logit `y' `x' `z', nolog
	AddInfo `startcol'
	est store Both
	local xafter = e(xcoef)	
	AddNote `x' `xbefore' `xafter'

	esttab Xonly Zonly Both using ssr`stub'.rtf, /// 
		pr2 b(%8.2f) scalar(VarYstar SDYstar) nodepvars nonumbers replace ///
		noconstant obslast nobase z  ///
		title(Table LPM: Logistic Regressions - `description') ///
		addnote(`=r(note)')
	
	********************************************
	* 2. OLS/LPM regressions
	display " "
	display as result "2. Linear Probability Models"
	display " "
	
	sum `y'
	local vary = r(Var)
	local sdy = r(Var)^.5

	reg `y' `x'
	estadd scalar VarY = `vary'
	estadd scalar SDY = `sdy'
	est store LPMXonly
	local xbefore = e(b)[1,`startcol']

	reg `y' `z'
	estadd scalar VarY = `vary'
	estadd scalar SDY = `sdy'
	est store LPMZonly

	reg `y' `x' `z'
	estadd scalar VarY = `vary'
	estadd scalar SDY = `sdy'
	est store LPMBoth
	local xafter = e(b)[1,`startcol']
	
	AddNote `x' `xbefore' `xafter'

	esttab LPMXonly LPMZonly LPMBoth using ssr`stub'.rtf, /// 
		r2 b(%8.2f) scalar(VarY SDY) nodepvars nonumbers append ///
		noconstant obslast nobase  ///
		mtitles (LPM-Xonly LPM-Zonly LPM-Both) ///
		title(Table LPM: Linear Probability Models - `description') ///
		addnote(`=r(note)')
		
	********************************************
	* 3. Y-standardization
	display " "
	display as result "3. Y-standardization"
	display " "
	
	est restore Xonly
	local xbefore = e(xcoef)/ e(SDYstar)
	est restore Both
	local xafter = e(xcoef)/ e(SDYstar)
	AddNote `x' `xbefore' `xafter'

	esttab Xonly Zonly Both using ssr`stub'.rtf, /// 
		main(bStdY %8.2f) nodepvars nonumbers append ///
		noconstant obslast nobase z  ///
		title(Table Ystd: Y Standardized Coefficients - `description') ///
		addnote(`=r(note)')
		
	* Confirm with listcoef
	est restore Xonly
	listcoef, std
	est restore Zonly
	listcoef, std
	est restore Both
	listcoef, std

	********************************************
	* 4. Do it yourself KHB
	display " "
	display as result  "4. KHB"
	display " "
	
	*** Do NOT show z-residual coefficients ***
	*** Z-residuals will be created as temporary variables
	*** and not shown in the table.
	if "`zshow'" == "" {
		fvrevar `z'
		foreach zvar of varlist `r(varlist)' {
			quietly sum `zvar'
			if r(min) != r(max) {
				quietly reg `zvar' `x'
				tempvar zresid
				predict `zresid', residual
				local zresids `zresids' `zresid'
			}
		}
	}
	*** DO show z-residual coefficientS ***
	*** Z-residuals will be created as permanent variables
	*** and shown in the table.
	else if "`zshow'" != "" {
		fvrevar `z', list
		foreach zvar of varlist `r(varlist)' {
			quietly reg `zvar' `x'
			predict `zvar'_resid, residual
			local zresids `zresids' `zvar'_resid
		}
	}
	
	* KHB Reduced model
	quietly logit `y' `x' `zresids', nolog
	AddInfo `startcol'
	est store Reduced
	local xbefore = e(xcoef)

	* KHB Full model
	logit `y' `x' `z', nolog
	AddInfo `startcol'
	est store Full
	local xafter = e(xcoef)
	
	AddNote `x' `xbefore' `xafter'
	
	* Modified Note is added to the table when the user has opted to not show the Z-residuals
	if "`zshow'" == "" {
		local znote ""The z-residual variables are included in the Reduced Model but coefficients are not shown" "`=r(note)'""
	}
	else if "`zshow'" != "" {
		local znote `=r(note)'
		local zresids
	}

	esttab Reduced Full using ssr`stub'.rtf, /// 
		pr2 b(%8.2f) scalar(VarYstar SDYstar) nodepvars nonumbers append ///
		noconstant obslast nobase z  ///
		drop(`zresids') /// 		
		title(Table KHB: KHB Models - `description') ///
		addnote (`znote')

	* Confirm with KHB Program
	khb logit `y' `x' || `z', nolog

	********************************************
	* 5. Marginal Effects
	display " "
	display as result "5A. Conventional Marginal Effects"
	display " "
	
	foreach estimates in Xonly Zonly Both Reduced Full {
		est restore `estimates'
		margins, dydx(_all) post
		estadd scalar xcoef = e(b)[1,`startcol']
		est store `estimates'marg
	}
	
	est restore Xonlymarg
	local xbefore = e(xcoef)
	est restore Bothmarg
	local xafter = e(xcoef)
	AddNote `x' `xbefore' `xafter'
	
	esttab Xonlymarg Zonlymarg Bothmarg using ssr`stub'.rtf, /// 
		b(%8.2f) nodepvars nonumbers append ///
		noconstant obslast nobase z   ///
		mtitles (Xonly Zonly Both) ///
		title(Table Marg: Marginal Effects - `description') ///
		addnote(`=r(note)')

	display " "
	display as result "5B. KHB Marginal Effects"
	display " "
	
	est restore Reducedmarg
	local xbefore = e(xcoef)
	est restore Fullmarg
	local xafter = e(xcoef)
	AddNote `x' `xbefore' `xafter'	
	
	esttab Reducedmarg Fullmarg using ssr`stub'.rtf, /// 
		b(%8.2f) nodepvars nonumbers append ///
		noconstant obslast nobase z   ///
		drop(`zresids') /// 
		mtitles (Reduced Full) /// 
		title(Table KHBMarg: KHB Marginal Effects - `description') ///
		addnote ("The z-residual variables are included in the Reduced Model but coefficients are not shown" "`=r(note)'") 
	
	* Confirm with KHB Program APE option
	khb logit `y' `x' || `z', nolog ape
		

	********************************************	
end

program AddInfo
	* Adds information involving ystar, y-standardized coefficients
	* to the saved results for a model
	args startcol
	quietly {
		tempname ystar bStdY
		predict `ystar', xb
		quietly sum `ystar'
		local vary = r(Var) + _pi^2/3
		local sdy = `vary' ^ .5
		matrix `bStdY' = e(b) / `sdy'
		local xcoef = e(b)[1, `startcol']
		
		estadd scalar VarYstar = `vary'
		estadd scalar SDYstar = `sdy'
		estadd matrix bStdY = `bStdY'
		estadd scalar xcoef = `xcoef'
	}
end

	********************************************
program AddNote, rclass
	* Adds notes to tables indicating the numeric and percentage changes
	* in coefficients across nested models
	args x xbefore xafter

	local xdiff = `xafter' - `xbefore'
	local xpctchg = round(`xdiff'/ `xbefore', .001) * 100
	local xdiff = round(`xdiff', .01)
	
	if `xpctchg' < 0 {
		local xpctchg = abs(`xpctchg')
		local note "With controls added the coefficient for `x' decreases (in magnitude) by `xdiff', i.e. by `xpctchg'%"
	}
	else if `xpctchg' > 0 {
		local note "With controls added the coefficient for `x' increases by `xdiff', i.e. by `xpctchg'%"
	}
	else {
		local note "When controls are added the coefficient for `x' does not change"
	}
	return local note "`note'"
end