Hello... Here is my SAS program. I want to point out that it is not only the termination criteria values that differ from the SAS documentation, but their names also ( ABSXCONV versus ABSXTOL) prociml;   start ARL(hr); shift=5;    m=100;   t=2*m+1;   delta=2*hr[1]/t;    z=-hr[1]+.5*delta; do i=2 to t;   k=-hr[1]+(i-.5)*delta;   z=z//k;   end;         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<-(hr[2]*hr[3])then fun1=tau1-(1-hr[2])*hr[3];      if abs(tau1)<=(hr[2]*hr[3])then fun1=tau1/hr[2];     if tau1>(hr[2]*hr[3])then fun1=tau1+(1-hr[2])*hr[3];    if tau2<-(hr[2]*hr[3])then fun2=tau2-(1-hr[2])*hr[3];    if abs(tau2)<=(hr[2]*hr[3])then fun2=tau2/hr[2];    if tau2>(hr[2]*hr[3])then fun2=tau2+(1-hr[2])*hr[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,];   r=t(r);   avrl=inv((I(t)-r))*j(t,1,1);   bb=(v*avrl);   return(bb);   finish;    start ARL0(hr);   m1=100;   shift=0;   bb0=j(1,1,0.);   t1=2*m1+1;   delta1=2*hr[1]/t1;   z1=-hr[1]+.5*delta1;   do i1=2 to t1;   k1=-hr[1]+(i1-.5)*delta1;   z1=z1//k1;   end;   do i1=1 to t1;         zz1=0;                 do j1=1 to t1;                  tau11=z1[j1, ]+(delta1/2)-z1[i1,];                 tau21=z1[j1, ]-(delta1/2)-z1[i1,];   if tau11<-(hr[2]*hr[3])then fun11=tau11-(1-hr[2])*hr[3];   if abs(tau11)<=(hr[2]*hr[3])then fun11=tau11/hr[2];    if tau11>(hr[2]*hr[3])then fun11=tau11+(1-hr[2])*hr[3]; if tau21<-(hr[2]*hr[3])then fun21=tau21-(1-hr[2])*hr[3];    if abs(tau21)<=(hr[2]*hr[3])then fun21=tau21/hr[2];    if tau21>(hr[2]*hr[3])then fun21=tau21+(1-hr[2])*hr[3];            kone1=z1[i1,]+fun11;            two1=z1[i1,]+fun21;            pij1=probnorm(kone1-shift)-probnorm(ktwo1-shift);                 zz1=zz1//pij1;                 end;         if i1=1 then r1=zz1;          else          r1=r1||zz1;         end;   v1=j(t1,1,0);         v1[((t1+1)/2),]=1;        v1=t(v1);   r1=r1[2:t1+1,];   r1=t(r1);   avrl1=inv((I(t1)-r1))*j(t1,1,1);   bb0[1]=1000-(v1*avrl1);   return(bb0);   finish;   con = {   0    0    0,                     .    1    .    };   hr={2.015  0.354  1.989};    optn= j(1,11,.);optn[2]=5;optn[10]=1;optn[11]=0;   CALL NLPNMS(rc,hrres,"ARL",hr) blc=con opt=optn nlc="ARL0";   quit; ... View more

Thank you so much for your reply. I am using SAS 9.0. I want to check that the termination criteria defaults for the NLPNMS subroutine are the same as those of later SAS versions. In the SAS IML Non-linear Optimization document, I found that the default for the "absolute parameter convergence criterion(ABSXTOL)," is 1E-08 and However, when running my code I got this:                                                                                                                               Nelder-Mead Simplex Optimization                                                                                                                                                           COBYLA Algorithm by M.J.D. Powell (1992)                                                                                                                                                                                                                                                  Minimum Iterations                                     0                                         Maximum Iterations                                  1000                                         Maximum Function Calls                              3000                                         Iterations Reducing Constraint Violation               0                                         ABSFCONV Function Criterion                            0                                         FCONV Function Criterion                    2.220446E-16                                         FCONV2 Function Criterion                           1E-6                                         FSIZE Parameter                                        0                                         ABSXCONV Parameter Change Criterion               0.0001                                         XCONV Parameter Change Criterion                       0                                         XSIZE Parameter                                        0                                         ABSCONV Function Criterion                  -1.34078E154                                         Initial Simplex Size (INSTEP)                        0.5                                         Singularity Tolerance (SINGULAR)                    1E-8    Does this mean that the absolute parameter convergence criterion (which is now called  ABSXCONV) is 0.0001 not 1E-08?!. I didn't change the defaults at all, so why are the termination criteria different?? Is ABSXCONV the same as ABSXTOL?? I am really confused!! Aya ... View more

