*! version 1.3.1  02June98     STB-47 sbe25
program define permband
version 4.0
local varlist "req ex min(1) max(2)"
local if "opt"
local in "opt"
#delimit ;
local options "MOMent(int 1) ALpha(real 0.05) Symbol(str) Connect(str)
 noGRaph SEed(int 0) SPan(real 0) REps(int 50) T1title(str) *" ;
#delimit cr

* double and smooth not options
parse "`*'"
if `moment'>4 {
	di in red "Must specify moment (1, 2, 3 or 4)"
	exit 198
}
else local moment `moment'
if `alpha'<0 | `alpha'>1 {
	di in red "Significance level (alpha) must be between 0 and 100"
	exit 198
}
if "`symbol'"!="" {
	local sym symbol(`symbol')
}
else local sym symbol(idd)
if "`connect'"!="" {
	local con connect(`connect')
}
else local con connect(l..)
if `seed'>0 { set seed `seed' }
if `reps'<2 {
	di in red "invalid reps()"
	exit 198
}
if `reps'!=50 {
	di in blue "Cannot perform hypothesis test if reps!=50"
}
if `span'>0 {
	local span "span(`span')" 
	di in blue "Note: cannot perform hypothesis test unless" /*
          */ " use default span"
	local spspan=1
}
else {
	local span
	local spspan=0
}

parse "`varlist'", parse(" ")
local z `1'
local x `2'

