SAS IML constraining a called function

Question

How do I properly constrain this minimizing function? Mincvf(cvf1) should minimize cvf1 with respect to h and I want to set so that h>=0.4

proc iml;

EDIT kirjasto.basfraaka var "open";

read all var "open" into cp;
p=cp[1:150];

conh={0.4 . .,. . .,. . .};

m=nrow(p);
m2=38;
pi=constant("pi");
e=constant("e");


start Kmod(x,h,pi,e);
k=1/(h#(2#pi)##(1/2))#e##(-x##2/(2#h##2));
return (k);
finish;


start mhatx2 (m2,newp,h,pi,e);
t5=j(m2,1);                   /*mhatx omit x=t*/
do x=1 to m2;
i=T(1:m2);
temp1=x-i;
ue=Kmod(temp1,h,pi,e)#newp[i];
le=Kmod(temp1,h,pi,e);
t5[x]=(sum(ue)-ue[x])/(sum(le)-le[x]);
end;
return (t5);
finish;

Start CVF1(h) global (newp,pi,e,m2);
cv3=j(m2,1);                                            
cv3=1/m2#sum((newp-mhatx2(m2,newp,h,pi,e))##2);
return(cv3);
finish;


start mincvf(CVF1);
optn={0,0};
init=1;
call nlpqn(rc, res,"CVF1",init) blc="conh";
return (res);
finish;


start outer(p,m) global(newp);
wl=38;    /*window length*/
m1=m-wl;    /*last window begins at m-wl*/
newp=j(wl,1);
hyi=j(m1,1);
do x=1 to m1;
we=x+wl-1;             /*window end*/
w=T(x:we);            /*window*/
newp=p[w];
hyi[x]=mincvf(CVF1);
end;
return (hyi);
finish;


wl=38;    /*window length*/
m1=m-wl;    /*last window begins at m-wl*/

time=T(1:m1);

ttt=j(m1,1);
ttt=outer(p,m);
print time ttt p;

However I get lots of:

WARNING: Division by zero, result set to missing value.

 count     : number of occurrences is 2
 operation : / at line 1622 column 22
 operands  : _TEM1003, _TEM1006

_TEM1003      1 row       1 col     (numeric)

         .

_TEM1006      1 row       1 col     (numeric)

         0

 statement : ASSIGN at line 1622 column 1
 traceback : module MHATX2 at line 1622 column 1
             module CVF1 at line 1629 column 1
             module MINCVF at line 1634 column 1
             module OUTER at line 1651 column 1

Which happens because losing of precision when h approaches 0 and "le" in "mhatx2" approaches 0. At value h=0.4, le is ~0.08 so I just artificially picked that as a lower bound which is still precise enough.

Also this output of "outer" subroutine, ttt which is vector of h fitted for the rolling windows, still provides values below the constraint 0.4. Why?

Answer 1

I have solved loss of precision issues previously by simply applying a multiplication transformation to the input... Multiply it by 10,000 or whatever is necessary, and then revert the transformation at the end.

Not sure if it will work in your situation, but it may be worth a shot.

Answer 2

This way it works, had to put that option and constrain vector both into the input parentheses:

Now I get no division by 0 warning. The previously miss-specified-due-loss-of-precision point's are now not specified at all and the value is substituted by 0.14 but the error isn't likely big.

start mincvf(CVF1);

con={0.14 . .,. . .,. . .};

optn={0,0};

init=1;

call nlpqn(rc, res,"CVF1",init,optn,con);

return (res);

finish;