### 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.

### 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 error error handling errors Excel file i/o floating network GAUSS Engine graphics GUI hotkeys installation Java API license licensing linux loading data loops matrices matrix matrix manipulation Maxlik MaxLikMT Memory optimization Optmum output PQG graphics procs RAM random numbers string functions strings structures threading Time Series writing data### Recent Questions

### Features

### Time Series 2.0 MT

### Industry Solutions

### Find out more now

### Time Series MT 2.1

### Find out more now

# Resources

# A looping problem

Dear Gaussians,

I would like to estimate a 10-separate probit models within a simple loop. Staying with the same structure, I am getting results from Proc Iml in SAS, but Gauss terminates my program and goes off. Here is the programs in Gauss and SAS with proc Iml;

iter=10;

bstart= zeros(24,iter);

bprbt = zeros(24,iter);

for iter(1,iter,1);

bb = invpd(x'x)*x'y[.,iter];

bstart[.,iter]= bb;

{b0,fmean,grad,covp,retcode} = maxlik( "mydata",vars,&loglik,bstart[.,iter] );

call maxprt(b0,fmean,grad,covp,retcode);

proc loglik(betas,xx);

local y,x,xb,ll;

y = xx[.,iter];

x = xx[.,11:nyx];

xb= x*betas[1:nparam];

ll= y.*lncdfn(xb) + (1-y).*ln(1-cdfn(xb));

retp(ll);

endp;

bprbt[.,iter]=b0;

endfor;

SAS proc Iml;

prbest = j(kx,kd,**0**);

bstart = j(kx,kd,**0**);

* start doing loop;

do j = **1**to kd;

dd = d[,j];

bols = inv(x`*x)*x`*dd;

bstart[,j] = bols;

start LogLik(param) global (dd,x);

k = ncol(x);

n = nrow(x);

xb = x*param[**1**:k];

f = sum ( (dd=**1**) # log( probnorm(xb) ) + (dd=**0**) # log( probnorm(-xb) ) );

return ( f );

finish;

opt = {**3**,

**2**};

call nlpnra(rc, beta, "LogLik", bstart[,j], opt);

prbest[,j] = beta`;

call nlpfdd(f,g,h,"LogLik",beta);

end;

I would greatly appreciate for any valuable comments and suggestions.

with best wishes.

## 5 Answers

- What version of GAUSS are you running? (i.e. GAUSS 10 for Mac or GAUSS 14 for Windows 64-bit, etc)
- Where does
`y`come from in this line?bb = invpd(x’x)*x’y[.,iter];

- The variable
`nparam`is not defined in this code snippet. Where is it defined? - It will be very helpful, I think, to remove the definition of the
`proc``loglik`from inside the loop. This will make the code easier to understand.

Hi,

I am running Gauss 14 for windows 64-bit. "y" is a n by 10 matrix, 0/1(matrix), and "nparam" is the number of covariates defined previously. The program does not work when I remove the proc loglik outside the loop and if I remove "proc", then no need to compute the program with maxlik procedure, instead I should write a loop program with gradient and hessian, an usual newton algorithm, to compute my parameters. or as I did:

mm=seqa(1,1,10);

if mm == mm[1,.];

proc loglik1(betas,xx);

local y,x,b,xb;

y = xx[.,1];

x = xx[.,11:nyx];

b = betas[1:nparam];

xb= x*b;

retp( y.*lncdfn(xb)+(1-y).*ln(1-cdfn(xb)) );

endp;

elseif mm == mm[2,.];

proc loglik2(betas,xx);

local y,x,b,xb;

y = xx[.,2];

x = xx[.,11:nyx];

b = betas[1:nparam];

xb= x*b;

retp( y.*lncdfn(xb)+(1-y).*ln(1-cdfn(xb)) );

endp;

and so on.

I think your understanding of GAUSS proc's is not quite correct. In GAUSS, procedures are defined at compile time (i.e. when GAUSS turns the text in a file into a program that can be executed). The `if` and `for` statements, etc happen at run-time.

Since the proc's are already defined at compile time, whether they are inside of a particular `if` clause does not make any difference. For example:

x = 0; if x == 0; proc (1) = foo(x); retp(x.*x); endp; else; proc (1) = foo2(x); retp(2.*x); endp; endif; print foo(4); print foo2(4);

Will return exactly the same output as:

x = 0; if x == 0; else; endif; print foo(4); print foo2(4); proc (1) = foo2(x); retp(2.*x); endp; proc (1) = foo(x); retp(x.*x); endp;

This is because, even though GAUSS is flexible enough to let you define a procedure most anywhere, the procedure will not be changed during run-time.

I think you can accomplish what you are trying to do with an array of function pointers. I will have to look over the code more. I will post something later.

Here is an example using an array of function pointers so that you can pass in a different function for each iteration from the loop. In this toy example, the procedure `callProcs` would be like `maxlik` in your example.

//Create an array of function pointers fn_pointer = &square|&half|&cube; //Pass in a different element of the //function pointer array on each iteration for i(1, 3, 1); print callProcs(fn_pointer[i], 4); endfor; proc (0) = callProcs(local_fn_pointer, x); local local_fn_pointer:proc; print local_fn_pointer(x); endp; proc (1) = square(x); retp(x.*x); endp; proc (1) = cube(x); retp(x.*x.*x); endp; proc (1) = half(x); retp(x./2); endp;

Dear All,

I am very grateful to all of you for your kind help and guidance. I finally come up with results after your suggestions as follows:

bstart= zeros(kz,kd);

bprbt = zeros(kz,kd);

b0 = zeros(kz,kd);

covb = zeros(kz,kd*kz);

lls = zeros(1,kd);

i=0;

do while i<=10;

bb = invpd(z'z)*z'd[.,i]; /*starting values from linear prob*/

bstart[.,i]= bb;

i = i+1;

{b0,fmean,grad,covp,retcode} = maxlik( "mydata1",varbdz,&loglik,bstart[.,i] );

call maxprt(b0,fmean,grad,covp,retcode);

bprbt[.,i]=b0;

covb[.,(i-1).*kz+1:(i).*kz] = covp;

lls[1,i] = nd.*fmean;

endo;

proc loglik(betas,dz);

local d,z,zb,ll;

d = dz[.,i];

z = dz[.,11:kdz];

zb = z*betas[1:kz];

ll = d.*lncdfn(zb) + (1-d).*ln(1-cdfn(zb));

retp(ll);

endp;

It produces all results but with an error indication (**G0058** : Index out of range [*maxlik.src*]). I will also try other suggestions. Thanks.