tempvar touse X smz newz smnewz maxz minz zp miss
quietly {
	mark `touse' `if' `in'
	markout `touse' `z' `x'
	count if `touse'
	local cnt = _result(1)
	if `cnt'<5 { noisily error 2001 }
	if "`x'"!=""{
		gen  `X'=`x' if `touse'
		_crcslbl `X' `x'
	}
	else {
		gen  `X'=_n if `touse'
		lab var `X' "n"
	}
	gen `zp'=`z'^`moment' if `touse'
	running `zp' `X' if `touse', `span' gen(`smz') nograph double

/* Do permutations. */
	tempname big
	scalar `big'=9e29
	gen `maxz'=-`big' if `touse'
	gen `minz'=`big' if `touse'
	local i 1
	while `i' <= `reps' {
		sort `X'
		noi di in gr "." _cont
		permpr `newz'=`zp'
		running `newz' `X' if `touse', /*
		 */ `span' gen(`smnewz') nograph double
		replace `maxz'=max(`maxz',`smnewz') if `touse'
		replace `minz'=min(`minz',`smnewz') if `touse'
		drop `newz' `smnewz'
		local i = `i' + 1
	}
*	if "`smooth'"!="" {
		local sp=0.5*$S_2	/* half span actually used */
		running `maxz' `X' if `touse',span(`sp') nograph gen(`smnewz')
		replace `maxz'=`smnewz'
		drop `smnewz'
		running `minz' `X' if `touse',span(`sp') nograph gen(`smnewz')
		replace `minz'=`smnewz'
*	}
	local mp=0
	local sp=0
	count if `touse'
	local nobs=_result(1)
	if `moment'<4 {
/*
	Significance test for non-randomness
*/
		if `moment'==1 {
			local mp=-0.439+(0.0448*log(`nobs'))
			local sp=0.362+(-0.0306*log(`nobs'))
		}
		else if `moment'==2 {
			local mp=-0.350+(0.0349*log(`nobs'))
			local sp=0.284+(-0.0225*log(`nobs'))
		}
		else if `moment'==3 {
			local mp=-0.350+(0.0384*log(`nobs'))
			local sp=0.306+(-0.0240*log(`nobs'))
		}
	}
	local critv=`mp'+invnorm(1-`alpha')*`sp'
	if `critv'<0 { local critv=0 }
	count if `touse' & (`smz'<`minz' | `smz'>`maxz')
	local outof=_result(1)/`nobs'
	local p=normprob((`mp'-`outof')/`sp')
	if `outof'==0 { local p=1 }
	local y=substr("`1'",1,5)
	cap drop s`moment'_`y'
	cap drop l`moment'_`y'
	cap drop u`moment'_`y'
	rename `smz' s`moment'_`y'
	rename `maxz' u`moment'_`y'
	rename `minz' l`moment'_`y'
	if `moment'==1 { lab var s`moment'_`y' "Smooth of `z'" }
	else { lab var s`moment'_`y' "Smooth of `z'^`moment'" }
	lab var l`moment'_`y' "Lower boundary"
	lab var u`moment'_`y' "Upper boundary"
}
if (`nobs'<50 | `nobs'>1000) & `reps'==50  & `spspan'==0 { local star="*" }
else {local star=" "}
di _n _n in gr "  Variable |   Obs" _col(26) "Proportion " /*
 */ _col(45) "Critical value" _col(65) "P-value" _n /*
 */  "           |    " in blue "`star'" in gr _col(25) /*
 */ "outside band" _col(45) "(alpha = " `alpha' ")" _n /*
 */ " ----------+" _dup(59) "-"

if `moment'==1 { local lname=substr("`z'",1,8) }
else { 
	local tname=substr("`z'",1,6)
	local lname=substr("`tname'^`moment'",1,8)
}
if `reps'==50 & `spspan'==0 & `moment'<4 {
	local skip=10-length("`lname'") 
	#delimit ;
	di in gr _skip(`skip') "`lname' |" in ye
		%6.0f `nobs' "        "
		%6.0g `outof' "            "
		%6.0g `critv' "           "
		%8.3f `p' 
		;
		#delimit cr
}
else {  
	local skip=10-length("`lname'") 
	#delimit ;
	di in gr _skip(`skip') "`lname' |" in ye
		%6.0f `nobs' "        "
		%6.0g `outof' "            "
		;
	#delimit cr
}

if (`nobs'<50 | `nobs'>1000) & `reps'==50 & `spspan'==0 { noi di _n in blue /*
 */ "* - Extrapolation used for test when n<50 or n>1000"}

if "`graph'"!="nograph" {
	graph s`moment'_`y' l`moment'_`y' u`moment'_`y' `X', `sym' `con' /*
         */ sort `options' t1title("`t1title'")
}
end

*! version 1.0.1 10Jul98.
program define permpr
* added lines to ignore missing
        version 4.0
        local varlist "req new max(1)"
        local exp "req noprefix"
        local options "Seed(int 0)"
        tempvar new
        parse "`*'"
        rename `varlist' `new'
                                        /* identify the sample  */
        tempvar rhs u order
	if `seed'>0 { set seed `seed' }
	quietly {
		sort `exp'
		gen long `order'=_n
		gen `u'=uniform()
		replace `u'=. if `exp'==.
		sort `u'
		replace `new'=`exp'[`order']
	}
        label var `new' "permuted `exp'"
        rename `new' `varlist'
end
*! version 2.0.0 PS/PR 10Jun97.
*!! Based on Peter Sasieni's "running4.ado" emailed 09Jun97
program define running
version 4.0
local varlist "req ex min(1) max(2)"
local if "opt"
local in "opt"
local weight "aweight" 
#delimit ;
local options "noGraph GEN(string) Knn(str) Double Repeat(int 1) Mean LOGit
 SPan(real 0) TItle(string) Symbol(str) Connect(str) PEn(str) CI GENSe(str) 
 GENB(str) TWice *";
#delimit cr
parse "`*'"
if "`gen'"!="" { confirm new var `gen' }
if "`genb'"!="" { 
	if "`twice'`mean'`logit'"!="" {
		di in red "genb not available with mean, twice or logit"
		exit 198
	}
	confirm new var `genb' 
}
if "`gense'"!="" { 
	confirm new var `gense' 
	if `repeat'>1 | "`double'`logit'`twice'"!="" {
		di in red "gense not available with repeat>1, twice or logit"
		exit 198
	}
}
if "`ci'"!="" {
	if `repeat'>1 | "`double'`logit'`twice'"!="" {
		di in red "ci not available with repeat>1, twice or logit"
		exit 198
	}
}
if `span'!=0 & "`knn'"!="" {
	di in red "cannot specify both span and knn"
	exit 198
}
if `span'<0 | `span'>2-("`mean'"!="") {
	di in red "span must be between 0 and " 2-("`mean'"!="")
	exit 198
}
if "`knn'"=="" {
	local knn 0 /* previous default */
}
else {
	cap confirm var `knn'
	if _rc {
		confirm integer num `knn'
	}
	else {
		unabbrev `knn'
		local kvec $S_1
	}
}
if "`double'"!="" { local repeat=2*`repeat' }
local nrep `repeat'
if `repeat'>7 {
	local repeat 7
	noisily di in bl "[repeat set to 7]"
}
parse "`varlist'", parse(" ")
local y `1'
local x `2'
tempvar Y X smooth sy rsy touse kr kl
quietly {
	if "`weight'" != "" { 
		tempvar sw w wt
		gen  `w' `exp'
		replace `w'=. if `w'<=0
		local exp "[`weight'=`w']"
	}
	mark `touse' `if' `in'
	markout `touse' `y' `x' `w'
	count if `touse'
	local cnt = _result(1)
	local IN "in 1/`cnt'"
	if `cnt'<5 { noisily error 2001 }
	if `span'!=0 {
		local knn=(`cnt'*`span'-1)/2 }
	else {
		if `knn'==0 {local knn = .5*`cnt'^0.8 }
		local span=(2*`knn'+1)/`cnt'
	}
	if "`kvec'" !="" { local kk  int(`kvec'/sqrt(`repeat') +.5)}
	else {
		local kk = int(`knn' /sqrt(`repeat') +.5)
		if `kk'<=0 {
			noi di in red "Span too small. Increase span or knn."
			exit 2002
		}
	}
	if `knn'<=1 & "`ci'"!="" {nois di "[ci not produced when knn=1]"}
	sum `y' `exp' if `touse'
	local ycen=_result(3)
	gen  `Y'=`y'-`ycen' if `touse'
	_crcslbl `Y' `y'
	if "`x'"!=""{
		gen  `X'=`x' if `touse'
		_crcslbl `X' `x'
	}
	else {
		gen  `X'=_n if `touse'
		lab var `X' "n"
	}
	sort `X'
	if "`weight'" != "" { 
		gen double `sw'=.
		replace `w'=. if !`touse'
		sum `w'
		noi di in gr "(sum of wgt is " _result(1)*_result(3) ")"
		replace `w' = `w'/_result(3)
		gen `wt' = `w'
		local ks "(`sw'[_n+`kr']-cond(_n>`kl',`sw'[_n-`kl'],0))"
		local d3 "`sw'[_N]"
		local wX "`w'*"
		local wtX "`wt'*"
		local wn "`wt'`nk'*"
		local g_swby "by `X':replace `sw'=sum(`w') "   
		local g_wt   "by `X':replace `wt'=`sw'[_N]/_N  "
		local gen_swt "replace `sw'=sum(`wt') `IN'"
	}
	else {
		local ks "(`kr'+`kl')"
		local d3 "_N"
	}
	gen  `smooth'=`Y'
	lab var `smooth' "Smooth fit"
	gen double `sy'=cond(_n<`cnt',(`X'==`X'[_n-1] |`X'==`X'[_n+1]), /*
		*/           `X'==`X'[_n-1]) `IN'
	gen  `rsy'=sum(`sy') `IN'
	local ties=`rsy'[`cnt']
	local xcen=(`X'[1]+`X'[`cnt'])/2
	replace `X'=`X'-`xcen' `IN'
	sort `X'
	if `ties'>0 {
		tempvar Yt 
		`g_swby'	/* sw must be defined before d3 is used */
		`g_wt'
		by `X':replace `sy'=sum(`wX'`Y') if `touse'  
		by `X':replace `smooth'=`sy'[_N]/`d3' if `touse' 
		if "`twice'"!=""{ gen `Yt'=`smooth' }
	}
	else {
		local Yt "`Y'" 
	}
        `gen_swt'
	gen int `kl' = min(`kk'+1,_n)
	gen int `kr' = min(`kk',`cnt'-_n)
	if "`mean'"=="" {
		local mean "line"
		tempvar sxy rsxy alpha beta sx sxx rsx rsxx 
		gen `alpha'= .
		gen `beta'= . 
		gen double `sxy'=.
		gen  `rsxy'=.
		gen double `sx'=sum(`wtX'`X')
		gen double `sxx'=sum(`wtX'((`X')^2))
		diff `rsx' `sx' `kr' `kl' `ks' `cnt' gen
		diff `rsxx' `sxx' `kr' `kl' `ks' `cnt' gen
		drop `sx' `sxx'
	}
	while `repeat'>0 {
		replace `sy'=sum(`wtX'`smooth') `IN'
		diff `rsy' `sy' `kr' `kl' `ks' `cnt' replace
		if "`mean'"=="line" {
		     replace `sxy'=sum(`wtX'`X'*`smooth') `IN'
		     diff `rsxy' `sxy' `kr' `kl' `ks' `cnt' replace
		     replace `beta'=(`rsxy'-`rsx'*`rsy')/(`rsxx'-(`rsx')^2) `IN'
		     replace `alpha'=`rsy'-`beta'*`rsx' `IN'
		     replace `smooth'=cond((`rsxx'-`rsx'*`rsx') > max(0, /*
			*/ `rsxx'/1e8), `alpha'+`beta'*`X', `rsy') `IN'
		}
		else {
		     replace `smooth'=`rsy' `IN'
		}
		local repeat=`repeat'-1
	}
	if "`twice'"!="" {
	   tempvar res
	   gen `res' = `Yt'-`smooth'
	   local repeat `nrep'
	   while `repeat'>0 {
		replace `sy'=sum(`wtX'`res') `IN'
		diff `rsy' `sy' `kr' `kl' `ks' `cnt' replace
		if "`mean'"=="line" {
		     replace `sxy'=sum(`wtX'`X'*`res') `IN'
		     diff `rsxy' `sxy' `kr' `kl' `ks' `cnt' replace
		     replace `beta'=(`rsxy'-`rsx'*`rsy')/(`rsxx'-(`rsx')^2) `IN'
		     replace `alpha'=`rsy'-`beta'*`rsx' `IN'
		     replace `res'=cond((`rsxx'-`rsx'*`rsx') > max(0, /*
			*/ `rsxx'/1e8), `alpha'+`beta'*`X', `rsy') `IN' 
		}
		else {
		     replace `res'=`rsy' `IN'
		}
		local repeat=`repeat'-1
	   }
	   replace `smooth'=`smooth'+`res' `IN'
	}
	if "`mean'"=="line" {drop `sxy' `rsxy'}
	if `ties'>0 {
		by `X':replace `sy'=sum(`smooth')  if `touse'
		by `X':replace `smooth'=`sy'[_N]/_N   if `touse'
	}
	if "`ci'"!="" | "`gense'"!="" {
		tempvar se sigma2
		gen double `se'=(`Y'-`smooth')^2
		if `ties'>0 {
			`g_swby'
			by `X':replace `se'=sum(`wX'`se')  if `touse'
			by `X':replace `se'=`se'[_N]/`d3'  if `touse'
		}
		`gen_swt'
		replace `se'=sum(`wtX'`se') `IN'
