topic Re: Number of iterations in NLPNM suproutine in SAS/IML Software and Matrix Computations

Number of iterations in NLPNM suproutine

RaniaMamdouh — Sat, 30 Nov 2013 20:45:25 GMT

Hello all;

How can i increase the number of iterations in proc IML, NLPNMS subroutine?

thanks

Re: Number of iterations in NLPNM suproutine

Rick_SAS — Mon, 02 Dec 2013 13:38:45 GMT

The doc for the termination criteria is here: SAS/IML(R) 12.3 User's Guide

tc[1] specifies the maximum number of iterations in the optimization process. To get, say, 5000 iterations, use the following:

tc = j(10,1,.); /* set all to default (missing) */

tc[1] = 5000; /* override default maximum iterations */

Re: Number of iterations in NLPNM suproutine

RaniaMamdouh — Mon, 02 Dec 2013 20:36:59 GMT

Thanks for your reply

but when i apply your solution also i find the same following error

ERROR: NMSIMP Optimization cannot be completed.

ERROR: NMSIMP needs more than 1000 iterations or 3000 function calls.

what can i do in this proplem?

Re: Number of iterations in NLPNM suproutine

Hutch_sas — Mon, 02 Dec 2013 20:52:30 GMT

Looks like you need to increase the second control parameter, the maximum number of function calls, i.e.

tc[2] = 10000;

Re: Number of iterations in NLPNM suproutine

Rick_SAS — Mon, 02 Dec 2013 20:54:58 GMT

That is a difficult question to answer, as many books and papers have been written on the topic. A few things to try:

1) Look at the return code to see why it is failing.

2) Look at the iteration history to see whether the optimization is jumping around, diverging to infinity, stuck in a valley/ridge, or some other situation.

3) Try a better initial guess.

4) Relax some of the convergence criteria. If the objective function is flat (see #1 and #2), you might not be able to converge by using the default values.

5) Try a different algorithm.

Re: Number of iterations in NLPNM suproutine

RaniaMamdouh — Sun, 08 Dec 2013 21:00:24 GMT

I tried these things, but i think that the problem can be in writing the routine. I will briefly describe my problem, i minimize the "cost" and i have 4 linear constrains on "x's" and three non-linear constrains on (arl), and this is the last part of my program which i think that something wrong in it:

tc = j(10,1,.); /* set all to default (missing) */

tc[1] = 5000; /* override default maximum iterations */

con = { 0. 0. 0. 1.,

1. . . 20.};

x={0.1 0.7 3.5 2};

optn= j(1,11,.);optn[1]=0;optn[2]= 5;optn[10]=3;optn[11]=0;

CALL nlpnms(rc,xres,"Cost",x,optn,con,tc) nlc="arl";

quit;

please tell me if there is wrong.

thanks

Re: Number of iterations in NLPNM suproutine

Rick_SAS — Mon, 09 Dec 2013 16:21:39 GMT

I don't see any linear constraints for x. All I see are ranges. You need four additional rows for the CON matrix. See SAS/IML(R) 12.3 User's Guide

Re: Number of iterations in NLPNM suproutine

RaniaMamdouh — Mon, 09 Dec 2013 17:23:05 GMT

yes, i mean that ranges for x's not linear constrains. Is there any other wrong?

Re: Number of iterations in NLPNM suproutine

Hutch_sas — Mon, 09 Dec 2013 17:44:27 GMT

You might try overriding the second termination criteria, set tc[2] = 10000 and see if that affects your results. it appears that you are running up against a limit on objective function calls.

Re: Number of iterations in NLPNM suproutine

Rick_SAS — Mon, 09 Dec 2013 17:50:25 GMT

Can you show the ARL module for nonlinear contstraints?

Re: Number of iterations in NLPNM suproutine

RaniaMamdouh — Mon, 09 Dec 2013 19:48:50 GMT

i'm already tried to set tc[2]=10000 but nothing improved. Is the place of "tc" in the CALL line right or not?

Re: Number of iterations in NLPNM suproutine

