Proc Optmodel - saving solutions as dataset.

[DNix1]

Self contained code.  First part generates data and is working.  Last part is for Proc Optmodel.   The model is running on array of 10 rows.  Goal is to solve row by row and save solutions.  Once working, the rows will increase to 1000 or more.  

Focus of problem is saving data from Proc Optmodel

Output needed is ki[i] , V[i] with minimization of objective function f.

ki involved iterative guesses until solution is obtained. 

Vi is an explicit Value determined for each iteration of ki.

the objective function depends on ki and Vi for each iteration.  

I'm getter correct solutions for each row, but I cannot save the output or pull up the array for solutions.  I have tried numerous approaches but always get errors. 

Code below. 

options nodate ;

proc datasets lib=work kill noprint;

run;

%let n=10;

%let fb=0.89;

%let CLcr = 100;

%let ka=2.202;

%let clnr=0.364;

%let weight=70;

%let dose=400;

data covar; *Create listing of covariates;

  do i=1 to &n;

  num=i;

  weight=50+rand("uniform")*50; *Weight range 50-100 kg;

  CLCR=20+rand("uniform")*10; *CLcr uniform distribtuion 20-30 ml/min;

  output;

  end;

  run;

*/;

data param;

  weight=&weight;

  clcr=&clcr;

  * clearance renal;

  CLr=CLCR*0.643/90; *covariate adjusted CLr;

  vCLr=(1.68*CLr)**2; *1.68 is the RSD for CLr;

  CLrphi=sqrt(vCLr+CLr**2); *Phi - used to calculate log parameters;

  LogCLrmu=log(CLr**2/clrphi); *Mean for simulations;

  logclrsig=sqrt(log(clrphi**2/CLr**2)); *sd for logCLtmu;

  * V;

  V=48.226/72*weight;

  vV=(0.696*V)**2; * 0.696 is the cov;

  Vphi=sqrt(vV+V**2);

  logVmu=log(V**2/Vphi); *Mean for log of V;

  logVsig=sqrt(log(Vphi**2/V**2)); *sd for log of Vc

  run;

data s1;  

    set param;

length _TYPE_ $4 _NAME_ $7 ;

    _TYPE_='MEAN'; _NAME_="logCLr";

    logclr=logclrmu;

logv=logvmu;

keep _TYPE_ _NAME_ logclr logv;

data s4;

set param;

length _TYPE_ $4 _NAME_ $7 ;

    _TYPE_='COV'; _NAME_="logclr";

    logclr=logclrsig**2;

logv=0.5*logclrsig*logVsig;

keep _TYPE_ _NAME_ logclr logv;

data s5;

    set param;

length _TYPE_ $4 _NAME_ $7 ;

    _TYPE_='COV'; _NAME_="logv";

    logclr=0.5*logclrsig*logVsig ;

logv=logvsig**2;

keep _TYPE_ _NAME_ logclr logv;

run;

data s (type=cov);

    set s4 s5 s1;

run;

Title2 'Parmeter Means and Correlation Matrix';

proc simnormal data=s outsim=SimMVN 

               numreal = &n           /* number of realizations = size of sample */

               seed = 12345;          /* random number seed */

   var logclr logv ;

   run;

data f_parm;

    set SimMVN ;

CLr=exp(logclr);

V=exp(logv);

CLs=clr+&clnr;

ke=CLs/V;

Dose=&dose;

ka=&ka;

run;

title2 'Patameter Values';

proc means data=f_parm;

    var cls v;

run;

proc univariate data=f_parm;

    var cls v;

run;

proc sgplot data=flu_parm;

    histogram cls;

run;

proc sgplot data=flu_parm;

    histogram v;

run;

data sim2;

   set SimMVN ;

   CLr=exp(logclr);

V=exp(logv);

CLs=clr+&clnr;

ke=CLs/V;

Dose=&dose;

ka=&ka;

