CMLMT error trap

If these errors are happening inside of your likelihood function, you could check these inputs to see if they are complex or Nans. For example, if we change the likelihood procedure in cmlmt1.e to artificially create a complex number before it is handed off to cdfn which does not accept complex numbers like this:

proc lpr(struct PV p, struct DS d, ind);
   local s2,b0,b,y,x,yh,u,res,g1,g2;

   struct modelResults mm;

   b0 = pvUnpack(p,"b0");
   b = pvUnpack(p,"b");
   s2 = pvUnpack(p,"variance");

   y = d[1].dataMatrix;
   x = d[2].dataMatrix;

   yh = b0 + x * b;
   res = y - yh;
   u = y[.,1] ./= 0;

   if ind[1];
       mm.function = u.*lnpdfmvn(res,s2) + (1-u).*(ln(cdfnc(yh/sqrt(s2))));
   endif;

   if ind[2];
      //make negative so 'yh' will be complex
      s2 = -s2;
      yh = yh/sqrt(s2);
      //stop program if 'yh' is complex
      if iscplx(yh);
         stop;
      endif;
      g1 = ((res~x.*res)/s2)~((res.*res/s2)-1)/(2*s2);
      g2 = ( -( ones(rows(x),1)~x )/sqrt(s2) )~(yh/(2*s2));
      g2 = (pdfn(yh)./cdfnc(yh)).*g2;
      mm.gradient = u.*g1 + (1-u).*g2;
   endif;

   retp(mm);

endp;

Now with this new likelihood procedure, the example file cmlmt1.e will catch the problem before it causes an error and the program will stop. Once we have done this, there is still one problem. You don't want the entire program to stop.

You can accomplish this by running this file with the run command using the -r flag. For example, if you made a main file which we will call cmlmtdriver.gss with these contents:

for i(1, 3, 1);
   print "starting iteration " i;
   run -r cmlmt1.e;
endfor;
print "main file still running";

If you have any questions about how to adapt these techniques to your simulation program, please ask.

Other useful functions
Besides just checking whether or not a matrix is complex with iscplx, you can check whether or not a matrix has a complex portion that is larger than a certain tolerance with the hasimag function. This will allow you to retain only the real portion of the matrix if the complex portion is 0 to machine precision.

The following GAUSS functions may also be helpful in this context: ismiss, isinfnanmiss, scalmiss.

There is also the trap function which can be used with chol, invpd, or solpd to trap in the case that the input matrix is not positive definite. trap can also be used with the slash operator and inv to trap cases in which the input matrix is singular, or with the eigenvalue functions to trap the case in which the eigenvalues cannot be calculated.

Thanks!

My code indicates that "cholup" and "choldn" sometimes give "error G0043 : Not implemented for complex matrices." The code from CMLMT.src has the following:
oldt = trapchk(1);
trap 1,1;
...
code involving cholup/choldn
...
trap oldt,1;

It appears that error trapping does not here for "cholup" and "choldn."

You are correct, cholup and choldn do not support the trap command.

One thing to know is that if a complex matrix is sent to cholup or choldn, they do perform a check on the complex portion of the input to see if it is significant like the hasimag function.

It looks like all the calls to cholup and choldn in CMLMT take the hessian as an input. I think the first step I would take is to add an iscplx check to any outputs from your likelihood function and add a stop if the output was complex.

If the programs take a while to run, I might also add an if iscplx() check to some of the other variables being passed into choldn or cholup, particularly g0 and v2. Or if you just want to be done with it, you could take a line like this:

h1 = cholup(h1,g1/v3-g0/v3+ak*dx/v3);

and change it to something like this:

local cholup_input;
cholup_input = g1/v3-g0/v3+ak*dx/v3;
if iscplx(h1) or iscplx(cholup_input);
   errorlogat "exiting CMLMT due to unexpected complex matrix";
   stop;
endif;
h1 = cholup(h1,cholup_input);