RaniaMamdouh — Mon, 09 Dec 2013 19:52:31 GMT

yes sure; it's as follows:

start arl(x);

bb0=j(3,1,0.);

m=100;

shift=0;

* L is ucl;

*x[1] is lamda, x[2] is L, x[3] is k;

* t is number of regions;

t=2*m+1;

* delta is width of each region;

delta=2*x[2]/t;

*now we are constructing the vector z;

*vector z represents different locations for the chart statistics;

*these are the midpoints of the t intervals;

z=-x[2]+.5*delta;

do i=2 to t;

k=-x[2]+(i-.5)*delta;

z=z//k;

end;

* now we are calculating the incontrol transitional probability matrix r;

do i=1 to t;

zz=0;

do j=1 to t;

tau1=z[j, ]+(delta/2)-z[i, ];

tau2=z[j, ]-(delta/2)-z[i, ];

if tau1<-(x[1]*x[3]) then fun1=tau1-(1-x[1])*x[3];

if abs(tau1)<=(x[1]*x[3])then fun1=tau1/x[1];

if tau1>(x[1]*x[3]) then fun1=tau1+(1-x[1])*x[3];

if tau2<-(x[1]*x[3]) then fun2=tau2-(1-x[1])*x[3];

if abs(tau2)<=(x[1]*x[3])then fun2=tau2/x[1];

if tau2>(x[1]*x[3]) then fun2=tau2+(1-x[1])*x[3];

kone=z[i,]+fun1;

ktwo=z[i,]+fun2;

pij=probnorm(kone-shift)-probnorm(ktwo-shift);

zz=zz//pij;

end;

if i=1 then r=zz;

else

r=r||zz;

end;

v=j(t,1,0);

v[((t+1)/2), ]=1;

v=t(v);

r=r[2:t+1,];

* now matrix r is t*t with the row of all zeros removed;

r=t(r);

*j(t,1,1) is tx1 vector of the value 1;

arl=inv((I(t)-r))*j(t,1,1);

arl0=v*arl;

bb0[1]=arl0-500;

m=100;

shift=1;

* L is ucl;

*x[1] is lamda, x[2] is L, x[3] is k;

* t is number of regions;

t=2*m+1;

* delta is width of each region;

delta=2*x[2]/t;

*now we are constructing the vector z;

*vector z represents different locations for the chart statistics;

*these are the midpoints of the t intervals;

z=-x[2]+.5*delta;

do i=2 to t;

k=-x[2]+(i-.5)*delta;

z=z//k;

end;

* now we are calculating the incontrol transitional probability matrix r;

do i=1 to t;

zz=0;

do j=1 to t;

tau1=z[j, ]+(delta/2)-z[i, ];

tau2=z[j, ]-(delta/2)-z[i, ];

if tau1<-(x[1]*x[3]) then fun1=tau1-(1-x[1])*x[3];

if abs(tau1)<=(x[1]*x[3])then fun1=tau1/x[1];

if tau1>(x[1]*x[3]) then fun1=tau1+(1-x[1])*x[3];

if tau2<-(x[1]*x[3]) then fun2=tau2-(1-x[1])*x[3];

if abs(tau2)<=(x[1]*x[3])then fun2=tau2/x[1];

if tau2>(x[1]*x[3]) then fun2=tau2+(1-x[1])*x[3];

kone=z[i,]+fun1;

ktwo=z[i,]+fun2;

pij=probnorm(kone-shift)-probnorm(ktwo-shift);

zz=zz//pij;

end;

if i=1 then r=zz;

else

r=r||zz;

end;

v=j(t,1,0);

v[((t+1)/2), ]=1;

v=t(v);

r=r[2:t+1,];

* now matrix r is t*t with the row of all zeros removed;

r=t(r);

*j(t,1,1) is tx1 vector of the value 1;

arl=inv((I(t)-r))*j(t,1,1);

arl1=v*arl;

bb0[2]=12-arl1;

m=100;

shift=0;

* L is ucl;

*x[1] is lamda, x[2] is L, x[3] is k, x[4] is n;