length Condition $12;

    Condition='400 mg Q 24 h';

  dose=400; tau=24;

    if CLs < 0.388 or CLs > 3.14 or V < 6.6 or v > 167 then delete;

   C_0=dose*ka/(v*(ka-ke))*(exp(-ke*0)/(1-exp(-ke*24))-exp(-ka*0)/(1-exp(-ka*24)));

   C_2=dose*ka/(v*(ka-ke))*(exp(-ke*2)/(1-exp(-ke*24))-exp(-ka*2)/(1-exp(-ka*24)));

   C_4=dose*ka/(v*(ka-ke))*(exp(-ke*4)/(1-exp(-ke*24))-exp(-ka*4)/(1-exp(-ka*24)));

   C_6=dose*ka/(v*(ka-ke))*(exp(-ke*6)/(1-exp(-ke*24))-exp(-ka*6)/(1-exp(-ka*24)));

   C_8=dose*ka/(v*(ka-ke))*(exp(-ke*8)/(1-exp(-ke*24))-exp(-ka*8)/(1-exp(-ka*24)));

   C_10=dose*ka/(v*(ka-ke))*(exp(-ke*10)/(1-exp(-ke*24))-exp(-ka*10)/(1-exp(-ka*24)));

   C_12=dose*ka/(v*(ka-ke))*(exp(-ke*12)/(1-exp(-ke*24))-exp(-ka*12)/(1-exp(-ka*24)));

   C_16=dose*ka/(v*(ka-ke))*(exp(-ke*16)/(1-exp(-ke*24))-exp(-ka*16)/(1-exp(-ka*24)));

   C_20=dose*ka/(v*(ka-ke))*(exp(-ke*20)/(1-exp(-ke*24))-exp(-ka*20)/(1-exp(-ka*24)));

   C_22=dose*ka/(v*(ka-ke))*(exp(-ke*22)/(1-exp(-ke*24))-exp(-ka*22)/(1-exp(-ka*24)));

   C_23=dose*ka/(v*(ka-ke))*(exp(-ke*23)/(1-exp(-ke*24))-exp(-ka*23)/(1-exp(-ka*24)));

   C_24=dose*ka/(v*(ka-ke))*(exp(-ke*24)/(1-exp(-ke*24))-exp(-ka*24)/(1-exp(-ka*24)));

  AUC24_sim=dose/cls;

  format c_0 best4. c_2 best4. c_4 best4. c_6 best4. c_8 best4. c_10 best4. c_12 best4. c_16 best4. c_20 best4. c_22 best4. c_23 best4. c_24 best4.;

  format auc24_sim best4.;

  run;

   proc means data=sim2 n median min max;

      var C_0 c_2 c_4 c_6 c_8 c_10 c_12 c_16 c_20 c_22 c_23 c_24;

  run;

   proc corr data=sim2;

      var cls v;

  run;

 proc optmodel;

       set OBS;

   num init {OBS};

   num ke {OBS};

   num V {OBS};

   num ka {OBS};

   num dose {OBS};

   num tau{OBS};

   num c_23{OBS};

   num auc24_sim {OBS};

   read data sim2 into OBS=[_N_] ke V ka tau c_23 auc24_sim;

   num simno;

   print ke v ka tau c_23 auc24_sim ;

   var ki init 0.01;

   var Vi;

       Vi = Dose[simno]*ka[simno]/(C_23[simno]*(ka[simno]-ki))*(exp(-ki*23)/(1-exp(-ki*tau[simno]))-exp(-ka[simno]*23)/(1-exp(-ka[simno]*tau[simno])));

   min f1=(ke[simno]-ki)**2/(ke[simno]*0.5)**2 + (V[simno]-Vi)**2/(V[simno]*0.5)**2;

   do simno=obs;

      put simno;

  solve;

  print simno ki vi;

  end;

*create data outki from [simno]=test;

quit;

[Ksharp · Accepted solution]

If you want to BY group processing with PROC OPTMODEL:
https://support.sas.com/kb/42/332.html

