Statalist


[Date Prev][Date Next][Thread Prev][Thread Next][Date Index][Thread Index]

st: RE: How to convert a do-file (that uses ml functions) into an ado file


From   "Nick Cox" <n.j.cox@durham.ac.uk>
To   <statalist@hsphsun2.harvard.edu>
Subject   st: RE: How to convert a do-file (that uses ml functions) into an ado file
Date   Thu, 18 Dec 2008 19:20:10 -0000

The overarching tip here is to learn about -set trace on- so that you
can see for yourself where Stata exits. 

David Elliott spotted one bug. Here are some more: 

1. The reference to `corr^2' will not work. You mean `corr'^2 -- in fact
you should really mean (`corr')^2 which will work as desired if `corr'
is ever negative. 

2. Your -syntax- statement allows -if- and -in- but you ignore any such
information within the program. The implication is that the program will
accept -if- and -in- conditions but ignore them, i.e. apply the program
to as many as observations as are non-missing on all the variables. 

3. Your -syntax- statement specifies -byable(recall)- but any attempt to
call your program under -by:- will fail. See the documentation for why. 

Here are some inefficiencies, minor but nevertheless not best style. 

5. You reciprocate by using powering with a power of -1. I don't think
there is any intelligence in Stata that catches that as a special case.
In fact Stata will call up a general powering routine and pass it -1 as
a power. It is faster and better just to divide. Thus foo^-1 is better
as 1/foo. 

Incidentally, I would always parenthesise that power, so foo^(-1) not
foo^-1; unless you've memorised Stata's precedence rules, it's quicker
to add the parentheses than to look up the rules, and the code is that
bit clearer to any reader who knows less Stata than you do. 

6. At various points you -summarize- when you only want the sum. In that
circumstance, you should use the -meanonly- option. This is documented
in 

SJ-7-3  st0135  . . . . . . . . . . . Stata tip 50: Efficient use of
summarize
        . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .  N.
J. Cox
        Q3/07   SJ 7(3):438--439                                 (no
commands)
        tip on using the meanonly option of summarize to more
        quickly determine the mean and other statistics of a
        large dataset

I think it's best to define programs separately, not nested. 

There may be other bugs, but that's all I saw. 

Nick 
n.j.cox@durham.ac.uk 

Tiago V. Pereira

I would like to convert a do-file into a ado file, but could not figure
out why my implementation outputs the following error:

. one_more_test a b c d e f
'b2' found where number expected
(error occurred while loading one_more_test.ado)
r(7);

Extremely grateful for any tip!

*/---------------- DATASET--START----------------------------
input a b c d e f
12 50 50 25 68 98
458 55 52 53 63 52
125 65 5 258 55 4
145 88 7 88 55 44
end
*/---------------- DATASET--END----------------------------

*/--------------- PROGRAM--START---------------------------

capture program drop one_more_test
program define one_more_test, rclass byable(recall)
version 8.0
syntax varlist(min=6 max=6 default=none numeric) [if] [in][,]
tokenize `varlist'


	tempvar b1 b2 v11 v22 v12 _wz1 _wz2 _z1wz1 _z2wz2

	qui gene `b1'=ln((`5'/ `4')/(`2'/`1'))
	qui gene `b2'=ln((`6'/ `4')/(`3'/`1'))
	qui gene `v11'=(`1'^-1) +(`2'^-1)+(`4'^-1) +(`5'^-1)
	qui gene `v22'=(`1'^-1) +(`4'^-1)+(`3'^-1) +(`6'^-1)
	qui gene `v12'=(`1'^-1) +(`4'^-1)
	qui gene `_wz1'=`v11'^-1
	qui gene `_wz2'=`v22'^-1
	qui gene `_z1wz1'=`b1'*`_wz1'
	qui gene `_z2wz2'=`b2'*`_wz2'


	summarize `_z1wz1'
	local muz1  = r(sum)

	summarize `_wz1'
	local swz1  = r(sum)
	local muz1  = `muz1'/`swz1'
	local tauz1 =`swz1'^-1

	summarize `_z2wz2'
	local muz2  = r(sum)
	summarize `_wz2'
	local swz2  = r(sum)
	local muz2  = `muz2'/`swz2'
	local tauz2 = `swz2'^-1
	local lambda= `muz1'/`muz2'
	local tauz2= ln(`tauz2')

	program LL
	args lnl lambda mub2 logtb2

	quietly {
		tempvar vb2 vb1 corr cov L
		scalar taub2=exp(`logtb2')
		scalar mub1=`mub2'*`lambda'
		gen double `vb1'=(`v11'+`lambda'^2*taub2)
		gen double `vb2'=( taub2 +`v22')
		gen double `cov'=(`lambda'*taub2 +`v12')
		gen double `corr'= `cov'/sqrt(`vb2'*`vb1')
		gen double `L'= ln(`vb1'*`vb2'*(1-`corr^2')) + /*
					*/ ( (`b1' - mub1 )^2/`vb1' + /*
					*/ (`b2' -`mub2')^2/`vb2' - /*
					*/ 2*`corr'*(b1 - mub1)*(`b2' -
`mub2')/sqrt(`vb2'*`vb1') )/* */
/(1-`corr'^2)
		replace `lnl'=-0.5*`L' - 0.5*ln(2*_pi)

		}

	end

	ml model lf LL (lambda:) (b2:) (logvar2:)
	ml init lambda:_cons=`lambda' b2:_cons=`muz2'
logvar2:_cons=`tauz2' ml check
	ml search
	ml maximize, difficult iter(30)

end

*/--------------- PROGRAM--END---------------------------

*/---------------- ORIGIANL DO-FILE--START---------------

gen b1=log(( e/ d)/( b/a))
gen b2=log(( f/ d)/( c/a))
gen V11=1/a +1/b+1/d +1/e
gen V22=1/a +1/d+1/c +1/f
gen V12=1/a +1/d
gen _wz1=1/V11
gen _wz2=1/V22
gen _z1wz1=b1*_wz1
gen _z2wz2=b2*_wz2
summarize _z1wz1
local muz1=r(sum)
summarize _wz1
local swz1=r(sum)
local muz1=`muz1'/`swz1'
local tauz1=1/`swz1'
summarize _z2wz2
local muz2=r(sum)
summarize _wz2
local swz2=r(sum)
local muz2=`muz2'/`swz2'
local tauz2=1/`swz2'
local lambda=`muz1'/`muz2'
drop _wz1 _wz2 _z1wz1 _z2wz2
local tauz2=log(`tauz2')
program drop _all

program define LL
args lnl lambda mub2 logtb2
quietly {
scalar taub2=exp(`logtb2')
scalar mub1=`mub2'*`lambda'
gen double vb1=`lambda'^2*taub2 + V11
gen double vb2=taub2 + V22
gen double cov=`lambda'*taub2 + V12
gen double corr=cov/sqrt(vb2*vb1)
gen double L=log(vb1*vb2*(1-corr^2)) + /*
*/ ( (b1 - mub1 )^2/vb1 + /*
*/ (b2 -`mub2')^2/vb2 - /*
*/ 2*corr*(b1 - mub1)*(b2 - `mub2')/sqrt(vb2*vb1) )/*
*/ /(1-corr^2)
replace `lnl'=-0.5*L - 0.5*log(2*_pi)
drop vb2 vb1 corr cov L
}
end
ml model lf LL (lambda:) (b2:) (logvar2:)
ml init lambda:_cons=`lambda' b2:_cons=`muz2' logvar2:_cons=`tauz2' ml
check
ml search
ml maximize, difficult iter(30)
*/------------------------ ORIGINAL--DO-FILE--END-------


*
*   For searches and help try:
*   http://www.stata.com/help.cgi?search
*   http://www.stata.com/support/statalist/faq
*   http://www.ats.ucla.edu/stat/stata/



© Copyright 1996–2014 StataCorp LP   |   Terms of use   |   Privacy   |   Contact us   |   What's new   |   Site index