* t is number of regions;

t=2*m+1;

* delta is width of each region;

delta=2*x[2]/t;

*now we are constructing the vector z;

*vector z represents different locations for the chart statistics;

*these are the midpoints of the t intervals;

z=-x[2]+.5*delta;

do i=2 to t;

k=-x[2]+(i-.5)*delta;

z=z//k;

end;

* now we are calculating the incontrol transitional probability matrix r;

do i=1 to t;

zz=0;

do j=1 to t;

tau1=z[j, ]+(delta/2)-z[i, ];

tau2=z[j, ]-(delta/2)-z[i, ];

if tau1<-(x[1]*x[3]) then fun1=tau1-(1-x[1])*x[3];

if abs(tau1)<=(x[1]*x[3])then fun1=tau1/x[1];

if tau1>(x[1]*x[3]) then fun1=tau1+(1-x[1])*x[3];

if tau2<-(x[1]*x[3]) then fun2=tau2-(1-x[1])*x[3];

if abs(tau2)<=(x[1]*x[3])then fun2=tau2/x[1];

if tau2>(x[1]*x[3]) then fun2=tau2+(1-x[1])*x[3];

kone=z[i,]+fun1;

ktwo=z[i,]+fun2;

pij=probnorm(kone-shift)-probnorm(ktwo-shift);

zz=zz//pij;

end;

if i=1 then r=zz;

else

r=r||zz;

end;

v=j(t,1,0);

v[((t+1)/2), ]=1;

v=t(v);

r=r[2:t+1,];

* now matrix r is t*t with the row of all zeros removed;

r=t(r);

*j(t,1,1) is tx1 vector of the value 1;

arl=inv((I(t)-r))*j(t,1,1);

arl0=v*arl;

m=100;

shift=4;

* L is ucl;

*x[1] is lamda, x[2] is L, x[3] is k;

* t is number of regions;

t=2*m+1;

* delta is width of each region;

delta=2*x[2]/t;

*now we are constructing the vector z;

*vector z represents different locations for the chart statistics;

*these are the midpoints of the t intervals;

z=-x[2]+.5*delta;

do i=2 to t;

k=-x[2]+(i-.5)*delta;

z=z//k;

end;

* now we are calculating the incontrol transitional probability matrix r;

do i=1 to t;

zz=0;

do j=1 to t;

tau1=z[j, ]+(delta/2)-z[i, ];

tau2=z[j, ]-(delta/2)-z[i, ];

if tau1<-(x[1]*x[3]) then fun1=tau1-(1-x[1])*x[3];

if abs(tau1)<=(x[1]*x[3])then fun1=tau1/x[1];

if tau1>(x[1]*x[3]) then fun1=tau1+(1-x[1])*x[3];

if tau2<-(x[1]*x[3]) then fun2=tau2-(1-x[1])*x[3];

if abs(tau2)<=(x[1]*x[3])then fun2=tau2/x[1];

if tau2>(x[1]*x[3]) then fun2=tau2+(1-x[1])*x[3];

kone=z[i,]+fun1;

ktwo=z[i,]+fun2;

pij=probnorm(kone-shift)-probnorm(ktwo-shift);

zz=zz//pij;

end;

if i=1 then r=zz;

else

r=r||zz;

end;

v=j(t,1,0);

v[((t+1)/2), ]=1;

v=t(v);

r=r[2:t+1,];

* now matrix r is t*t with the row of all zeros removed;

r=t(r);

*j(t,1,1) is tx1 vector of the value 1;

arl=inv((I(t)-r))*j(t,1,1);

arl1=v*arl;

c0=10;

c1=100;

gam1=1;

gam2=1;

O=0.01;

W=25;

a=0.5;

b=0.1;

F=50;

E=0.05;

T0=0;

T1=2;

T2=0;

h4=(x[4]*E)+(gam1*T1)+(gam2*(T2**2));