Or using PROC APPEND as I showed here:
https://communities.sas.com/t5/Mathematical-Optimization/Redistributing-numbers/m-p/963804#M4302

[DNix1]

Thank very much.  This works like a charm, and is greatly appreciated after two weeks of struggling. 

[Ksharp]

OK. Try this one.

But I noticed you have some error when running your code.

1874! Dose[simno]*ka[simno]/(C_23[simno]*(ka[simno]-ki))*(exp(-ki*23)/(1-exp(-ki*tau[simno]))-exp(-ka[simno]*23)/(1-exp(-ka[sim
1874! no]*tau[simno])));
ERROR: 符号“simno”在第 1874 行第 18 列没有值。

options nodate ;
 
proc datasets lib=work kill noprint;
run;
 
%let n=10;
%let fb=0.89;
%let CLcr = 100;
%let ka=2.202;
%let clnr=0.364;
%let weight=70;
%let dose=400;
 
data covar; *Create listing of covariates;
  do i=1 to &n;
  num=i;
  weight=50+rand("uniform")*50; *Weight range 50-100 kg;
  CLCR=20+rand("uniform")*10; *CLcr uniform distribtuion 20-30 ml/min;
  output;
  end;
  run;
*/;
 
 
data param;
  
  weight=&weight;
  clcr=&clcr;
  
  * clearance renal;
 
  CLr=CLCR*0.643/90; *covariate adjusted CLr;
  vCLr=(1.68*CLr)**2; *1.68 is the RSD for CLr;
  CLrphi=sqrt(vCLr+CLr**2); *Phi - used to calculate log parameters;
  LogCLrmu=log(CLr**2/clrphi); *Mean for simulations;
  logclrsig=sqrt(log(clrphi**2/CLr**2)); *sd for logCLtmu;
 
  * V;
  V=48.226/72*weight;
  vV=(0.696*V)**2; * 0.696 is the cov;
  Vphi=sqrt(vV+V**2);
  logVmu=log(V**2/Vphi); *Mean for log of V;
  logVsig=sqrt(log(Vphi**2/V**2)); *sd for log of Vc
  run;
 
data s1;  
    set param;
length _TYPE_ $4 _NAME_ $7 ;
    _TYPE_='MEAN'; _NAME_="logCLr";
    logclr=logclrmu;
logv=logvmu;
keep _TYPE_ _NAME_ logclr logv;
 
data s4;
set param;
length _TYPE_ $4 _NAME_ $7 ;
    _TYPE_='COV'; _NAME_="logclr";
    logclr=logclrsig**2;
logv=0.5*logclrsig*logVsig;
keep _TYPE_ _NAME_ logclr logv;
 
data s5;
    set param;
length _TYPE_ $4 _NAME_ $7 ;
    _TYPE_='COV'; _NAME_="logv";
    logclr=0.5*logclrsig*logVsig ;
logv=logvsig**2;
keep _TYPE_ _NAME_ logclr logv;
run;
 
data s (type=cov);
    set s4 s5 s1;
run;
 
Title2 'Parmeter Means and Correlation Matrix';
 
proc simnormal data=s outsim=SimMVN 
               numreal = &n           /* number of realizations = size of sample */
               seed = 12345;          /* random number seed */
   var logclr logv ;
   run;
 
 
data f_parm;
    set SimMVN ;
CLr=exp(logclr);
V=exp(logv);
CLs=clr+&clnr;
ke=CLs/V;
Dose=&dose;
ka=&ka;
run;
 
title2 'Patameter Values';
 
proc means data=f_parm;
    var cls v;
run;
 
proc univariate data=f_parm;
    var cls v;
run;
 
proc sgplot data=flu_parm;
    histogram cls;
run;
 
proc sgplot data=flu_parm;
    histogram v;
run;
 
data sim2;
   set SimMVN ;
   CLr=exp(logclr);
