The obstacle is here: rho_var = partialderiv*Sigma*t(partialderiv);
It is quadratic form .
data partialderiv ;
set Covparms end=last;
retain a b d g h i;
if _n_=1 then a=Estimate;
if _n_=2 then b=Estimate;
if _n_=3 then c=Estimate;
if _n_=4 then g=Estimate;
if _n_=5 then h=Estimate;
if _n_=6 then i=Estimate;
if last then do;
rhohat = (b+h)/sqrt((a+g)*(c+i));
df_da = -0.5*((b+h)*(c+i))/sqrt(((a+g)*(c+i))**3);
df_db = 1/sqrt((a+g)*(c+i));
df_dc = -0.5*(b+h)*(a+g)/sqrt(((a+g)*(c+i))**3);
df_dg = -0.5*(b+h)*(c+i)/sqrt(((a+g)*(c+i))**3);
df_dh = 1/sqrt((a+g)*(c+i));
df_di = -0.5*(b+h)*(a+g)/sqrt(((a+g)*(c+i))**3);
output;
end;
keep rhohat df_da df_db df_dc df_dg df_dh df_di; /*<--*/
run;
/***********************/
/***********************/
data _null_;
set partialderiv;
temp=catx(' ',df_da, df_db, df_dc, df_dg, df_dh, df_di);
call symputx('value',temp);
stop;
run;
data cov;
set Asycov;
array x{6} _temporary_ (&value);
array c{*} CovP1-CovP6;
i+1;
do j=1 to dim(c);
c{j}=c{j}*x{i}*x{j};
end;
drop i j;
run;
proc sql;
select sum(sum(CovP1, CovP2, CovP3 ,CovP4 ,CovP5, CovP6)) into : rho_var
from cov;
quit;
data want;
set partialderiv;
rho_SE = sqrt(&rho_var);
l95=rhohat-1.96*rho_SE;
u95=rhohat+1.96*rho_SE;
p = (1-cdf("normal",abs(rhohat),0,rho_SE))*2;
label
rhohat='Estimate';
rho_SE='Std Error';
l95='95% CI Lower Bound';
u95='95% CI UpperBound';
p='Pvalue' ;
run;
... View more