### Aptech Systems, Inc. Worldwide Headquarters

Address:

Aptech Systems, Inc.

2350 East Germann Road, Suite #21

Chandler, AZ 85286Phone: 360.886.7100

FAX: 360.886.8922**Ready to Get Started?**### Request Quote & Product Information

### Industry Solutions

### Products

### Resources

### Support

### Training & Events

Want more guidance while learning about the full functionality of GAUSS and its capabilities? Get in touch for in-person training or browse additional references below.

### Tutorials

Step-by-step, informative lessons for those who want to dive into GAUSS and achieve their goals, fast.

### What’s New in GAUSS 17

### Want to find out more?

### Have a Specific Question?

### Q&A: Register and Login

### Support Plans

Premier Support and Platinum Premier Support are annually renewable membership programs that provide you with important benefits including technical support, product maintenance, and substantial cost-saving features for your GAUSS System or the GAUSS Engine.

### User Forums

Join our community to see why our users are considered some of the most active and helpful in the industry!

### Where to Buy

Available across the globe, you can have access to GAUSS no matter where you are.

### Recent Tags

applications character vectors CML CMLMT Constrained Optimization datasets dates dlibrary dllcall econometrics Editor error error codes error handling errors Excel file i/o GAUSS Engine GAUSS Light graphics GUI hotkeys installation license licensing linux loading data loop loops matrix manipulation Maximum Likelihood Maxlik MaxLikMT Memory multidimensional array optimization Optmum output panel data PQG graphics procs random numbers strings structures threading### Recent Questions

### Features

### Time Series 2.0 MT

### Industry Solutions

### Find out more now

### Time Series MT 2.1

### Find out more now

### Find out more now

# Resources

# CMLMT error trap

Suppose I am running a Monte Carlo, and in each simulation, CMLMT is called. Most of the time, CMLML works find, but there are occasions when CML may fail, for example, due to "Not implemented for complex matrices." If this happens, I would like to discard this simulation and start a new simulation. Typically, I do not know what type of error CMLMT may create, as the simulated data may be ill-behaved in many different ways. Is there a decent way of "trapping" such that whenever CMLMT generates errors (and the simulation will stop if without "trapping"), the simulation loop will skip this round of simulation and continue to the next simulation?

## 3 Answers

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);

## Your Answer

## 3 Answers

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);