V=exp(logv);
CLs=clr+&clnr;
ke=CLs/V;
Dose=&dose;
ka=&ka;
length Condition $12;
    Condition='400 mg Q 24 h';
  dose=400; tau=24;
    if CLs < 0.388 or CLs > 3.14 or V < 6.6 or v > 167 then delete;
   C_0=dose*ka/(v*(ka-ke))*(exp(-ke*0)/(1-exp(-ke*24))-exp(-ka*0)/(1-exp(-ka*24)));
   C_2=dose*ka/(v*(ka-ke))*(exp(-ke*2)/(1-exp(-ke*24))-exp(-ka*2)/(1-exp(-ka*24)));
   C_4=dose*ka/(v*(ka-ke))*(exp(-ke*4)/(1-exp(-ke*24))-exp(-ka*4)/(1-exp(-ka*24)));
   C_6=dose*ka/(v*(ka-ke))*(exp(-ke*6)/(1-exp(-ke*24))-exp(-ka*6)/(1-exp(-ka*24)));
   C_8=dose*ka/(v*(ka-ke))*(exp(-ke*8)/(1-exp(-ke*24))-exp(-ka*8)/(1-exp(-ka*24)));
   C_10=dose*ka/(v*(ka-ke))*(exp(-ke*10)/(1-exp(-ke*24))-exp(-ka*10)/(1-exp(-ka*24)));
   C_12=dose*ka/(v*(ka-ke))*(exp(-ke*12)/(1-exp(-ke*24))-exp(-ka*12)/(1-exp(-ka*24)));
   C_16=dose*ka/(v*(ka-ke))*(exp(-ke*16)/(1-exp(-ke*24))-exp(-ka*16)/(1-exp(-ka*24)));
   C_20=dose*ka/(v*(ka-ke))*(exp(-ke*20)/(1-exp(-ke*24))-exp(-ka*20)/(1-exp(-ka*24)));
   C_22=dose*ka/(v*(ka-ke))*(exp(-ke*22)/(1-exp(-ke*24))-exp(-ka*22)/(1-exp(-ka*24)));
   C_23=dose*ka/(v*(ka-ke))*(exp(-ke*23)/(1-exp(-ke*24))-exp(-ka*23)/(1-exp(-ka*24)));
   C_24=dose*ka/(v*(ka-ke))*(exp(-ke*24)/(1-exp(-ke*24))-exp(-ka*24)/(1-exp(-ka*24)));
 
  AUC24_sim=dose/cls;
 
  format c_0 best4. c_2 best4. c_4 best4. c_6 best4. c_8 best4. c_10 best4. c_12 best4. c_16 best4. c_20 best4. c_22 best4. c_23 best4. c_24 best4.;
  format auc24_sim best4.;
  run;
 
 
   proc means data=sim2 n median min max;
      var C_0 c_2 c_4 c_6 c_8 c_10 c_12 c_16 c_20 c_22 c_23 c_24;
  run;
 
   proc corr data=sim2;
      var cls v;
  run;
 
 
 proc optmodel;
       set OBS;
   num init {OBS};
   num ke {OBS};
   num V {OBS};
   num ka {OBS};
   num dose {OBS};
   num tau{OBS};
   num c_23{OBS};
   num auc24_sim {OBS};
   read data sim2 into OBS=[_N_] ke V ka  dose tau c_23 auc24_sim;
 
   num simno;
          
   print ke v ka tau c_23 auc24_sim ;
   
   var ki init 0.01;
   var Vi;
       Vi = Dose[simno]*ka[simno]/(C_23[simno]*(ka[simno]-ki))*(exp(-ki*23)/(1-exp(-ki*tau[simno]))-exp(-ka[simno]*23)/(1-exp(-ka[simno]*tau[simno])));
   min f1=(ke[simno]-ki)**2/(ke[simno]*0.5)**2 + (V[simno]-Vi)**2/(V[simno]*0.5)**2;
   
 

 /* Declare parameters to store solution values */
         num ki_sol {OBS};
         num vi_sol {OBS};

