Dear list,
I ask for some help. (Trust me, I'm not just outsourcing my homework,
I spent quite some hours bugfixing, only to partial avail.)
I try to code up an EM algorithm, with alternating steps updating the
observation-wise mixing probabilities and maximizing the likelihood
with the mixing given. I think the first step should be
straightforward, but it requires me to have some extra arguments that
-moptimize- loads into the evaluator and maybe that is the cause of my
problems (an extra -st_view()-, with variable names kept track of
etc.). Otherwise my model stays really close to the manual's example
of a classical normal linear model with -moptimize-, so I'm puzzled
which step is invalid in Mata.
Also, if you could help me with any bugfixing options in Mata that
helps me check errors and trace the execution even if other codes call
the subroutines, I'd appreciate it. I have -set matadebug on-, -set
trace on- and -mata: mata set matalnum on-, yet as my code is called
by other code (namely -moptimize-), I barely get any error messages or
trace. And it is a waste of my time to copy parts of my code to
partial do-files just to track its run better.
Thank you for any help.
Laszlo
Here is my evaluator:
function mstep(transmorphic M,
real rowvector b,
real colvector lnf)
{
real colvector p1, p2, s1, s2
real colvector y
p1 = moptimize_util_xb(M, b, 1)
p2 = moptimize_util_xb(M, b, 2)
s1 = moptimize_util_xb(M, b, 3)
s2 = moptimize_util_xb(M, b, 4)
y = moptimize_util_depvar(M, 1)
st_view(mix,.,moptimize_util_userinfo(M,1),moptimize_util_userinfo(M,2))
lnf = ln(mix:*normalden(y:-p1, 0,
sqrt(s1))+(J(rows(mix),1,1):-mix):*normalden(y:-p2, 0, sqrt(s2)))
}
mata mosave mstep(), replace dir(PERSONAL)
Here is the call from the program (lhs, rhs, touse and mix are string
scalars from the call):
// M-step
M = moptimize_init()
moptimize_init_evaluator(M, &mstep())
moptimize_init_evaluatortype(M,"lf")
moptimize_init_touse(M, touse)
moptimize_init_ndepvars(M, 1)
moptimize_init_depvar(M, 1, lhs)
moptimize_init_eq_n(M, 4)
moptimize_init_eq_indepvars(M, 1, rhs)
moptimize_init_eq_cons(M, 1, "on")
moptimize_init_eq_indepvars(M, 2, rhs)
moptimize_init_eq_cons(M, 2, "on")
moptimize_init_eq_indepvars(M, 3, "")
moptimize_init_eq_cons(M, 1, "on")
moptimize_init_eq_indepvars(M, 4, "")
moptimize_init_eq_cons(M, 2, "on")
moptimize_init_eq_name(M, 1, "class1")
moptimize_init_eq_name(M, 2, "class2")
moptimize_init_eq_name(M, 3, "se1")
moptimize_init_eq_name(M, 4, "se2")
// moptimize_init_by(M, by)
moptimize_init_eq_coefs(M, 1, b1)
moptimize_init_eq_coefs(M, 2, b2)
moptimize_init_eq_coefs(M, 3, s1)
moptimize_init_eq_coefs(M, 4, s2)
moptimize_init_search(M, "on")
moptimize_init_search_random(M, "on")
moptimize_init_search_bounds(M, 3, (0,.))
moptimize_init_search_bounds(M, 4, (0,.))
moptimize_init_valueid(M, "LogL with given mixing")
moptimize_init_nuserinfo(M, 2)
moptimize_init_userinfo(M, 1, mix)
moptimize_init_userinfo(M, 2, touse)
moptimize_init_tracelevel(M,"hessian")
_moptimize(M)
I get the error message:
: _moptimize(M)
execution begins..
st_view(): 3204 matrix found where scalar required
mstep(): - function returned error [13] {103}
mopt__calluser_lf(): - function returned error {79}
opt__lf_calluser(): - function returned error {66}
deriv__call1user_v(): - function returned error {91}
_deriv__compute_value(): - function returned error {165}
_deriv(): - function returned error {187}
opt__eval_nr_lf(): - function returned error {130}
opt__eval(): - function returned error {98}
_optimize_evaluate(): - function returned error {126}
_mopt__evaluate(): - function returned error {76}
_moptimize_evaluate(): - function returned error {87}
_moptimize_search(): - function returned error {77}
_moptimize(): - function returned error {77}
<istmt>: - function returned error [1] {11}
(11 lines skipped)
--------------------------------------------------------------------------------
r(3204);
*
* 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/