BookmarkSubscribeRSS Feed
dpullarkat
Calcite | Level 5

I have a Proc Optmodel code which runs (on SAS 9.3), but the solved values of the variables which I have are beyond the upper and lower bounds I have specified in my var declaration. 

Any idea why this could be happening and what I can do to fix it?

 

As you can see from the code below, ps can been bounded between 0 and 1. But in the solved output, the values frequently go negative or above 1. I am having a similar issue with compliance, which also goes above 1 or negative.

 

PROC OPTMODEL; 

   SET <STRING> DRUG;
   SET <STRING> REGIMEN;   
   SET <NUMBER> MONTH;

   read data month into month=[month];   

   NUMBER dailydose{DRUG};
   READ DATA DRUG_consumption INTO DRUG=[DRUG] dailydose;

   NUMBER DOT{regimen};
   READ DATA REGIMEN_dot INTO regimen=[regimen] dot;    
   PRINT DOT;

   NUMBER UNITS{DRUG,MONTH};
   READ DATA units_data INTO [DRUG MONTH] UNITS;
   PRINT UNITS;

   number weight{regimen, drug};
   read data regimen_drug into [regimen drug] weight; 
  
	/*   identify if it is a partial year wrt to IMS data - check if there are 12 months of data*/
	num month_count=card(MONTH);
	print month_count;

	number launch_date{regimen};
	read data regimen_Launch_dates into [regimen] launch_date;
	print launch_date;

   /*DEFINE THE PARAMETERS*/

   var ps {REGIMEN, MONTH} init 0.3 >=0.00 <=.900,
   annualpatients init 10000 <=25000 >=1000, 
   compliance {drug} >=0.3 <=.95,
   dd{drug} >=4 <=6; 

   var factor >=0.2 <=1;

	for {d in drug} fix dd[d]=dailydose[d];
	for {d in drug} unfix dd["Ribavarin"];

	if month_count=12 then fix factor=1; 
	else unfix factor;


   impvar su_add{ d in drug, M IN MONTH}=units[d,m]/dd[d];   

   impvar monthly_patient_share{M IN MONTH}=sum{r in regimen}ps[r,M];
   impvar regimen_tot_ps{r in regimen}=sum{m in month} ps[r,M];
   impvar annual_ps{r in regimen}=regimen_tot_ps[r]/month_count;
   impvar avg_dot_1 {r in regimen}=annual_ps[r]*dot[r];
   impvar avg_dot_2=sum{r in regimen} avg_dot_1[r];
   impvar avg_dot=avg_dot_2;
   impvar avg_dot_months=avg_dot_2/(30/7);

/*   annual total treated patients*/
   impvar pat_tot=annualpatients*avg_dot_months;
   impvar total_patients=sum{r in regimen, m in month} (ps[r,M]*annualpatients/(avg_dot/(30/7)));
   impvar total_units=sum{m in month,d in drug} su_add[d,M];
   impvar monthly_units{m in month}=sum{d in drug} su_add[d,M];
   impvar ratio{m in month}=monthly_units[m]/total_units;
   impvar monthlypatients{m in month}=pat_tot*ratio[m]*factor;
   
/*   calculate total drug->regimen shares*/
   impvar drug_share{d in drug,M IN MONTH} = sum{r in regimen}ps[r,m]*weight[r,d];
/*   calculated units*/
    impvar calc_units{d in drug,M IN MONTH}=monthlypatients[m]*dd[d]*compliance[d]*30*drug_share[d,m];
/*	deviation*/
	impvar deviation{d in drug, M IN MONTH}=(units[D,M]-calc_units[d,m])/calc_units[d,m];

   /*SET THE CONSTRAINTS*/
   constraint  MONTHLY_ps {M IN MONTH}:  monthly_patient_share[m]=1;
   constraint  patient_match: total_patients=annualpatients;   
   constraint calc_units_dev{d in drug,M IN MONTH}: -0.05 <= deviation[d,m] <= 0.05;


   /*SET THE OBJECTIVE*/
   var error_m_d {d in drug,M IN MONTH};
   CON error_con {d in drug,M IN MONTH}: error_m_d[d,m] = (units[D,M]-calc_units[d,m])/dd[d];
   Min error  = sum{d in drug,M IN MONTH}error_m_d[D,M]**2;

/*min error{d in drug, m in month} = deviation[d,m];*/

/* expand;*/
solve with nlp;
create data pat_solver_shares from [regimen month] ps;
create data matching from [drug month] units calc_units deviation error_m_d; 
create data compliance from [drug] compliance DD;
create data dot from [regimen]dot;
create data annual_pats from annualpatients;
create data factor from factor;
 quit;
