*! version 6.0.0  06feb2000
program define lncure, eclass
	version 6
	if replay() {
		if `"`e(cmd)'"' != "lncure" { error 301 }
		syntax [, Level(integer $S_level) noCOEF noHEADer TR ]
	}
	else {
		syntax [varlist] [if] [in] [fweight pweight iweight] /*
			*/ [, noCOEF CLuster(varname) Dead(varname numeric)/*
			*/ DEBUG FROM(string) noHEADer MLMethod(string) /*
			*/ MLOpt(string) noCONStant/*
			*/ Level(integer $S_level) noLOG /*
			*/ OFFset(varname numeric) /*
			*/ Robust SCore(string) T0(varname numeric) SKIP TR *]

		tokenize `varlist'
		local t `1'
		mac shift
		local rhs `*'

		if "`cluster'"!="" {
			local cluopt cluster(`cluster')
		}
		if "`from'" != "" { local iniopt init(`from') }
		if "`mlmethod'" == "" { local mlmetho = "lf" }
		if "`offset'" !="" { local offopt = "offset(`offset')" }
		mlopts options, `options'

		if "`score'" != "" {
			local n : word count `score'
			if `n' != 2 {
				di in red /*
			*/ "score() invalid:  two new variable names required"
				exit 198
			}
			confirm new var `score'
			local scopt "score(`score')"
		}

		if "`weight'" != "" {
			tempvar wv
			qui gen double `wv' `exp'
			local w [`weight'=`wv']
		}

		tempvar touse
		mark `touse' `w' `if' `in'
		markout `touse' `t' `rhs' `dead' `t0'
		markout `touse' `cluster', strok

		if "`dead'" != "" {
			local sdead "`dead'"
			capture assert `dead'==0 | `dead'==1 if `touse'
			if _rc {
				tempvar mydead
				qui gen byte `mydead' = `dead'!=0 if `touse'
				local dead "mydead"
			}
		}
		else {
			tempvar dead
			qui gen byte `dead'=1
			local sdead 1
		}
		if "`t0'" == "" {
			local t0 0
		}
		capture assert `t0' < `t' if `touse'
		if _rc {
			di in red "`t0' >= `t' in some obs."
			exit 498
		}


		_rmcoll `rhs' `w' if `touse', `constant'
		local rhs "`r(varlist)'"
		global S_1

		global EREGd `dead'
		* global EREGt `t'
		global EREGt0 `t0'

		if "`log'"!="" { local nlog="*" }
		tempvar mysamp num den
		tempname b f V  g b0 bc0
		quietly {
			if "`weight'"=="aweight" | "`weight'"=="pweight" {
				tempvar wvn
				summ `wv' if `touse', meanonly
				gen double `wvn' = `wv'/r(mean)
				local wvngen 1
			}
			else if "`weight'"!="" {
				local wvn `wv'
			}
			else {
				local wvn 1
			}
			gen double `num' = `wvn'*ln(`t') if `touse' & `dead'
			replace `num' = sum(`num')
			global EREGa = `num'[_N]
			drop `num'
			if "`wvngen'"!=""  { drop `wvn' }
		}
		if "`constant'"!="" {
			local skip = "skip"
			`nlog' di in gr "Fitting full model:"
		}
		if "`rhs'" != "" & "`skip'"=="" {
			`nlog' di ""
			`nlog' di in gr "Fitting constant-only model:"
			ml model `mlmetho' lncur_lf (`t': `t'=, `offopt') /*
				*/ /ln_sig /ln_c `w' if `touse', init(_cons=1) /*
				*/ missing collin nopreserve wald(0) `mlopt' /*
				*/ max search(quietly) noout `log' `options'
			local cont continue
			`nlog' di ""
			`nlog' di in gr "Fitting full model:"
		}
		else {
			 local cont wald(1)
		}
		ml model `mlmetho' lncur_lf /*
			*/ (`t': `t'=`rhs' , `offopt' `constant') /*
			*/ /ln_sig /ln_c `w' if `touse', `cont' noout /*
			*/ `robust' `cluopt' `scopt' `iniopt' `mlopt' /*
			*/ missing collin nopreserve /*
			*/ max search(off) `log' `options'

		if "`e(wtype)'" != "" {
			est local wexp `"`exp'"'
		}
		est local title "accelerated failure-time form"
		est local predict lnorma_p
		est local cmd lnormal
		est local t0 "`t0'"
        	est local dead `sdead'

		global S_E_cmd lnormal
	}
	global EREGw
	global EREGd
	global EREGt
	global EREGt0
	global EREGa

	if "`coef'"=="" {
		if `"`tr'"'!=`""' {
			local hr `"eform(Tm. Ratio)"'
                }

		if "`header'" == "" {
			di _n in gr /*
			*/ "Log-normal regression -- entry time `e(t0)'"
			ml di, `hr' level(`level') first
		}
		else	est di, `hr' level(`level') first
		local bg "[ln_sig]_b[_cons]"
		local sg "[ln_sig]_se[_cons]"
		local iz = invnorm(1-(1-`level'/100)/2)
		local ll = `bg'-`sg'*`iz'
		local uu = `bg'+`sg'*`iz'
		local zz = `bg'/`sg'

		local bc "[ln_c]_b[_cons]"
		local sc "[ln_c]_se[_cons]"
		local lc = `bc'-`sc'*`iz'
		local uc = `bc'+`sc'*`iz'
		local zc = `bc'/`sc'


		est scalar sigma = exp(`bg')
		est scalar cure = exp(`bc')

                local se = sqrt(exp(`bg')*exp(`bg')*`sg'*`sg')
                local sec = sqrt(exp(`bc')*exp(`bc')*`sc'*`sc')

		di in gr   " /ln_sig |  "  in ye %9.0g `bg' /*
			*/ "  " %9.0g (`sg') " " %10.3f (`zz') /*
			*/ " " %7.3f (normprob(-abs(`zz'))*2) /*
			*/ _col(58) %9.0g (`ll') "   " %9.0g (`uu')
		di in gr   "   /ln_c |  "  in ye %9.0g `bc' /*
			*/ "  " %9.0g (`sc') " " %10.3f (`zc') /*
			*/ " " %7.3f (normprob(-abs(`zc'))*2) /*
			*/ _col(58) %9.0g (`lc') "   " %9.0g (`uc')

		di in gr _dup(78) "-"
		di in gr   "   sigma |  "  in ye %9.0g exp(`bg') /*
                        */ "  " %9.0g (`se') " "  /*
			*/ _col(58) %9.0g exp(`ll') "   " %9.0g exp(`uu')
		di in gr   "   cure  |  "  in ye %9.0g exp(`bc') /*
                        */ "  " %9.0g (`sec') " "  /*
			*/ _col(58) %9.0g exp(`lc') "   " %9.0g exp(`uc')

		di in gr _dup(78) "-"
	}
end