*! version 3.2.3  23aug2001
program define ksmapprox, sort
	version 7
	syntax varlist(min=2 max=2 numeric) [if] [in] [, /*
		*/ NUMpnts(int 0) /*
		*/ Line Weight noGraph GENerate(string) BWidth(real 0) /*
		*/ BY(varlist) LOWess LOgit T1title(string) /*
		*/ Symbol(string) Sort Connect(string)  * ]
	local gen "`generate'"
	local generat

	if "`gen'"~="" {
		confirm new var `gen'
	}

	local tmpwgt "`weight'"
	local weight 
	marksample touse
	local weight "`tmpwgt'"
	local tmpwgt
	
	if "`lowess'"!="" {
		local weight "weight"
		local line "line"
	}

	local sort "sort"	/* Force -sort- */

	tokenize `varlist'
	local Y `1'
	local X `2'

	tempvar grp
	if "`by'"~="" {
		sort `touse' `by' 
		qui by `touse' `by': gen long `grp'=1 if _n==1 & `touse'
		qui replace `grp'=sum(`grp') if `touse'
		local n = `grp'[_N]
					/* check more than one groups */
		if `n'<2 {
			di in red "by() variable takes on only one value"
			exit 198
		}
	}
	else {
		local n=1
		qui gen `grp'=1
	}
/*
	Set up complete.
*/
	local bw = `bwidth'
	if `bw' <= 0 | `bw' >= 1 { local bw .8 }

	local smlist
	local lsymbol
	local lconnect
	local i=1

	while `i'<=`n' {
		tempvar subuse
		qui gen byte `subuse'=`touse' & `grp'==`i'
		tempvar Y1 X1 smooth`i'
		quietly {
			gen `Y1'=`Y' if `subuse'
			gen `X1'=`X' if `subuse'
			summ `X1'
			local cnt = r(N)
			if `cnt'<3 { noisily error 2001 }
			if r(max)-r(min)<1e-30 {
				di in red "Range of `2' is too small"
				exit 499
			}
			_crcslbl `Y1' `1'
			_crcslbl `X1' `2'
/*
	smopt=0 (unwtd mean), 1(unwtd line), 2(wtd mean), 3(wtd line).
*/
			local smopt = ("`line'"=="line") /* 
				*/ + cond("`weight'"=="weight",2,0)
			local estpnts = `numpnts'
			if (`estpnts' <= 0 | `estpnts' >= `cnt') /*
			        */ { local estpnts `cnt' }

	Ksmapprox `Y1' `X1' `cnt' `bw' `smooth`i'' `smopt' `subuse' `estpnts'

			if "`logit'" == "logit" {
				local adj = 1/`cnt'
				local small = 0.0001
				replace `smooth`i'' = `adj' if `smooth`i'' /*
					*/ <`small' & `touse'
				replace `smooth`i'' = 1-`adj'  /*
					*/ if `smooth`i''>(1-`small') & `touse'
				replace `smooth`i'' = ln(`smooth`i''/(1- /*
					*/ `smooth`i'')) if `touse'
			}
		}

		local smlist `smlist' `smooth`i''
		local lsymbol "`lsymbol'i"
		local lconnect "`lconnect'l"
		
		capture drop `subuse' `Y1' `X1'
		local i=`i'+1
	}

/*
	Graph (if required)
*/
	if "`graph'" == "" { 
		if `"`t1title'"' ==""{
			if "`lowess'"!="" {
				local t1title "Lowess smoother"
			}
			else if "`line'" == "line" {
				local t1title "Running line smoother"
			}
			else 	local t1title "Running mean smoother"
			local t1title `"`t1title', bandwidth = `bw'"'
		}
		if "`logit'"=="" {
			tempvar sm
			qui egen `sm'=rmax(`smlist')
			if "`by'"~="" {sort `by'}

			local symbol = cond("`symbol'"=="", /*
				*/ "oi", "`symbol'")
			local connect = cond("`connect'"=="", /*
				*/ ".l", "`connect'")
			gr `Y' `sm' `X' if `touse', `options' /*
				*/ t1(`"`t1title'"') s(`symbol') /*
				*/ c(`connect') `sort' by(`by')
		}
		else {
			tempvar sm
			qui egen `sm'=rmax(`smlist')
			if "`by'"~="" {sort `by'}

			local symbol = cond("`symbol'"=="", /*
				*/ "i", "`symbol'")
			local connect = cond("`connect'"=="", /*
				*/ "l", "`connect'")
			gr `sm' `X' if `touse', `options' /*
				*/ t1(`"`t1title'"') s(`symbol') /*
				*/ c(`connect') `sort' by(`by')
		}
	}

	if "`gen'"~="" {
		qui egen `gen'=rmax(`smlist')
	}

end


program define Ksmapprox
	args Y X count bwidth Ys smopt touse estpnts
	tempvar W
	capture drop `Ys'
/*
	Args: 1=Y, 2=X, 3=count, 4=bandwidth, 5=smoothed Y,
	6=0 (unwtd mean), 1(unwtd line), 2(wtd mean), 3(wtd line).
	NOTE: bandwidth is expressed as a fraction of the sample size,
	not a number of sample values.
	7=number of estimation points
*/

	tempvar revuse 
	gen byte `revuse' = cond(`touse',1,2)
	sort `revuse' `X'
	drop `revuse'
	/* touse sample is now 1/`count' */

/*
	Using symmetric neighbourhoods.
	`k' is half-bandwidth. Initialise using points 1,...,k+1.
	Delta is distance between current Y and furthest neighbour in the
	interval. `W' is tricube weight function within the interval.
*/
	

	gen `W' = 1 in 1/`count'
	gen `Ys' = .
	local wtd = (`smopt'>1)
	local mean = (`smopt'==0 | `smopt'==2)
	local k = int((`bwidth'*`count'-0.5)/2)
	local i 0
	local j 0
	while `j'<`estpnts' {
		local j = `j' + 1
		local i = round(`count'*`j'/(`estpnts'+1),1)
		local x0 = `X'[`i']
		local imk = max(1, `i'-`k')
		local ipk = min(`i'+`k', `count')
		if `wtd' {
			local delta = 1.0001* /*
				*/ max(`X'[`ipk']-`x0', `x0'-`X'[`imk'])
			replace `W' = cond(`delta'>1e-30, /*
				*/ (1-(abs(`X'-`x0')/`delta')^3)^3 , 1) /*
				*/ in `imk'/`ipk'
		}
		if `mean' {
			summ `Y' [aw=`W'] in `imk'/`ipk'
			replace `Ys' = r(mean) in `i'
		}
		else {
			reg `Y' `X' [aw=`W'] in `imk'/`ipk'
			replace `Ys' = _b[_cons]+_b[`X']*`x0' in `i'
		}
	}
end