*		Adjustment for estimating alpha and beta made in se not sigma2
*
		diff `sigma2' `se' `kr' `kl' `ks' `cnt' gen
		if "`mean'"=="line" {
			replace `se'=sqrt(`sigma2'* /*
			*/ cond((`rsxx'-`rsx'*`rsx') > max(0,`rsxx'/1e8), /*
			*/  ((`kl'+`kr')/(`kl'+`kr'-2)) /*
			*/  *(1+(`X'-`rsx')^2/(`rsxx'-(`rsx')^2))/`ks' , /*
			*/ 1/(`kl'+`kr'-1) ))  `IN'  
		}
		else {
			replace `se'=sqrt(`sigma2'/(`kl'+`kr'-1)) `IN'
		}
		if `ties'>0 {
			`g_swby'
			by `X':replace `se'=sum((`se')^2)   if `touse'
			by `X':replace `se'=sqrt(`se'[_N]/_N) if `touse'  
		}
	}
	replace `X'=`X'+`xcen'	 `IN'
	replace `Y'=`Y'+`ycen'	 `IN'
	replace `smooth'=`smooth'+`ycen' `IN'	
	if "`mean'"=="mean" &" `kvec'"=="" {
 		local k2=`cnt'+1-`kk'
                replace `smooth'=. in 1/`kk'
                replace `smooth'=. in `k2'/l
	}
	if "`logit'"!="" {
		sum `y' if `touse'
		if _result(5)<-5e-9 | _result(6)>1+5e-9 {
			noi di "[Yvar out of range [0,1]. logit not available]"
		} 
		else {
		sum `y' if `y'>0 & `y'<1 &`touse'
		local aL = min(1/(min(`span',1)*`cnt'+1),_result(5))
		local aU = max(1-1/(min(`span',1)*`cnt'+1),_result(6))
		replace `smooth'=log(cond(`smooth'< `aL', `aL'/(1-`aL'), /*
		 */  cond(`smooth'>`aU',`aU'/(1-`aU'),`smooth'/(1-`smooth')) /*
		 */	))  if `smooth' !=. `IN'
		local laL = log(`aL'/(1-`aL'))
		local laU = log(`aU'/(1-`aU'))
		local adj=(`laU'-`laL')/25
		replace `Y'=cond(`y'<`aL', `laL'-(1.3+uniform())*`adj', /*
		 */  	    cond(`y'>`aU', `laU'+(1.3+uniform())*`adj', /*
		 */    		log(`y'/(1-`y')) )) `IN'
		}
	}
* GRAPH (If required)
	if "`graph'"!="nograph" {
		if "`ci'"!="" {
			local df=`kl'+`kr'-("`mean'"=="line")
			local t=invt(`df',${S_level}/100)
			tempvar l u
			gen `l'=`smooth'-`t'*`se'
			gen `u'=`smooth'+`t'*`se'
			local pp "44"
		}
		if "`title'" ==""{
			local title "Running `mean' smoother"
		}
		if "`symbol'" ==""{ local symbol  "oiii" }
		if "`connect'"==""{ local connect ".lll" }
		if "`pen'"    ==""{ local pen     "23`pp'" }
		noisily graph `Y' `smooth' `l' `u' `X' `exp', `options' /*
		 */ ti(`title') s(`symbol') c(`connect') pen(`pen')
	}
	if "`gen'"  !="" {rename `smooth' `gen'}
	if "`gense'"!="" {rename `se' `gense'}
	if "`genb'" !="" {rename `beta' `genb'}
}
global S_1 `knn'
global S_2 `span'
end

prog def diff 
	local new `1'
	local old `2'
	local kr `3'
	local kl `4'
	local ks `5'
	local cnt `6'
	local com `7'
   	`com' `new' =(`old'[_n+`kr'] - cond(_n>`kl',`old'[_n-`kl'],0))/`ks'  /*
	*/	in 1/`cnt'
end