Dear all,
1. Thank you very much all of you who answered my r(2000) question.
2. Observe, however, that I didn't ask for help with keeping the loop
running. I knew how to deal with that -- as you can see in the code
attached, which I am sending in full this time, and which I didn't send
in full before because it didn't seem relevant to the question to me.
Again, my question was: how can it be that a few rounds of the loop
produce a "no observations" error, but most don't, given that the data
are exactly the same in all rounds.
3. To those interested on what I am doing, I am using the ml set of
commands to minimize a quite convoluted loss function with both linear
and non linear restrictions. The linear constraints are handled with the
constraint define command. The non-linear constraints are handled by
making the evaluator add a punishment to the result for those vectors of
parameters outside the permitted space. I repeat the optimization many
times because (a) sometimes it does not converge (about 3-5% of cases),
and (b) because it is a convoluted loss function, and thus I want to get
some assurance that what I am finding is a global and not a local maximum.
4. It would be better that Stata had an optimization command not
specific for ml, that accepted constraints (ideally, both linear and non
linear), but that is not the case (nl does not accepts constraints of
any type; besides, it is not really powerful compared to ml -- with all
its excellent commands like ml check, ml search, the subroutines
available to write evaluators, the four optimizations algorithms it can
use, etc.). Using ml to run an optimization of something that is not a
maximum likelihood function has its challenges, but the pay-off is great
because the ml group of commands is so incredible good . . .
However, given that the code has been already developed, I really think
Stata should offer a general purpose optimization command, which takes
advantage of all the tools included in ml, but that does not require you
to "trick" it in order for it to work. Moreover, if you are not doing
maximum likelihood, the syntax of ml model is not really what you want –
you want to be able to simply indicate what the vector of parameters is,
avoiding all the stuff about equations (which, if you are not doing ml,
are completely artificial) that ml requires. And it should be possible
to include non-linear as well as linear constraints. Rival packages like
SPSS do have an easy-to-use general purpose command for optimization,
which admits both linear and non-linear constraints (but which I think
is not as powerful as ml). And, of course, Matlab also has such a thing ..
5. The full code of the do file is attached.
I hope I haven't violated the list's etiquette code by posting such a
long message, or in some other way -- but I am sure you will let me know
if I did.
Thank you again for your help.
Best,
Pablo
--
Pablo A. Mitnik
University of Wisconsin-Madison (http://www.wisc.edu/ )
Department of Sociology ( http://www.ssc.wisc.edu/soc/)
Center on Wisconsin Strategy (http://www.cows.org/ )
1180 Observatory Drive
Room 7114A
Madison, WI 53706
TEL (608) 2621839
E-mail: [email protected]
log using model6p_lf_rep_search.log, replace
set more off
clear
/*Define range of run*/
local list_of_ind = 3
local n = 20
foreach ind of numlist `list_of_ind' {
/*Data to be used in the optimization*/
/* for layouts of matrices see file layouts.xls */
use "/project/cows/career_ladders/Ladders6/competitors"
keep if ourind==`ind'
sort startyear below_or_outside edlevel
mkmat competitors, matrix(C)
/*mat list C*/
drop _all
use "/project/cows/career_ladders/Ladders6/openings"
keep if ourind==`ind'
sort startyear
mkmat openings, matrix(Y)
/*mat list Y*/
drop _all
use "/project/cows/career_ladders/Ladders6/mobilities"
keep if ourind==`ind'
sort startyear below_or_outside edlevel
mkmat mob_rate, matrix(Z)
foreach rep of numlist 1/`n'{
drop _all
svmat double Z
rename Z1 mob_rate
/*Defines the linear constraint; "[theta#]_cons" denotes the constant in equation # (see below)*/
constraint define 1 [theta1]_cons + [theta2]_cons + [theta3]_cons + [theta4]_cons = 1
/*Runs optimization; model6p_lf calculates loss-function*/
capture noisily {
ml model lf model6p_lf (alpha:mob_rate=) (beta:) (theta1:) (theta2:) (theta3:) (theta4:), technique (nr bfgs dfp bhhh) constraints(1)
if (`rep'==1) {
ml check
}
ml search alpha 0 1 beta 0 70 theta1 0 1 theta2 0 1 theta3 0 1 theta4 0 1, repeat(30)
ml maximize
/*transforms column vector with mobility rate hat for one rep into a variable, renames, saves it*/
svmat double MRH, names(mrh)
rename mrh1 mrh_ind`ind'_rep`rep'
drop mob_rate
save mrh_ind`ind'_rep`rep', replace
drop _all
/*defines matrix with loss function and parameters from ereturns*/
mat par_ind`ind'rep`rep'= e(ll), e(b), [alpha]_se[_cons], [beta]_se[_cons], [theta1]_se[_cons], [theta2]_se[_cons], [theta3]_se[_cons], [theta4]_se[_cons]
/*transforms above matrix into variables, adds variable for rep #, renames other vars, and saves*/
svmat double par_ind`ind'rep`rep'
rename par_ind`ind'rep`rep'1 loss_function
rename par_ind`ind'rep`rep'2 alpha
rename par_ind`ind'rep`rep'3 beta
rename par_ind`ind'rep`rep'4 theta1
rename par_ind`ind'rep`rep'5 theta2
rename par_ind`ind'rep`rep'6 theta3
rename par_ind`ind'rep`rep'7 theta4
rename par_ind`ind'rep`rep'8 alpha_se
rename par_ind`ind'rep`rep'9 beta_se
rename par_ind`ind'rep`rep'10 theta1_se
rename par_ind`ind'rep`rep'11 theta2_se
rename par_ind`ind'rep`rep'12 theta3_se
rename par_ind`ind'rep`rep'13 theta4_se
generate rep =`rep'
order rep
save par_ind`ind'_rep`rep', replace
} /*ends capture noisily*/
if _rc!=0 {
display ""
display ""
capture noisily display "ERROR #"
capture noisily display _rc
display ""
display" "
display" "
}
} /*ends repetition loop*/
/*appends all reps of parameteres for one industry, deletes working files*/
use par_ind`ind'_rep1.dta, clear
foreach rep of numlist 2/`n' {
capture noisily append using par_ind`ind'_rep`rep'.dta
}
save par_ind`ind'.dta, replace
foreach rep of numlist 1/`n' {
capture noisily erase par_ind`ind'_rep`rep'.dta
}
/*merges all reps of mob_rates_hat for one industry, deletes working files*/
use mrh_ind`ind'_rep1.dta, clear
foreach rep of numlist 2/`n' {
capture noisily merge using mrh_ind`ind'_rep`rep'.dta
capture noisily drop _merge
}
save mrh_ind`ind'.dta, replace
foreach rep of numlist 1/`n' {
capture noisily erase mrh_ind`ind'_rep`rep'.dta
}
} /*ends industry loop*/
log close
end