cofor {b in OBS} do;
  simno = b;
  solve;
  ki_sol[b] = ki.sol;
  vi_sol[b] = vi.sol;
end;
 
create data outki from [simno]=OBS ki_sol vi_sol;

quit;
   

[Ksharp]

I think you should use IMPAR statement to get rid of that ERROR information:

options nodate ;
 
proc datasets lib=work kill noprint;
run;
 
%let n=10;
%let fb=0.89;
%let CLcr = 100;
%let ka=2.202;
%let clnr=0.364;
%let weight=70;
%let dose=400;
 
data covar; *Create listing of covariates;
  do i=1 to &n;
  num=i;
  weight=50+rand("uniform")*50; *Weight range 50-100 kg;
  CLCR=20+rand("uniform")*10; *CLcr uniform distribtuion 20-30 ml/min;
  output;
  end;
  run;
*/;
 
 
data param;
  
  weight=&weight;
  clcr=&clcr;
  
  * clearance renal;
 
  CLr=CLCR*0.643/90; *covariate adjusted CLr;
  vCLr=(1.68*CLr)**2; *1.68 is the RSD for CLr;
  CLrphi=sqrt(vCLr+CLr**2); *Phi - used to calculate log parameters;
  LogCLrmu=log(CLr**2/clrphi); *Mean for simulations;
  logclrsig=sqrt(log(clrphi**2/CLr**2)); *sd for logCLtmu;
 
  * V;
  V=48.226/72*weight;
  vV=(0.696*V)**2; * 0.696 is the cov;
  Vphi=sqrt(vV+V**2);
  logVmu=log(V**2/Vphi); *Mean for log of V;
  logVsig=sqrt(log(Vphi**2/V**2)); *sd for log of Vc
  run;
 
data s1;  
    set param;
length _TYPE_ $4 _NAME_ $7 ;
    _TYPE_='MEAN'; _NAME_="logCLr";
    logclr=logclrmu;
logv=logvmu;
keep _TYPE_ _NAME_ logclr logv;
 
data s4;
set param;
length _TYPE_ $4 _NAME_ $7 ;
    _TYPE_='COV'; _NAME_="logclr";
    logclr=logclrsig**2;
logv=0.5*logclrsig*logVsig;
keep _TYPE_ _NAME_ logclr logv;
 
data s5;
    set param;
length _TYPE_ $4 _NAME_ $7 ;
    _TYPE_='COV'; _NAME_="logv";
    logclr=0.5*logclrsig*logVsig ;
logv=logvsig**2;
keep _TYPE_ _NAME_ logclr logv;
run;
 
data s (type=cov);
    set s4 s5 s1;
run;
 
Title2 'Parmeter Means and Correlation Matrix';
 
proc simnormal data=s outsim=SimMVN 
               numreal = &n           /* number of realizations = size of sample */
               seed = 12345;          /* random number seed */
   var logclr logv ;
   run;
 
 
data f_parm;
    set SimMVN ;
CLr=exp(logclr);
V=exp(logv);
CLs=clr+&clnr;
ke=CLs/V;
Dose=&dose;
ka=&ka;
run;
 
title2 'Patameter Values';
 
proc means data=f_parm;
    var cls v;
run;
 
proc univariate data=f_parm;
    var cls v;
run;
 
proc sgplot data=flu_parm;
    histogram cls;
run;
 
proc sgplot data=flu_parm;
    histogram v;
run;
 
data sim2;
   set SimMVN ;
   CLr=exp(logclr);