h1=((arl1-0.5)*(O*(F+(c1*T0*(gam1-1)))-(2*arl0*(c0+O*(((arl1-0.5)*a)+((arl1-0.5)*b*x[4])+W)+(c1*(-1+h4*O-O*((x[4]*E)+T1+T2)))))))/(2*O*arl0);

h2=((-2*(arl1-0.5)*(F+(c1*T0*(gam1-1))+(arl0*(a+(b*x[4]))*(1+(h4*O))))))/(O*arl0);

h3=(-1/(2*(O**2)*arl0))*((2*F)+(2*c0*T0*(gam1-1))-(a*T0*O)-(2*(arl1-0.5)*a*T0*O)-(2*c1*h4*T0*O)-(b*x[4]*T0*O)-(2*(arl1-0.5)*b*x[4]*T0*O)-(2*T0*W*O)+(2*F*O*((x[4]*E)+T1+T2))+(a*T0*gam1*O)+(2*(arl1-0.5)*a*T0*gam1*O)+(2*c1*h4*T0*gam1*O)+(b*x[4]*T0*gam1*O)+(2*(arl1-0.5)*b*x[4]*T0*gam1*O)+(2*T0*W*gam1*O)-(a*h4*T0*(O**2))-(b*x[4]*h4*T0*(O**2))+(a*h4*T0*gam1*(O**2))+(a*x[4]*h4*T0*gam1*(O**2))+(2*arl0*(a+b*x[4])*(1+h4*O)*(1+O*((x[4]*E)+T1+T2))));

h=(-h2+SQRT((h2**2)-(4*h1*h3)))/(2*h1);

s=(1/(O*h))-0.5;

r1=(c0/O)+(c1*((x[4]*E)+(h*(arl1-0.5))+(gam1*T1)+(gam2*T2))+((s*F)/arl0)+W);

r2=((a+(b*x[4]))/h)*((1/O)+(x[4]*E)+(h*(arl1-0.5))+(gam1*T1)+(gam2*T2));

r3=(1/O)+((1-gam1)*s*T0)/arl0+(x[4]*E)+(h*(arl1-0.5))+T1+T2;

C2=(r1+r2)/r3;

print C2;

bb0[3]=(1.005*13.355)-C2;

return (bb0);

finish arl;

Re: Number of iterations in NLPNM suproutine

Rick_SAS — Mon, 09 Dec 2013 20:46:07 GMT

Ugh. Too complicated. I don't see how to verify that the ARL module is correct.

Can you solve the unconstrained problem? What about with only the first or only the first 2 constraints?

Re: Number of iterations in NLPNM suproutine

RaniaMamdouh — Tue, 10 Dec 2013 18:01:39 GMT

if i put only the third constrain i can face the same problem, i.e. the problem isn't in the number of constrains?

Re: Number of iterations in NLPNM suproutine

Rick_SAS — Tue, 10 Dec 2013 18:19:32 GMT

Is it possible that the problem is in the (very complicated) nonlinear constraints? Can you show that there is a feasible solution that satisfies the constraints? If not, you will encounter this kind of non-convergence.

To give a simple example, the problem min F(x,y) s.t. y> -x^2 and y>x^2 +1 has no solution because no (x,y) pair is feasible.

Re: Number of iterations in NLPNM suproutine

RaniaMamdouh — Thu, 26 Dec 2013 17:38:54 GMT

I'm sorry for the delay in reply because I have seen your comment immediately.

yes, there is a feasible solution satisfies the non-linear constrains. Is the place of "tc" in the CALL line right or not?

thanks

Re: Number of iterations in NLPNM suproutine

Rick_SAS — Sat, 28 Dec 2013 13:13:44 GMT

1) Give us a feasible solution

2) Yes, the tc arg is in the right place, but you can always use the alternative keyword-pair" syntax when in doubt:

CALL nlpnms(...) tc=tc nlc="arl";

Re: Number of iterations in NLPNM suproutine

RaniaMamdouh — Sat, 28 Dec 2013 17:03:49 GMT

A feasible solution is : Cost= 17.058 at n=1, h=0.28238, Lammda =0.13199, , L=0.82698, K=2.858