**********************************************************************

Attached below are the following ado files:

	biprob.ado
	suprob.ado
	bipr_l2.ado
	bipr_lf.ado

And help files:

	biprob.hlp

**********************************************************************
biprob.ado

Instructions for use:  Cut out the following program and put it in a
file called "biprob.ado".  Be sure that the file is saved as a plain
text (ASCII) file and that it has a hard return at the end of the
last line in the file. 

Put biprob.ado in your C:\ADO directory (you may have to create this
directory) or your Stata working directory (i.e., your current directory).

Note:  Do not put it in the C:\STATA\ADO directory.  Only official Stata
ado and hlp files should be placed in this directory.

For more information, see [2] ado in Vol. 1 of the Reference manuals.

---------------------------------cut here---start-of-biprob.ado-------
*! version 1.6.0  28aug1997
program define biprob
	version 5.0
        local options "Level(integer $S_level)"
        if (substr("`1'",1,1)=="," | "`*'"=="") {
                if "$S_E_cmd"~="biprob" {
                        error 301
                }
                parse "`*'"
                if "`level'" != "" {
                        if `level' < 10 | `level' > 99 {
                                local level = 95
                        }
                }
        }
        else {
		local varlist   "req ex min(2)"
		local if	"opt"
		local in	"opt"
		local options	"ROBust SCore(string) TRace noLOg"
		local options	"`options' Level(integer $S_level)"
		local options	"`options' noCHI CLuster(string) *"
		parse "`*'"

		if "`score'" != "" {
			if "`robust'" == "" {
				local robust "robust"
			}
			local n : word count `score'
			if `n' != 2 {
				noi di in red "need 2 arguments in score()"
				exit 198
			}
			local sc1nam : word 1 of `score'
			local sc2nam : word 2 of `score'
			confirm new var `sc1nam'
			confirm new var `sc2nam'
		}
		parse "`varlist'", parse(" ")
		local dep1 "`1'"
		local dep2 "`2'"
		mac shift 2
		local ind1 "`*'"
		local ind2 "`*'"
		eq `dep1' : `dep1' `ind1'
		eq `dep2' : `dep2' `ind2'
		local eq1 "`dep1'"
		local eq2 "`dep2'"
		
		quietly {
			tempvar touse
			mark `touse' `if' `in' 
			markout `touse' `dep1' `ind1' `dep2' 

			if "`cluster'" != "" {
				unabbrev `cluster'
				if $S_2 > 1 {
					noi di in red /*
				*/ "only one varname allowed in cluster()"
					exit 198
				}
				local cluster "$S_1"
				local clarg ", cluster(`cluster')"
				local robust "robust"
				confirm var `cluster'
				markout `touse' `cluster', strok
			}

			global SC_1
			global SC_2
			global SC_3
			global S_y1 "`dep1'"
			global S_y2 "`dep2'"
			global S_X1 "`ind1'"
			global S_X2 "`ind2'"

			if "`chi'" == "" {
				/* Fit constant only model */
				probit `dep1' if `touse'
				tempname b0 b1 b2 b3 eq01 eq02
				capture mat `b1' = get(_b)
			        if _rc {
					count if `touse' `mif'
					local n1 = _result(1)
					count if `dep1'~=0
						*/ & `touse' `mif'
					local n2 = _result(1)
					noi BadProb `n1' /*
						*/ `n2' `b1' 0 `dep1'
				}
				mat coleq `b1' = `eq01'
				probit `dep2' if `touse'
				capture mat `b2' = get(_b)
				if _rc {
					count if `touse'
					local n1 = _result(1)
					count if `dep2'~=0 & `touse'
					local n2 = _result(1)
					noi BadProb `n1' /*
						*/ `n2' `b2' 0 `dep2'
				}
				mat coleq `b2' = `eq02'
				mat `b3' = (0)
				mat colnames `b3' = athrho:_cons
				mat `b0' = `b1',`b2'
				mat `b0' = `b0',`b3'
				eq athrho:
				eq `eq01' : `dep1'
				eq `eq02' : `dep2'

				ml begin
				tempvar mltouse 
				ml sample `mltouse' if `touse', noauto
				ml function bipr_l2
				ml method deriv2
				ml model b = `eq01' `eq02' athrho, /*
					*/ depv(110) from(`b0')
				if "`log'" == "nolog" {
					local prel "qui"
				}
				else    local prel "noi"
				`prel' di 
				`prel' di in gr "Fitting constant only model"
				`prel' ml max f V, `trace' `options' 
				local ll0 = scalar(f)
				eq drop `eq01'
				eq drop `eq02'
				local llarg "lf0(`ll0') pr2"
			}
			else {
				local llarg ""
			}

			/* Fit full model */
		
			probit `dep1' `ind1' if `touse'
			tempname b0 b1 b2 b3
			capture mat `b1' = get(_b)
			if _rc {
				count if `touse' `mif'
				local n1 = _result(1)
				count if `dep1'~=0 & `touse' `mif'
				local n2 = _result(1)
				noi BadProb `n1' `n2' `b1' 1 `dep1'
			}
			mat coleq `b1' = `dep1'
			probit `dep2' `ind2'  if `touse'
			capture mat `b2' = get(_b)
			if _rc {
				count if `touse'
				local n1 = _result(1)
				count if `dep2'~=0 & `touse'
				local n2 = _result(1)
				noi BadProb `n1' `n2' `b1' 1 `dep2'
			}
			mat coleq `b2' = `dep2'
			mat `b3' = (0)
			mat colnames `b3' = athrho:_cons
			mat `b0' = `b1',`b2'
			mat `b0' = `b0',`b3'
			eq athrho:

			if "`robust'" != "" {
				tempvar s1 s2 s3
				gen double `s1' = .
				gen double `s2' = .
				gen double `s3' = .
				global SC_1 "`s1'"
				global SC_2 "`s2'"
				global SC_3 "`s3'"
			}
			ml begin
			tempvar mltouse 
			ml sample `mltouse' if `touse', noauto
			ml function bipr_l2
			ml method deriv2
			ml model b = `eq1' `eq2' athrho, depv(110) from(`b0')
			if "`log'" == "nolog" {
				local prel "qui"
			}
			else    local prel "noi"
			`prel' di 
			`prel' di in green "Fitting full model"
			`prel' ml max f V, `trace' `options'
			local title "Bivariate probit regression"
			ml post biprob, `llarg' title(`title')
			if "`robust'" != "" {
				_robust `s1' `s2' `s3' `clarg'
				if "`score'" != "" {
					rename `s1' `sc1nam'
					rename `s2' `sc2nam'
				}
			}
		}
		if "`cluster'" != "" {
			global S_E_cvn "`cluster'"
		}
		global S_E_cmd "biprob"
	}
	ml mlout biprob, level(`level') 
	birpt `level'
        di in gr "     rho |  " in ye %9.5f $S_1 _c
	di in ye _col(3) %9.5f $S_2 _col(37) %9.5f $S_3 /*
                */ %12.5f $S_4
	di in gr _dup(78) "-"



end

program define birpt
	tempname arho serho z t1 t2 y
	scalar `z' = invnorm((`1'+100)/200)
	scalar `arho' = [athrho]_cons
	scalar `serho' = [athrho]_se[_cons]
	zap_s
	global S_1 = (exp(2*`arho')-1)/(exp(2*`arho')+1)
        global S_2 = (1 - ($S_1)^2)*`serho'
	if $S_2 ==0 {
		global S_2 = .
		global S_3 = .
		global S_4 = .
		exit
	}
	scalar `t1' = [athrho]_cons-`z'*[athrho]_se[_cons]
        scalar `t2' = [athrho]_cons+`z'*[athrho]_se[_cons]
        scalar `y' = exp(`t1')
        global S_3 = (`y' - 1/`y')/(`y' + 1/`y')
        scalar `y' = exp(`t2')
        global S_4 = (`y' - 1/`y')/(`y' + 1/`y')
end


program define BadProb
	local ntot = `1'
	local none = `2'
	local mat "`3'"
	local mod = `4'
	local depv "`5'"

	
	if `mod' {
		local trail "in full model for `depv'"
	}
	else    local trail "in constant only model for `depv'"

	if `ntot' == 0 {
		noi di in red "error in obtaining starting values `trail'"
		noi di in red "no observations"
		exit 2000
	}

	if `ntot' == `none' | `none' == 0 {
		noi di in red "error in obtaining starting values `trail'"
		noi di in red "outcome does not vary"
		exit 198
	}
	noi di in red "error in obtaining starting values `trail'"
	exit 198
end
---------------------------------cut here---end-of-biprob.ado-------
**********************************************************************
suprob.ado

Instructions for use:  Cut out the following program and put it in a
file called "suprob.ado".  Be sure that the file is saved as a plain
text (ASCII) file and that it has a hard return at the end of the
last line in the file. 

Put suprob.ado in your C:\ADO directory (you may have to create this
directory) or your Stata working directory (i.e., your current directory).

Note:  Do not put it in the C:\STATA\ADO directory.  Only official Stata
ado and hlp files should be placed in this directory.

For more information, see [2] ado in Vol. 1 of the Reference manuals.

---------------------------------cut here---start-of-suprob.ado-------
*! version 1.6.0  28aug1997
program define suprob
	version 5.0
        local options "Level(integer $S_level)"
        if (substr("`1'",1,1)=="," | "`*'"=="") {
                if "$S_E_cmd"~="suprob" {
                        error 301
                }
                parse "`*'"
                if "`level'" != "" {
                        if `level' < 10 | `level' > 99 {
                                local level = 95
                        }
                }
        }
        else {
		parse "`*'", parse(",")
		local eqns	"`1'"
		mac shift 
		local if	"opt"
		local in	"opt"
		local options	"ROBust SCore(string) TRace noLOg"
		local options	"`options' Level(integer $S_level)"
		local options	"`options' noCHI CLuster(string) NESTED *"
		parse "`*'"

		if "`score'" != "" {
			if "`robust'" == "" {
				local robust "robust"
			}
			local n : word count `score'
			if `n' != 2 {
				noi di in red "need 2 arguments in score()"
				exit 198
			}
			local sc1nam : word 1 of `score'
			local sc2nam : word 2 of `score'
			confirm new var `sc1nam'
			confirm new var `sc2nam'
		}

		parse "`eqns'", parse(" ")
		local eq1	"`1'"
		local eq2	"`2'"

		eq ? `eq1'
		local dep1 "$S_3"
		parse "$S_1", parse(" ")
		mac shift
		local ind1 "`*'"

		eq ? `eq2'
		local dep2 "$S_3"
		parse "$S_1", parse(" ")
		mac shift
		local ind2 "`*'"

		tempvar touse
		mark `touse' `if' `in' 
		
		if "`cluster'" != "" {
			unabbrev `cluster'
			if $S_2 > 1 {
				noi di in red /*
				*/ "only one varname allowed in cluster()"
				exit 198
			}
			local cluster "$S_1"
			local clarg ", cluster(`cluster')"
			local robust "robust"
			markout `touse' `cluster', strok
		}

		quietly {
			global SC_4 "suprob"
			local fcn "bipr_l2"
			local mth "deriv2"
			local tle "Seemingly unrelated probit regression"
				
			if "`nested'" != "" {
				if "`robust'" != "" {
					noi di in red /*
*/ "robust, cluster(), and score() options not available with nested model"
					exit 198
				}
				summ `dep1' `if' `in', meanonly
				local n1 = _result(1)
				summ `dep2' `if' `in', meanonly
				local n2 = _result(1)

				if `n1' > `n2' {
					local tmp "`dep1'"
					local dep1 "`dep2'"
					local dep2 "`tmp'"
					local tmp "`ind1'"
					local ind1 "`ind2'"
					local ind2 "`tmp'"
				}
				tempvar marker
				gen `marker' = 1 if `dep1'==.
				local mif "& `marker'~=1" 
				replace `dep1' = 0 if `dep1'==.
				global SC_4 "nested"
				local fcn "bipr_lf"
				local mth "lf"
				local tle "Nested probit regression"
			}
			else {
				global SC_4 "suprob"
				local fcn "bipr_l2"
				local mth "deriv2"
				local tle /*
				*/ "Seemingly unrelated probit regression"
			}
			markout `touse' `dep1' `ind1' `dep2' `ind2'

			capture {

				global SC_1
				global SC_2
				global SC_3
				global S_y1 "`dep1'"
				global S_y2 "`dep2'"
				global S_X1 "`ind1'"
				global S_X2 "`ind2'"

				if "`chi'" == "" {

					/* Fit constant only model */

					probit `dep1' if `touse' `mif'
					tempname b0 b1 b2 b3 eq01 eq02
					capture mat `b1' = get(_b)
					if _rc {
						count if `touse' `mif'
						local n1 = _result(1)
						count if `dep1'~=0 
							*/ & `touse' `mif'
						local n2 = _result(1)
						noi BadProb `n1' /*
							*/ `n2' `b1' 0 `dep1'
					}
					mat coleq `b1' = `eq01'
					probit `dep2' if `touse' 
					capture mat `b2' = get(_b)
					if _rc {
						count if `touse' 
						local n1 = _result(1)
						count if `dep2'~=0 & `touse' 
						local n2 = _result(1)
						noi BadProb `n1' /*
							*/ `n2' `b2' 0 `dep2'
					}
					mat coleq `b2' = `eq02'
					mat `b3' = (0)
					mat colnames `b3' = athrho:_cons
					mat `b0' = `b1',`b2'
					mat `b0' = `b0',`b3'
					eq athrho:
					eq `eq01' : `dep1'
					eq `eq02' : `dep2'

					ml begin
					tempvar touse 
					ml sample `touse' `dep1' `dep2' /*
					*/ `if' `in'
					ml function `fcn'
					ml method `mth'
					ml model b = `eq01' `eq02' athrho, /*
						*/ depv(110) from(`b0')
					if "`log'" == "nolog" {
						local prel "qui"
					}
					else    local prel "noi"
					`prel' di 
					`prel' di in gr /*
					*/ "Fitting constant only model"
					`prel' ml max f V, `trace' `options' 
					local ll0 = scalar(f)
					eq drop `eq01'
					eq drop `eq02'
					local llarg "lf0(`ll0') pr2"
				}
				else {
					local llarg ""
				}
				/* Fit full model */
			
				probit `dep1' `ind1' if `touse' `mif'
				tempname b0 b1 b2 b3
				capture mat `b1' = get(_b)
				if _rc {
					count if `touse' `mif'
					local n1 = _result(1)
					count if `dep1'~=0 & `touse' `mif'
					local n2 = _result(1)
					noi BadProb `n1' `n2' `b1' 1 `dep1'
				}
				mat coleq `b1' = `dep1'
				probit `dep2' `ind2'  if `touse' 
				capture mat `b2' = get(_b)
				if _rc {
					count if `touse' 
					local n1 = _result(1)
					count if `dep2'~=0 & `touse' 
					local n2 = _result(1)
					noi BadProb `n1' `n2' `b1' 1 `dep2'
				}
				mat coleq `b2' = `dep2'
				mat `b3' = (0)
				mat colnames `b3' = athrho:_cons
				mat `b0' = `b1',`b2'
				mat `b0' = `b0',`b3'
				eq athrho:

				if "`robust'" != "" {
					tempvar s1 s2 s3
					gen double `s1' = .
					gen double `s2' = .
					gen double `s3' = .
					global SC_1 "`s1'"
					global SC_2 "`s2'"
					global SC_3 "`s3'"
				}
				ml begin
				tempvar touse 
				ml sample `touse' `dep1' `dep2' `ind1' /*
				*/ `ind2' `if' `in'
				ml function `fcn'  
				ml method `mth'     
			noi ml model b = `eq1' `eq2' athrho, depv(110) /*
				*/ from(`b0')
				if "`log'" == "nolog" {
					local prel "qui"
				}
				else    local prel "noi"
				`prel' di 
				`prel' di in green "Fitting full model"
				`prel' ml max f V, `trace' `options'
				ml post suprob, `llarg' title(`tle')
				if "`robust'" != "" {
					_robust `s1' `s2' `s3' `clarg'
					if "`score'" != "" {
						rename `s1' `sc1nam'
						rename `s2' `sc2nam'
					}
				}
			}
			local rc = _rc
                        if "$SC_4" == "nested" {
                                replace `dep1' = . if `marker'==1
                        }
			/*
                        eq drop `eq01'
                        eq drop `eq02'
			*/
			if `rc' {
				exit `rc'
			}
		}
		if "`cluster'" != "" {
			global S_E_cvn "`cluster'"
		}
		global S_E_cmd "suprob"
	}
	di
	ml mlout suprob, level(`level') 
        surpt `level'
        di in gr "     rho |  " in ye %9.5f $S_1 _c
        di in ye _col(3) %9.5f $S_2 _col(37) %9.5f $S_3 /*
                */ %12.5f $S_4
        di in gr _dup(78) "-"
end

program define surpt
        tempname arho serho z t1 t2 y
        scalar `z' = invnorm((`1'+100)/200)
        scalar `arho' = [athrho]_cons
        scalar `serho' = [athrho]_se[_cons]
	zap_s
        global S_1 = (exp(2*`arho')-1)/(exp(2*`arho')+1)
        global S_2 = (1 - ($S_1)^2)*`serho'
	if $S_2 == 0 {
		global S_2 = .
		global S_3 = .
		global S_4 = .
		exit
	}
        scalar `t1' = [athrho]_cons-`z'*[athrho]_se[_cons]
        scalar `t2' = [athrho]_cons+`z'*[athrho]_se[_cons]
	scalar `y' = exp(`t1')
        global S_3 = (`y' - 1/`y')/(`y' + 1/`y')
	scalar `y' = exp(`t2')
        global S_4 = (`y' - 1/`y')/(`y' + 1/`y')
end

program define BadProb
	local ntot = `1'
	local none = `2'
	local mat "`3'"
	local mod = `4'
	local depv "`5'"

	
	if `mod' {
		local trail "in full model for `depv'"
	}
	else    local trail "in constant only model for `depv'"

	if `ntot' == 0 {
		noi di in red "error in obtaining starting values `trail'"
		noi di in red "no observations"
		exit 2000
	}

	if `ntot' == `none' | `none' == 0 {
		noi di in red "error in obtaining starting values `trail'"
		noi di in red "outcome does not vary"
		exit 198
	}
	noi di in red "error in obtaining starting values `trail'"
	exit 198
end
---------------------------------cut here---end-of-suprob.ado-------
**********************************************************************
bipr_l2.ado

Instructions for use:  Cut out the following program and put it in a
file called "bipr_l2.ado".  Be sure that the file is saved as a plain
text (ASCII) file and that it has a hard return at the end of the
last line in the file. 

Put bipr_l2.ado in your C:\ADO directory (you may have to create this
directory) or your Stata working directory (i.e., your current directory).

Note:  Do not put it in the C:\STATA\ADO directory.  Only official Stata
ado and hlp files should be placed in this directory.

For more information, see [2] ado in Vol. 1 of the Reference manuals.

---------------------------------cut here---start-of-bipr_l2.ado-------
*! version 1.1.0  17jun1997
program define bipr_l2
	version 5.0
	local b "`1'"   /* beta */
	local f "`2'"   /* value of L */
	local g "`3'"   /* gradient */
	local h "`4'"   /* Hessian */

	mac shift 4
	local options "LASTIT FIRSTIT FAST(string)"
	parse "`*'"

	if "`lastit'"!="" { exit }

	tempvar qi1 qi2 wi1 wi2 gi1 gi2 rhoi phi2 t1 t2

	qui gen byte `qi1' = 2*($S_y1!=0) - 1 if $S_mlwgt==1
	qui gen byte `qi2' = 2*($S_y2!=0) - 1 if $S_mlwgt==1

	tempname b1 b2 rho
	matsplit `b' `b1' `b2' `rho'  /* Also defines S_L1 S_L2 */

	mat score double `wi1' = `b1' if $S_mlwgt==1
	mat score double `wi2' = `b2' if $S_mlwgt==1

	qui replace `wi1' = `qi1'*`wi1' if $S_mlwgt==1
	qui replace `wi2' = `qi2'*`wi2' if $S_mlwgt==1

	qui gen double `rhoi' = `qi1'*`qi2'*`rho' if $S_mlwgt==1

	/* Calculate log likelihood */

	qui gen double `phi2' = binorm(`wi1',`wi2',`rhoi') if $S_mlwgt==1
	qui gen double `t1' = log(`phi2') if $S_mlwgt==1
	qui summ `t1', meanonly
	scalar `f' = _result(18)
	if `fast' == 0  { 
		exit 
	}
	
	drop `t1'

	local nr = rowsof(matrix(`b'))
	local nc = colsof(matrix(`b'))
	matrix `g' = J(1,`nc',0)
	matrix `h' = J(`nc',`nc',0)

	qui gen double `gi1' = normd(`wi1')*/*
		*/ normprob((`wi2'-`rhoi'*`wi1')/sqrt(1-`rhoi'*`rhoi')) /*
		*/ if $S_mlwgt==1
	qui gen double `gi2' = normd(`wi2')*/*
		*/ normprob((`wi1'-`rhoi'*`wi2')/sqrt(1-`rhoi'*`rhoi')) /*
		*/ if $S_mlwgt==1

	/* Calculate gradient */

	qui gen double `t1' = `qi1'*`gi1'/`phi2' if $S_mlwgt==1
	if "$SC_1" != "" {
		qui replace $SC_1 = `t1' if $S_mlwgt==1
	}
	
	qui gen double `t2' = .
	local i 1
	while `i' <= $S_L1 {
		local xx : word `i' of $S_X1
		qui replace `t2' = `t1'*`xx' if $S_mlwgt == 1
		qui summ `t2', meanonly
		mat `g'[1,`i'] = _result(18)
		local i = `i'+1
	}
	qui summ `t1', meanonly
	mat `g'[1,`i'] = _result(18)
	local j = `i'+1

	qui replace `t1' = `qi2'*`gi2'/`phi2' if $S_mlwgt == 1
	if "$SC_2" != "" {
		qui replace $SC_2 = `t1' if $S_mlwgt
	}
	
	local i 1
	while `i' <= $S_L2 {
		local xx : word `i' of $S_X2
		qui replace `t2' = `t1'*`xx' if $S_mlwgt == 1
		qui summ `t2', meanonly
		mat `g'[1,`j'] = _result(18)
		local j = `j'+1
		local i = `i'+1
	}
	qui summ `t1', meanonly
	mat `g'[1,`j'] = _result(18)
	local j = `j'+1

	qui replace `t1' = exp(-.5*(`wi1'^2+`wi2'^2-2*`rhoi'*`wi1'*`wi2')/ /*
	*/ (1-`rhoi'^2))/(2*_pi*sqrt(1-`rhoi'^2)) if $S_mlwgt == 1

	qui replace `t2' = `qi1'*`qi2'*`t1'/`phi2' if $S_mlwgt == 1
	if "$SC_3" != "" {
		qui replace $SC_3 = `t2' if $S_mlwgt == 1
	}
	
	qui summ `t2', meanonly
	mat `g'[1,`j'] = _result(18)
	
	/* Calculate Hessian */

	tempvar di vi1 vi2

	qui gen double `di' = 1/sqrt(1-`rhoi'^2) if $S_mlwgt == 1
	qui gen double `vi1' =  `di'*(`wi2'-`rhoi'*`wi1') if $S_mlwgt == 1
	qui gen double `vi2' =  `di'*(`wi1'-`rhoi'*`wi2') if $S_mlwgt == 1

	qui replace `gi1' = normd(`wi1')*normprob(`vi1') if $S_mlwgt == 1
	qui replace `gi2' = normd(`wi2')*normprob(`vi2') if $S_mlwgt == 1

	qui replace `t1' = `di'*normd(`wi1')*normd(`vi1') if $S_mlwgt == 1

	local i 1
	while `i' <= $S_L1+1 {
		local xi : word `i' of $S_X1
		if "`xi'" == "" {
			local xi = 1 
		}
		
		local j = `i'
		while `j' <= $S_L1+1 {
			local xj : word `j' of $S_X1
			if "`xj'" == "" {
				local xj = 1 
			}
			qui replace `t2' = `xi'*`xj'*(-`wi1'*`gi1'/`phi2' - /*
				*/ `rhoi'*`t1'/`phi2' - `gi1'^2/`phi2'^2) /*
				*/ if $S_mlwgt == 1
			qui summ `t2', meanonly
			mat `h'[`i',`j'] = _result(18)
			mat `h'[`j',`i'] = _result(18)
			local j = `j'+1
		}

		local k 1
		while `k' <= $S_L2+1 {
			local xk : word `k' of $S_X2
			if "`xk'" == "" {
				local xk = 1 
			}
			qui replace `t2' = `qi1'*`qi2'*`xi'*`xk'* /*
				*/ (`t1'/`phi2' - `gi1'*`gi2'/`phi2'^2) /*
				*/ if $S_mlwgt == 1
			qui summ `t2', meanonly
			mat `h'[`i',`j'] = _result(18)
			mat `h'[`j',`i'] = _result(18)
			local k = `k' + 1
			local j = `j'+1
		}

		qui replace `t2' = `qi2'*`xi'*`t1'/`phi2'* /*
			*/ (`rhoi'*`di'*`vi1' - `wi1'-`gi1'/`phi2') /*
			*/ if $S_mlwgt == 1
		qui summ `t2', meanonly
		mat `h'[`i',`j'] = _result(18)
		mat `h'[`j',`i'] = _result(18)
		local i = `i'+1
	}

	local k = $S_L1+2
	local i 1
	while `i' <= $S_L2+1 {
		local xi : word `i' of $S_X2
		if "`xi'" == "" {
			local xi = 1 
		}
		local m = `k'
		local j = `i'
		while `j' <= $S_L2+1 {
			local xj : word `j' of $S_X2
			if "`xj'" == "" {
				local xj = 1 
			}
			qui replace `t2' = `xi'*`xj'*(-`wi2'*`gi2'/`phi2' - /*
				*/ `rhoi'*`t1'/`phi2' - `gi2'^2/`phi2'^2) /*
				*/ if $S_mlwgt == 1
			qui summ `t2', meanonly
			mat `h'[`k',`m'] = _result(18)
			mat `h'[`m',`k'] = _result(18)
			local m = `m'+1	
			local j = `j'+1
		}

		qui replace `t2' = `qi1'*`xi'*`t1'/`phi2'*/*
			*/ (`rhoi'*`di'*`vi2' - `wi2'-`gi2'/`phi2') /*
			*/ if $S_mlwgt == 1
		qui summ `t2', meanonly
		mat `h'[`k',`m'] = _result(18)
		mat `h'[`m',`k'] = _result(18)
		local i = `i'+1
		local k = `k'+1
	}

	local k = $S_L1+$S_L2+3
	qui replace `t2' = `t1'/`phi2'*(`di'^2*`rhoi'*(1-`di'^2*(`wi1'^2+/*
	*/ `wi2'^2-2*`rhoi'*`wi1'*`wi2'))+`di'*`wi1'*`wi2'-`t1'/`phi2') /*
	*/ if $S_mlwgt == 1
	qui summ `t2', meanonly
	mat `h'[`k',`k'] = _result(18)

	tempname minus
	scalar `minus' = -1
	mat `h' = `h'*`minus'
end

program define matsplit
	local mat "`1'"
	local mat1 "`2'"
	local mat2 "`3'"
	local mat3 "`4'"

	local j : coleq(`mat')
	local k : word count `j'
	local km1 = `k'-1
	parse "`j'", parse(" ")
	local i 2
	while "``i''" == "`1'" {
		local i = `i'+1
	}
	local im1 = `i'-1
	mat `mat1' = `mat'[1,1..`im1']
	mat `mat2' = `mat'[1,`i'..`km1']
	scalar `mat3' = `mat'[1,`k']
* Convert to atanh
	scalar `mat3' = (exp(2*`mat3')-1) / (exp(2*`mat3')+1)
	global S_L1 = `im1'-1
	global S_L2 = `km1'-`i'
end
---------------------------------cut here---end-of-bipr_l2.ado-------
**********************************************************************
bipr_lf.ado

Instructions for use:  Cut out the following program and put it in a
file called "bipr_lf.ado".  Be sure that the file is saved as a plain
text (ASCII) file and that it has a hard return at the end of the
last line in the file. 

Put bipr_lf.ado in your C:\ADO directory (you may have to create this
directory) or your Stata working directory (i.e., your current directory).

Note:  Do not put it in the C:\STATA\ADO directory.  Only official Stata
ado and hlp files should be placed in this directory.

For more information, see [2] ado in Vol. 1 of the Reference manuals.

---------------------------------cut here---start-of-bipr_lf.ado-------
*! version 1.1.0  17jun1997
program define bipr_lf
	version 5.0
	local ll  "`1'"  /* Log l */
	local xb1 "`2'"  /* x1*b1 */
	local xb2 "`3'"  /* x2*b2 */
	local rr  "`4'"  /* rho */

	
	local rr = (exp(2*`rr')-1) / (exp(2*`rr')+1)

	quietly {
		replace `ll' = binorm(`xb1',`xb2',`rr')   /*
			*/ if $S_y1==1 & $S_y2==1
		replace `ll' = binorm(-`xb1',`xb2',-`rr') /*
			*/ if $S_y1==0 & $S_y2==1
		replace `ll' = 1-normprob(`xb2')  if $S_y2==0
		replace `ll'= log(`ll')
	}
end
---------------------------------cut here---end-of-bipr_lf.ado-------
**********************************************************************
biprob.hlp

Instructions for use:  Cut out the following program and put it in a
file called "biprob.hlp"  Be sure that the file is saved as a plain
text (ASCII) file and that it has a hard return at the end of the
last line in the file. 

Put biprob.hlp in your C:\ADO directory (you may have to create this
directory) or your Stata working directory (i.e., your current directory).

Note:  Do not put it in the C:\STATA\ADO directory.  Only official Stata
ado and hlp files should be placed in this directory.

For more information, see [2] ado in Vol. 1 of the Reference manuals.

---------------------------------cut here---start-of-biprob.hlp------
.-
help for ^biprob^                                   (STB:       ^biprob^)
.-

Estimate bivariate probit models
- --------------------------------

     ^biprob^ depvar1 depvar2 [varlist] [^if^ exp] [^in^ range] [^,^ ^nochi^
                 ^l^evel^(^#^)^ ^sc^ore^(^score1 score2^)^ ^tr^ace ^rob^ust
                 ^cl^uster^(^clustervar^)^ ^nolo^g ]

^biprob^ shares the features of all estimation commands; see help @est@.

To reset problem-size limits, see help @matsize@.


Description
- -----------

^biprob^ estimates bivariate probit models where you have two dependent
binary variables and a common list of independent variables.  This is a 
convenient program as results can be obtained using the ^suprob^ command.


Options
- -------

^nochi^ specifies that a constant only model should not be fit.  For some 
    problems, the constant only model can take many iterations to converge.
    Skipping the fit of the constant only model means that there will be no
    chi-squared or pseudo r-squared values reported in the output. 

^level(^#^)^ specifies the confidence level, in percent, for the confidence 
    intervals of the coefficients; see help @level@.

^score(score1 score2)^ specifies that the scores for a robust estimate of 
    variance should be saved in the new variables ^score1^ and ^score2^.  The
    ^score1^ variable will contain the scores for the bivariate estimation and
    the ^score2^ variable will contain the scores for the correlation
    parameter equation.  This option implies specification of the ^score^
    option.

^trace^ specifies that the trace of the iterations should be presented.

^robust^ specifies that the sandwich estimate of variance should be calculated
    and the robust standard errors presented for the estimation.

^cluster(clustervar)^ specifies that the robust estimate of variance should be
    adjusted for clustering on the groups specified by the values of the 
    ^clustervar^ variable.

^nolog^ specifies that the iteration log should be suppresed.

Examples
- --------

	. ^biprob priv vote ln_inc ptaxes years^
	. ^biprob priv vote ln_inc ptaxes years, robust^


Also see
- --------

    STB:  ^suprob^, if installed

Author
- ------

    James W. Hardin
    Stata Corporation
    jwh@@stata.com
---------------------------------cut here---end-of-biprob.hlp------

<end>