Thank you so much, but are these functions available in SAS 9.0?? A. Emara ... View more

Thank u so much , but actually i have a constraint which restricts my minimum value, I am not free to choose any value for my parameters so that the objective function is minimized. I need arlo>=200. The syntax gives answers that seem reasonable but they take very long time. Hutch suggested putting x in a global clause. Would this help and how could this be done? Could i fix the random generations through different optimization iterations so that at every guess of the solution , the generated data that is used is the same?? I hope I have clarified my point. ... View more

Yes, this is very close to the truth. I am trying to find the optimal parameters vector (hr) where the objective function represent the average of a certain function over 100,000 runs . This function counts the number of generated values until a certain metric "ti" exceeds the parameter hr[1]. I have one constraint which takes a form very close to the objective function form. I added the comments to the program Really appreciate your help. A. Emara prociml; start ARL global(hr);                     shift=0.5;     * p is number of variables;    p=2;    * nruns is the number of simulation runs;     nruns=100000;   arlv=0; * starting the simulation; do reps=1 to nruns;                     rl=0;   *initialy z is set as a px1 vector of zeros;    z=0//j(p-1,1,0);    *ti carries the sum of squares of z, i.e, each element of z is squared and then we take the sum;   ti=ssq(z);    do until (ti>hr[1]);    *generating p univariate standard normal variables;    x=normal(repeat(-1,p,1));    *adding the shift for the out of control case;   sh=shift//j(p-1,1,0);   x=x+sh;    *computing the norm of the error;   ei=sqrt(ssq(x-z));      if ei>hr[3] then w=(ei-(1-hr[2])*hr[3])/ei;      else w=hr[2];      z=(1-w)*z+w*x;    ti=ssq(z);      rl=rl+1;         end;    arlv=arlv+rl;   end;   arlr=arlv/nruns;   bb=arlr;   return(bb);   finish;    start ARL0(hry)global(hr);   hr=hry;   shift=0;   bb0=j(1,1,0.);   p=2;   * nruns is the number of simulation runs;    nruns=100000;   arlv=0;    * starting the simulation;    do reps=1 to nruns;   rl=0;   *initialy z is set as a px1 vector Of zeros;    z=0//j(p-1,1,0);    *ti carries the sum of squares of Z, i.e, each element of z is squared and then we take the sum;    ti=ssq(z);    do until (ti>hr[1]);    *generating p univariate standard normal variables;     x=normal(repeat(-1,p,1));    *adding the shift for the out of control case;    sh=shift//j(p-1,1,0);   x=x+sh;   *computing the norm of the error;      ei=sqrt(ssq(x-z));      if ei>hr[3] then w=(ei-(1-hr[2])*hr[3])/ei;      else w=hr[2];      z=(1-w)*z+w*x; ti=ssq(z);     rl=rl+1;   end;   arlv=arlv+rl;    end;   arlr=arlv/nruns;   bb0[1]=(arlr)-200; return(bb0);   finish;   con = {   0    0    0,                         .    1    . }; *initial vector; hr={0.8856  0.17  3.85};   optn= j(1,11,.);optn[2]=5;optn[10]=1;optn[11]=0;   CALL NLPNMS(rc,hrres,"ARL",hr) blc=con opt=optn nlc="ARL0";   quit; ... View more

Hello Rick,I tried the following bounds for the inner integral [-10,10] . I got errors concerning matrix inverse and convergence. The log file is attached ... View more

Sorry Rrick for troubling you, I really appreciate your continious efforts. I tried your first suggestion "scaling the integral". As I understood, I should substitute with ZED=0 and QQ=0 in the integrand. However, this substitution yields h(0)=0. Am I correct? Aya ... View more

Sorry for repeating my question, but I don't know what do you mean by scaling the integral? I tried to understand the procedure you explained but I couldn't get the point. Thanks, Aya Emara ... View more

No, I don't know the antiderivative information. What I know is that the v*avrl is the average run length of a certain control chart. I am computing the average run length using the markov chain technique. The martrix "r" is the transitional probability matrix and I have t  transiet states. Aya ... View more