*! version 1  3sep1999 by Duolao Wang   (STB-52 sg123)

program define npshift, rclass
	version 6
	syntax varname [if] [in],BY(varname)[Level(int $S_level)]
	marksample touse
	markout `touse' `by'

	local x "`varlist'"
	tempname m n g1 g2 theta theta_l theta_u alpha
	tempfile grp1
	tempvar x1 x2 x2_1 

	preserve

	quietly {
		summarize `by' if `touse'
		if r(N) == 0 { noisily error 2000 }
		if r(min) == r(max) {
			di in red `"1 group found, 2 required"'
			exit 499
		}
		scalar `g1' = r(min)    
		scalar `g2' = r(max)    

		count if `by'!=`g1' & `by'!=`g2' & `touse'
		if r(N) != 0 {
			di in red `"more than 2 groups found, only 2 allowed"'
			exit 499
		}

		gen `x1'=`x' 
		keep if `by'==`g1' & `touse'
		keep `x1' 
		scalar `m'=_N
		save `grp1',replace

		restore,preserve
		gen `x2'=`x' 
		keep if `by'==`g2' & `touse'
		keep `x2'
		scalar `n'=_N

		cross using `grp1'
		gen `x2_1'=`x2'-`x1'

		sort `x2_1'
		centile `x2_1',centile(50)
		scalar `theta'=r(c_1)
}

		scalar `alpha'=(100-`level')/100
		scalar Calpha1=round(`m'*`n'/2-invnorm(1-`alpha'/2)*sqrt(`m'*`n'*(`m'+`n'+1)/12),1)
		scalar Calpha2=`m'*`n'-Calpha1+1

		scalar `theta_l'=`x2_1'[Calpha1]
		scalar `theta_u'=`x2_1'[Calpha2]

*formating

		local fmt : format `x1'

		if substr("`fmt'",-1,1)=="f" { 
			local ofmt="%9."+substr("`fmt'",-2,2)
			local lfmt="%-9."+substr("`fmt'",-2,2)
		}
		else if substr("`fmt'",-2,2)=="fc" {
			local ofmt = "%9." + substr("`fmt'",-3,3)
			local lfmt = "%-9." + substr("`fmt'",-3,3)
		}
		else	{local ofmt "%9.0g"
	 		 local lfmt "%-9.0g"
		}

#delimit ;
    	di in gr _n `"Hodges-Lehmann Estimates of Shift Parameters"' _n
         in gr _dup(65) `"-"'
	   in gr _n `"Point Estimate of Shift : Theta = Pop_2 - Pop_1 = "' _col(46) 
         in ye `lfmt' `theta'
	   in gr _n `"`level'% Confidence Interval for Theta: "' _col(44) 
         in ye `"["'  `lfmt' `theta_l' _col(52) 
         in re `","'  _col(54) 
         in ye `ofmt' `theta_u]' `"]"' _n
	   in gr _dup(65) `"-"';
#delimit cr


*save results

	return scalar  m      = `m'
	return scalar  n      = `n'
 	return scalar  theta    = `theta'
	return scalar  theta_l  = `theta_l' 
	return scalar  theta_u  = `theta_u'
 	return scalar  level    = `level'

* Double saves 

	global S_1 "`return(m)'"
	global S_2 "`return(n)'"
	global S_3 "`return(theta)'"
	global S_4 "`return(theta_u)'"
	global S_5 "`return(theta_l)'"
	global S_6 "`return(level)'"

end