*! 8mar2011
program bivar
version 11.1 

mata: DoIt()
end

version 11.1
mata:
void DoIt() {
	
	syntaxl(st_local("0"),(&(varlist="anything"),&(ifin="[ifin]")), (&(hi="[hi()]"),&(decimal="[d:ecimal]"),&(overall="[o:verall]"), &(moren="[m:oren]"),&(store="[s:tore]"),&(out="[out()]"),&(labs="[l&v()]")))
	landv(labs,nlb="",vlb="")
	varlist=varlist(varlist,"","",vtypes="")
	if (length(varlist)<2) errel("At least two variables must be specified")
	origby=varlist[1]
	varlist=varlist[2..cols(varlist)]
	vtypes=vtypes[2..cols(vtypes)]
	if (st_isstrvar(origby)) {
		stata(sprintf("encode %s, g(%s) label(%s_bivar)",origby,byvar=st_tempname(),byvar))
		st_global(byvar+"[tmp_varsub]",origby)
		}
	else byvar=origby
	(void) st_addvar("byte",touse=st_tempname())
	stata(sprintf("qui replace %s=1 %s",touse,ifin))
	stata(sprintf("qui replace %s=0 if mi(%s,%s)",touse,touse,byvar))
	if (strlen(ifin)) stata(sprintf("qui count %s & mi(%s)",ifin,byvar))
	else stata(sprintf("qui count if mi(%s)",byvar))
	if (st_numscalar("r(N)")) printf("%f observations excluded for missing %s\n",st_numscalar("r(N)"),origby)
	
	decimal=decimal?"1":"0"
	hi=strtoreal(hi+"0")
	byfreqs=freq(byvar,touse)'
	byfreqs=select(byfreqs,!colmissing(byfreqs[1,]))
	body=J(4,3+3*(cols(byfreqs)+overall),"")
	cb=cols(body)
	vcols=rangel(4,cb,3)
	body[1,4]=vardisp(origby,nlb)
	body[4,2..3]=("p","N")
	body[3,vcols]=concat(valdisp(byvar,byfreqs[1,],vlb)," =","r")',(overall?"Overall":J(1,0,""))
	body[4,vcols]=strofreal((byfreqs[2,],(overall?sum(byfreqs[2,]):J(1,0,.))))
	
	table=table_createM()
	asarray(table,"head",4) //exp
	asarray(table,"stub",3) //exp
	table_code(table,"vline","t",1,"alt") 
	table_code(table,"vline","c",rangel(4,cb,3),"none")
	table_code(table,"vline","c",rangel(5,cb,3),"none")
	table_code(table,"vline","c",cb,"thick")
	if (!moren) table_code(table,"padL","c",rangel(4,cb,3),0)
	table_code(table,"padR","c",rangel(4,cb,3),0)
	table_code(table,"padR","c",rangel(5,cb,3),0)
	table_code(table,"align","c",rangel(6,cb,3),"<")
	table_code(table,"hline","r",4,"thick")
	table_code(table,"hline","cell",(3,1\3,2),"none")
	table_code(table,"span","cell",(1,4\2,4),cb-3)
	table_code(table,"span","cell",pairs((3\4),rangel(4,cb,3)),3)
	table_code(table,"class","c",rangel(4,cb,3),"weaker")
	//table_code(table,"class","r",rangel(rb+2,rb+rows(cat),2),"weaker")
	asarray(table,"spchars",(
	"-"," ",""\
	">","&#x2192;","&#x2192;"\
	"V","&#x21D3;","&#x2193;"))
	
	vtypes=vtype(varlist,vtypes,"vtype",store)
	for (i=1;i<=length(varlist);i++) {
		var=varlist[i]
		vtype=vtypes[i]
		isstr=st_isstrvar(var)
		if (isstr & substr(vtype,1,3)!="cat") vtype="cat"
		if (anyof(("cat","catc","catr"),vtype)) {
			vtype2=substr(vtype,strlen(vtype),1)
			vtype="cat"
			}
		if (substr(vtype,1,2)=="di") {
			vtype2=substr(vtype,3)
			vtype="di"
			}
		if (!anyof(("cont","di","cat"),vtype)) continue
		
		if (isstr) {
			str=var
			stata("encode "+var+", gen("+(var=st_tempname())+")")
			st_global(var+"[tmp_varsub]",str)
			}
		
		if (vtype=="cat") table_code(table,"hline","r",rb,"clear") 
		body=body\(isstr?char(1)+"hidata1"+char(2)+"S-"+char(3):"")+vardisp(var,nlb),J(1,cb-1,"")
		rb=rows(body)
		
		if (vtype=="cat") {
			sub=toindices(("r","c","t"):==vtype2)
			(void) st_addvar("byte",catmark=st_tempname())
			stata(sprintf("qui replace %s=%s & !mi(%s)",catmark,touse,var)) //excluding missing cat values
			freqs=freq((var,byvar),catmark)
			freqs=select(freqs,!rowmissing(freqs[,1..2]))
			st_dropvar(catmark)
			stub=uniqrows(freqs[,1])
			freqs=colshape(freqs[,3],rows(uniqrows(freqs[,2]))) //know 3 cols returned
			denom=*(&rowsum(freqs),&colsum(freqs),&sum(freqs))[sub]
			rown=strofreal(rowsum(freqs))
			if (overall) {
				if (vtype2=="c") denom=denom,sum(freqs)
				freqs=freqs,rowsum(freqs)
				}
			cat=strofreal(100:*freqs:/denom,"%3."+decimal+"f")
			ns=moren?strofreal(colshape(freqs,1)):J(length(cat),1,"")
			sign=("%"+char(1)+"1"+char(1)+"2","% "+char(1)+"3","%")[sub]
			cat=colshape((ns,colshape(cat,1),J(length(cat),1,sign)),cb-3)
			if (isstr) cathead=valdisp(var,stub,"v")
			else cathead=concat(valdisp(var,stub,vlb)," =","r")
			cat=cathead,J(rows(cat),1,""),rown,cat
			table_code(table,"class","cell",(rangel(rb+1,rb+rows(cat),1),J(rows(cat),1,1)),"it")
			body=body\cat
			table_code(table,"hline","r",rb+rows(cat),"clear")
			}
		else {
			if (vtype=="cont") {
				for(c=1;c<=cols(byfreqs)+overall;c++) { 
					obs= c>cols(byfreqs)?.:obswhere(byvar+"=="+strofreal(byfreqs[1,c]))
					st_view(v,obs,var,touse)
					mv=meanvariance(v)
					body[rb,1+3*c]=moren?strofreal(rows(v)):""
					body[rb,2+3*c]=strofreal(mv[1],st_varformat(var)) //override var formats?
					body[rb,3+3*c]=char(177)+strofreal(sqrt(mv[2]),"%9.3g") //sd won't be meaningful for date/times...
					}
				}
			else { //di
				(void) st_addvar("byte",var=st_tempname()) //change var so it's used for chi2
				if (vtype2=="") vtype2="#V!=0 & #V<."
				stata("qui replace "+var+"="+subinstr(vtype2,"#V",varlist[i]))
				freqs=freq((var,byvar),touse)
				freqs=select(freqs,!rowmissing(freqs[,2]))
				prod=freqs[,1]:*freqs[,3]
				prod=colshape(freqs[,3],rows(uniqrows(freqs[,2])))
				if (overall) prod=prod,rowsum(prod)
	//			body[rb,vcols]=strofreal(prod[rows(prod),])
				body[rb,vcols]=moren?strofreal(prod[rows(prod),]):J(1,length(vcols),"")
				body[rb,vcols:+1]=strofreal(100*prod[rows(prod),]:/colsum(prod),"%3."+decimal+"f") 
				body[rb,vcols:+2]=J(1,length(vcols),"%") 
				}
			}
		
		if (vtype=="cont") {
			stata("qui oneway "+var+" "+byvar+" if "+touse)
			if (!length(st_numscalar("r(N)"))) body[rb,2..3]=(".",".")
			else body[rb,2..3]=strofreal(Ftail(st_numscalar("r(df_m)"),st_numscalar("r(df_r)"),st_numscalar("r(F)")),"%4.3f") ,strofreal(st_numscalar("r(N)"))
			}
		else {
			stata("qui ta "+var+" "+ byvar+" if "+touse+", chi2")
			body[rb,2..3]=strofreal(empty(st_numscalar("r(p)"),"scalar"),"%4.3f"), strofreal(empty(st_numscalar("r(N)"),"scalar"))
			}
		if (hi & strtoreal(body[rb,2])<=hi) table_code(table,"class","cell",(rb,2),"hidata1")
		
		}
	table_code(table,"hline","r",rows(body),"thick")
	asarray(table,"body",body)
	printf("\n")
	table_present(table,out)
	}	

end