V=exp(logv);
CLs=clr+&clnr;
ke=CLs/V;
Dose=&dose;
ka=&ka;
length Condition $12;
    Condition='400 mg Q 24 h';
  dose=400; tau=24;
    if CLs < 0.388 or CLs > 3.14 or V < 6.6 or v > 167 then delete;
   C_0=dose*ka/(v*(ka-ke))*(exp(-ke*0)/(1-exp(-ke*24))-exp(-ka*0)/(1-exp(-ka*24)));
   C_2=dose*ka/(v*(ka-ke))*(exp(-ke*2)/(1-exp(-ke*24))-exp(-ka*2)/(1-exp(-ka*24)));
   C_4=dose*ka/(v*(ka-ke))*(exp(-ke*4)/(1-exp(-ke*24))-exp(-ka*4)/(1-exp(-ka*24)));
   C_6=dose*ka/(v*(ka-ke))*(exp(-ke*6)/(1-exp(-ke*24))-exp(-ka*6)/(1-exp(-ka*24)));
   C_8=dose*ka/(v*(ka-ke))*(exp(-ke*8)/(1-exp(-ke*24))-exp(-ka*8)/(1-exp(-ka*24)));
   C_10=dose*ka/(v*(ka-ke))*(exp(-ke*10)/(1-exp(-ke*24))-exp(-ka*10)/(1-exp(-ka*24)));
   C_12=dose*ka/(v*(ka-ke))*(exp(-ke*12)/(1-exp(-ke*24))-exp(-ka*12)/(1-exp(-ka*24)));
   C_16=dose*ka/(v*(ka-ke))*(exp(-ke*16)/(1-exp(-ke*24))-exp(-ka*16)/(1-exp(-ka*24)));
   C_20=dose*ka/(v*(ka-ke))*(exp(-ke*20)/(1-exp(-ke*24))-exp(-ka*20)/(1-exp(-ka*24)));
   C_22=dose*ka/(v*(ka-ke))*(exp(-ke*22)/(1-exp(-ke*24))-exp(-ka*22)/(1-exp(-ka*24)));
   C_23=dose*ka/(v*(ka-ke))*(exp(-ke*23)/(1-exp(-ke*24))-exp(-ka*23)/(1-exp(-ka*24)));
   C_24=dose*ka/(v*(ka-ke))*(exp(-ke*24)/(1-exp(-ke*24))-exp(-ka*24)/(1-exp(-ka*24)));
 
  AUC24_sim=dose/cls;
 
  format c_0 best4. c_2 best4. c_4 best4. c_6 best4. c_8 best4. c_10 best4. c_12 best4. c_16 best4. c_20 best4. c_22 best4. c_23 best4. c_24 best4.;
  format auc24_sim best4.;
  run;
 
 
   proc means data=sim2 n median min max;
      var C_0 c_2 c_4 c_6 c_8 c_10 c_12 c_16 c_20 c_22 c_23 c_24;
  run;
 
   proc corr data=sim2;
      var cls v;
  run;
 
 
 proc optmodel;
       set OBS;
   num init {OBS};
   num ke {OBS};
   num V {OBS};
   num ka {OBS};
   num dose {OBS};
   num tau{OBS};
   num c_23{OBS};
   num auc24_sim {OBS};
   read data sim2 into OBS=[_N_] ke V ka  dose tau c_23 auc24_sim;
 
   num simno;
          
   print ke v ka tau c_23 auc24_sim ;
   
   var ki init 0.01;
   
   impvar   Vi = Dose[simno]*ka[simno]/(C_23[simno]*(ka[simno]-ki))*(exp(-ki*23)/(1-exp(-ki*tau[simno]))-exp(-ka[simno]*23)/(1-exp(-ka[simno]*tau[simno])));
   min f1=(ke[simno]-ki)**2/(ke[simno]*0.5)**2 + (V[simno]-Vi)**2/(V[simno]*0.5)**2;
   
 

 /* Declare parameters to store solution values */
         num ki_sol {OBS};
         num vi_sol {OBS};

cofor {b in OBS} do;
  simno = b;
  solve;
  ki_sol[b] = ki.sol;
  vi_sol[b] = vi.sol;
end;
 
create data outki from [simno]=OBS ki_sol vi_sol;

quit;