SAS IML约束一个被调用的函数



我如何正确地约束这个最小化函数?Mincvf(cvf1)应该最小化cvf1相对于h,我想设置使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;

但是我得到了很多:

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

这是因为当h接近0时失去了精度,而"mhatx2"中的"le"接近0。当h=0.4时,le = ~0.08所以我只是人为地选择了这个下限,它仍然足够精确。

同样,"外部"子程序的输出ttt是h的向量,它适合于滚动窗口,仍然提供低于约束0.4的值。为什么?

我以前通过简单地对输入应用乘法变换来解决精度损失问题…将其乘以10,000或任何必要的值,然后在最后恢复转换。

不确定它是否适用于你的情况,但它可能值得一试。

这样做的话,必须将该选项和约束向量都放入输入括号中:

现在我没有得到除0的警告。先前未指定的由于精度损失的点现在根本不指定,其值被0.14替换,但错误可能不大。

开始mincvf (CVF1);

con={0.14…,…,…};

optn = {0};

init = 1;

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

返回(res);

完成;

最新更新