7 REPLIES 7
RobPratt
SAS Super FREQ

Please also include the input data.

dpullarkat
Calcite | Level 5

Hi Robb,

 

Below is the code, and attached are the input data file and also the output data from the optmodel procedure.

 

Thanks,

Deepthi

 

%let folder=C:\Users\pullarkd\Documents\Work\Solvers SAS Trial\;
%let input_file=Masked Data.xlsx;

libname solver "&folder.";

/*read in data from excel input template*/

LIBNAME inputs EXCEL PATH="&folder.\&input_file.";

data units_data_raw;
set inputs.'Inputs-Units Data$'n;
run;

data drug_consumption_master;
set inputs.drugdosing;
where dailydose is not null;
run;

data regimen_dot_master;
set inputs.regimen_dot_master;
where dot is not null;
run;

data regimen_drug_master;
set inputs.mapping;
where weight is not null;
run;

data compliance_constraints;
set inputs.compliance_range;
where upper is not null;
run;

data regimen_key_drug;
set inputs.regimen_key_drug;
where drug is not null;
run;

LIBNAME inputs CLEAR;

proc sql;
create table units_data_3 as
select country, period as month, drug, sum(value) as units 
from units_data_raw
group by 1,2,3;
run;

proc sql;
create table month as
select distinct month
from units_data_3;
quit;

proc sql;
create table master as
select drug, month
from (select distinct drug from units_data_3) as a, month;
quit;

proc sql;
create table units_data as
select b.country, a.drug, a.month, case when units is null then 0 else units end as units
from master as a
left join units_data_3 as b
on a.drug=b.drug and a.month=b.month;
quit;

proc sql;
create table drug_consumption as
select * from drug_consumption_master
where drug in (select distinct drug from units_data);
quit;

proc sql;
create table regimen_drug as
select regimen, drug, weight 
from regimen_drug_master  
where regimen in
(select distinct regimen from
regimen_key_drug
where drug in (select distinct drug from units_data));
quit;

proc sql;
create table regimen_dot as 
select distinct regimen, dot 
from regimen_dot_master
where regimen in
(select distinct regimen from
regimen_key_drug
where drug in (select distinct drug from units_data));
quit;

proc sql;
create table launch_dates as
select drug, min(month) as launch_date 
from
(select drug, month, units
from units_data where units >0) as a
group by 1;
quit;

proc sql;
create table regimen_launch_dates as
select regimen, min(launch_date) as launch_date 
from
(select b.regimen, b.drug, a.launch_date, b.weight
from launch_dates as a
left join regimen_drug as b
on a.drug=b.drug) as a
where weight>0 and 
(drug in ("Drug3", "Drug4", "Drug5", "Drug6", "Drug7", "Drug8", "Drug9") or regimen="Reg1")
group by 1;
quit;

PROC OPTMODEL; 

   SET <STRING> DRUG;
   SET <STRING> REGIMEN;   
   SET <NUMBER> MONTH;

   read data month into month=[month];   

   NUMBER dailydose{DRUG};
   READ DATA DRUG_consumption INTO DRUG=[DRUG] dailydose;

   NUMBER DOT{regimen};
   READ DATA REGIMEN_dot INTO regimen=[regimen] dot;    
   PRINT DOT;

   NUMBER UNITS{DRUG,MONTH};
   READ DATA units_data INTO [DRUG MONTH] UNITS;
   PRINT UNITS;

   number weight{regimen, drug};
   read data regimen_drug into [regimen drug] weight; 
  
	/*identify if it is a partial year wrt to IMS data - check if there are 12 months of data*/
	num month_count=card(MONTH);
	print month_count;

	number launch_date{regimen};
	read data regimen_Launch_dates into [regimen] launch_date;	

   /*DEFINE THE PARAMETERS*/

   var ps {REGIMEN, MONTH} init 0.3 >=0.00 <=.900,
   annualpatients init 10000 <=25000 >=1000,    
   dd{drug} >=4 <=6; 

   var factor >=0.2 <=1;

   var compliance{drug};
   read data compliance_constraints into [drug] compliance.lb=lower compliance.ub=upper;

	for {d in drug} fix dd[d]=dailydose[d];
	for {d in drug} unfix dd["Drug2"];

	FOR {R IN REGIMEN,M IN MONTH} IF M<LAUNCH_DATE[r] THEN FIX PS[R,M]=0;

	if month_count=12 then fix factor=1; 
	else unfix factor;


   impvar su_add{ d in drug, M IN MONTH}=units[d,m]/dd[d];   

   impvar monthly_patient_share{M IN MONTH}=sum{r in regimen}ps[r,M];
   impvar regimen_tot_ps{r in regimen}=sum{m in month} ps[r,M];
   impvar annual_ps{r in regimen}=regimen_tot_ps[r]/month_count;
   impvar avg_dot_1 {r in regimen}=annual_ps[r]*dot[r];
   impvar avg_dot_2=sum{r in regimen} avg_dot_1[r];
   impvar avg_dot=avg_dot_2;
   impvar avg_dot_months=avg_dot_2/(30/7);

