Hi all,
I am stuggling to debug the code below. Specifically, the code runs without generating errors in the SAS log, but the outer-most do loop, i.e. do cutoff=-10 to 10, does not run correctly. In the calculation of the "theory_stay" matrix, the values generated in each iteration of the loop are the same as the those in the "theory_stay" matrix from the first iteration, i.e. cutoff=-10. If I should attach some data or give more background into what the code is doing, just let me know. I appreciate all feedback!
Here is the code:
proc iml;
use work.output1;
read all var {jun04 jul04 aug04 sep04 oct04 nov04 dec04 jan05 feb05 mar05 apr05 may05 jun05 jul05 aug05 sep05 oct05 nov05 dec05 jan06 feb06 mar06 apr06 may06 jun06 jul06 aug06 sep06 oct06 nov06 dec06 jan07 feb07 mar07 apr07 may07 jun07 jul07 aug07 sep07 oct07 nov07 dec07 jan08 feb08 mar08 apr08 may08 jun08 jul08 aug08 sep08 oct08 nov08 dec08 jan09 feb09 mar09 apr09 may09 jun09 jul09 aug09 sep09 oct09 nov09 dec09 jan10 feb10 mar10 apr10 may10 jun10 jul10 aug10 sep10 oct10 nov10 dec10 jan11 feb11 mar11 apr11 may11 jun11 jul11 aug11 sep11 oct11 nov11 dec11 jan12 feb12 mar12 apr12 may12 jun12 jul12 aug12 sep12 oct12 nov12 dec12 jan13 feb13 mar13 apr13 may13 jun13 jul13 aug13 sep13 oct13 nov13 dec13 jan14 feb14 mar14 apr14 may14 jun14 jul14 aug14 sep14 oct14 nov14 dec14 jan15 feb15 mar15 apr15 may15} into month_lb;
n_months=ncol(month_lb);
n_farms=nrow(month_lb);
use work.avc1;
read all var {jun04 jul04 aug04 sep04 oct04 nov04 dec04 jan05 feb05 mar05 apr05 may05 jun05 jul05 aug05 sep05 oct05 nov05 dec05 jan06 feb06 mar06 apr06 may06 jun06 jul06 aug06 sep06 oct06 nov06 dec06 jan07 feb07 mar07 apr07 may07 jun07 jul07 aug07 sep07 oct07 nov07 dec07 jan08 feb08 mar08 apr08 may08 jun08 jul08 aug08 sep08 oct08 nov08 dec08 jan09 feb09 mar09 apr09 may09 jun09 jul09 aug09 sep09 oct09 nov09 dec09 jan10 feb10 mar10 apr10 may10 jun10 jul10 aug10 sep10 oct10 nov10 dec10 jan11 feb11 mar11 apr11 may11 jun11 jul11 aug11 sep11 oct11 nov11 dec11 jan12 feb12 mar12 apr12 may12 jun12 jul12 aug12 sep12 oct12 nov12 dec12 jan13 feb13 mar13 apr13 may13 jun13 jul13 aug13 sep13 oct13 nov13 dec13 jan14 feb14 mar14 apr14 may14 jun14 jul14 aug14 sep14 oct14 nov14 dec14 jan15 feb15 mar15 apr15 may15} into avc;
n_months=ncol(avc);
n_farms=nrow(avc);
use work.prices1;
read all var {jun04 jul04 aug04 sep04 oct04 nov04 dec04 jan05 feb05 mar05 apr05 may05 jun05 jul05 aug05 sep05 oct05 nov05 dec05 jan06 feb06 mar06 apr06 may06 jun06 jul06 aug06 sep06 oct06 nov06 dec06 jan07 feb07 mar07 apr07 may07 jun07 jul07 aug07 sep07 oct07 nov07 dec07 jan08 feb08 mar08 apr08 may08 jun08 jul08 aug08 sep08 oct08 nov08 dec08 jan09 feb09 mar09 apr09 may09 jun09 jul09 aug09 sep09 oct09 nov09 dec09 jan10 feb10 mar10 apr10 may10 jun10 jul10 aug10 sep10 oct10 nov10 dec10 jan11 feb11 mar11 apr11 may11 jun11 jul11 aug11 sep11 oct11 nov11 dec11 jan12 feb12 mar12 apr12 may12 jun12 jul12 aug12 sep12 oct12 nov12 dec12 jan13 feb13 mar13 apr13 may13 jun13 jul13 aug13 sep13 oct13 nov13 dec13 jan14 feb14 mar14 apr14 may14 jun14 jul14 aug14 sep14 oct14 nov14 dec14 jan15 feb15 mar15 apr15 may15} into price;
n_months=ncol(price);
n_farms=nrow(price);
call randseed(datetime());
u=j(n_farms,n_months,.);
mu=0;
sigma= 4.0233929;
call randgen(u,"normal",mu,sigma);
stoch_avc=j(n_farms,n_months,0);
month_cwt=j(n_farms,n_months,0);
month_prof=j(n_farms,n_months,0);
profits=j(n_farms,n_months,0);
outputs=j(n_farms,n_months,0);
vprofpercwt=j(n_farms,n_months,0);
real_exit=j(n_farms,n_months,.);
real_stay=j(n_farms,n_months,.);
theory_exit=j(n_farms,n_months,0);
theory_stay=j(n_farms,n_months,0);
stay_acc=j(n_farms,n_months,.);
exit_acc=j(n_farms,n_months,.);
exit_output=j(21,2,.);
stay_output=j(21,2,.);
stay_rule=j(1,2,.);
exit_rule=j(1,2,.);
do cutoff=-10 to 10;
do i=1 to n_farms;
do j=7 to n_months;
if avc[i,j]>0 then stoch_avc[i,j]=(avc[i,j]+u[i,j]);
month_cwt[i,j]=month_lb[i,j]/100;
if month_cwt[i,j]>0 then month_prof[i,j]=(price[i,j]-stoch_avc[i,j])*month_cwt[i,j];
outputs[i,j]=(month_cwt[i,j-1]+month_cwt[i,j-2]+month_cwt[i,j-3]+month_cwt[i,j-4]+month_cwt[i,j-5]+month_cwt[i,j-6]);
profits[i,j]=(month_prof[i,j-1]+month_prof[i,j-2]+month_prof[i,j-3]+month_prof[i,j-4]+month_prof[i,j-5]+month_prof[i,j-6]);
if outputs[i,j]>0 then vprofpercwt[i,j]=profits[i,j]/outputs[i,j];
if month_cwt[i,j]>0 then real_stay[i,j]=1;
if vprofpercwt[i,j]>=cutoff then theory_stay[i,j]=1;
if real_stay[i,j]=. then theory_stay[i,j]=.;
if real_stay[i,j]^=. then stay_acc[i,j]=(real_stay[i,j]=theory_stay[i,j]);
if ((month_cwt[i,j]=0)*(month_cwt[i,j-1]^=0))=1 then real_exit[i,j]=1;
if ((month_cwt[i,j]=0)*(month_cwt[i,j-1]=0))=1 then real_exit[i,j]=.;
if vprofpercwt[i,j]<cutoff then theory_exit[i,j]=1;
if real_exit[i,j]=. then theory_exit[i,j]=.;
if real_exit[i,j]^=. then exit_acc[i,j]=(real_exit[i,j]=theory_exit[i,j]);
end;
end;
r=cutoff+11;
*creates outlier fence for exit;
exit=exit_acc[:];
exit_output[r,1]=exit;
exit_output[r,2]=cutoff;
if exit_output[r,1]<.95 then exit_output[r,1]=.;
exit_max=exit_output[><, ];
if exit_output[r,1]<exit_max then exit_output[r,1]=.;
if exit_output[r,1]>0 then exit_rule[1,1]=exit_output[r,1];
exit_rule[1,2]=exit_output[r,2];
exit_outlier=exit_rule[1,1];
exit_outlier_percent=exit_rule[1,2];
*creates outlier fence for stay;
stay=stay_acc[:];
stay_output[r,1]=stay;
stay_output[r,2]=cutoff;
if stay_output[r,1]<.95 then stay_output[r,1]=.;
stay_max=stay_output[><, ];
if stay_output[r,1]<stay_max then stay_output[r,1]=.;
if stay_output[r,1]>0 then stay_rule[1,1]=stay_output[r,1];
stay_rule[1,2]=stay_output[r,2];
stay_outlier=stay_rule[1,1];
stay_outlier_percent=stay_rule[1,2];
end;
print exit_output stay_output;
print stay_rule exit_rule;
quit;