/*   annual total treated patients*/
   impvar pat_tot=annualpatients*avg_dot_months;
   impvar total_patients=sum{r in regimen, m in month} (ps[r,M]*annualpatients/(avg_dot/(30/7)));
   impvar total_units=sum{m in month,d in drug} su_add[d,M];
   impvar monthly_units{m in month}=sum{d in drug} su_add[d,M];
   impvar ratio{m in month}=monthly_units[m]/total_units;
   impvar monthlypatients{m in month}=pat_tot*ratio[m]*factor;
   
/*   calculate total drug->regimen shares*/
   impvar drug_share{d in drug,M IN MONTH} = sum{r in regimen}ps[r,m]*weight[r,d];
/*   calculated units*/
    impvar calc_units{d in drug,M IN MONTH}=monthlypatients[m]*dd[d]*compliance[d]*30*drug_share[d,m];
/*	deviation*/
	impvar deviation{d in drug, M IN MONTH}=(units[D,M]-calc_units[d,m])/calc_units[d,m];

   /*SET THE CONSTRAINTS*/
   constraint  MONTHLY_ps {M IN MONTH}:  monthly_patient_share[m]=1;
   constraint  patient_match: total_patients=annualpatients;   
   constraint calc_units_dev{d in drug,M IN MONTH}: -0.05 <= deviation[d,m] <= 0.05;

   /*SET THE OBJECTIVE*/
   var error_m_d {d in drug,M IN MONTH};
   CON error_con {d in drug,M IN MONTH}: error_m_d[d,m] = (units[D,M]-calc_units[d,m])/dd[d];
   Min error  = sum{d in drug,M IN MONTH}error_m_d[D,M]**2;

/*min error{d in drug, m in month} = deviation[d,m];*/

/* expand;*/
solve with nlp/maxiter=100;
create data pat_solver_shares from [regimen month] ps;
create data matching from [drug month] units calc_units deviation error_m_d; 
create data compliance from [drug] compliance DD;
create data dot from [regimen]dot;
create data annual_pats from annualpatients;
create data factor from factor;
 quit;
RobPratt
SAS Super FREQ

You have stopped the solver early by setting MAXITER=100, and the resulting Infeasibility value is positive, meaning that some bounds or constraints are violated.  Are you sure that the problem has a feasible solution?

dpullarkat
Calcite | Level 5

Hi Robb,

 

Sorry, I had forgotten to remove the MAXITER option from the code I shared (it was for testing purposes).

 

But the issue is still present even if the optmodel ran 5000 iterations, which is the maximum limit (or can I increase it using the MAXITER?).

 

The problem does not have a 'fully' feasible solution in all scenarios, but it is more important that the variables are within the constraints, even if the objective value is somewhat large. How does the proc optmodel algorithm work? Does is use all possible values of the variables, even those outside the constraints?

 

If I reduce the scale of the objective, by mutiplying with e.g. 10e-4, would it work?

 

Thanks,

Deepthi

 

 

wezhou
SAS Employee

There are two main SAS nlp solvers available. One is IP (interior point solver), and the other is AS (active set solver). Currently, constraints, including variable bounds are relaxed, in solver IP, in order to permit larger steps to enhance the performance. When a feasible point is not found, some of the constraints can be violated, including bound constraints. On the other hand, solver AS will honor variable bounds. So if you prefer the algorithm stays inside the box bounds, for now, you might should use solver AS.

dpullarkat
Calcite | Level 5

Thanks for helping out!

 

Is there any good material that can help me in understanding the mathematical basis for this? I was not able to find the information I wanted from support.sas

 

Thanks again,

Deepthi

sas-innovate-2024.png

Join us for SAS Innovate April 16-19 at the Aria in Las Vegas. Bring the team and save big with our group pricing for a limited time only.

Pre-conference courses and tutorials are filling up fast and are always a sellout. Register today to reserve your seat.

 

Register now!

Multiple Linear Regression in SAS

Learn how to run multiple linear regression models with and without interactions, presented by SAS user Alex Chaplin.

Find more tutorials on the SAS Users YouTube channel.

Discussion stats
  • 7 replies
  • 1873 views
  • 0 likes
  • 3 in conversation