topic Re: Solving an equation , non linear in Mathematical Optimization, Discrete-Event Simulation, and OR

Solving an equation , non linear

Thomas_mp — Sat, 20 Jun 2020 20:02:55 GMT

Good morning,

I am trying to solve an equation of the type a1 = c1+ (b1-X)*c1/(1+X) + (b2-X)*c2/(1+X)^2 + (b3-X)*c3/(1+X)^3 +… (bt-X)*ct/X*(1/(1+X)^t-1)

Note that the first (t-1) terms are of the same form (bj-Xj)*cj/(1+X)^j-1 , but the last term is different.

I need to find the value of X ( 0<X<1). In the case of multiple solutions, chose the smaller non-negative solution.

In my case I have 25 variables.

ai
12 bi
12 ci

The data looks like:

Input a b1 b2 b3 b4 b5 b6 b7 b8 b9 b10 b11 b12 c1 c2 c3 c4 c5 c6 c7 c8 c9 c10 c11 c12

cards

84.3 0.33 0.31 0.29 0.28 0.26 0.25 0.24 0.22 0.21 0.19 0.18 0.17 11.1 13.4 16.1 19.1 22.4 26.2 30.4 34.9 39.8 45.1 50.7 56.5

91.5 0.332 0.310 0.291 0.277 0.263 0.249 0.235 0.221 0.207 0.194 0.180 0.166 11.076 13.422 16.074 19.061 22.433 26.201 30.368 34.927 39.859 45.132 50.700 56.502

I can, by trial and error, find the value of X, this is: solve the equation in Excel

84.3= 11.1 +(0.33-X)*11.1/(1+X) + (0.31-X)*13.4/(1+X)^2 + (0.29-X)*16.1/(1+X)^3+ (0.28-X)*19.1/(1+X)^4 + ............. + (0.18-X)*50.7/(1+X)^11 + (0.17-X)*56.5/X*[1/(1+X)^11]

In this case the solution is X=0.068, this X produces a value for a=84.3.

Doing the same, the solution for the second row of data (which starts with a value of a=91.5) is X = 0.0638

But I have several thousands of observations … so I do not how to do this in Excel by trial and error (function IRR does not work in excel because X is both in the numerator and in the denominator).

I am attaching the excel file with the data and the solution.

Thank you !!!

Tomas

Re: Solving an equation , non linear

RobPratt — Tue, 23 Jun 2020 14:18:37 GMT

Here's one way to do it with PROC OPTMODEL by using a DO loop to call the solver once per observation:

%let n = 12;
data indata;
   input a b1-b&n c1-c&n @@;
   datalines;
84.3 0.33  0.31  0.29  0.28  0.26  0.25  0.24  0.22  0.21  0.19  0.18  0.17  
     11.1  13.4  16.1  19.1  22.4  26.2  30.4  34.9  39.8  45.1  50.7  56.5
91.5 0.332 0.310 0.291 0.277 0.263 0.249 0.235 0.221 0.207 0.194 0.180 0.166 
     11.076 13.422 16.074 19.061 22.433 26.201 30.368 34.927 39.859 45.132 50.700 56.502
;

proc optmodel printlevel=0;
   num n = &n;
   set OBS;
   num a {OBS};
   num b {OBS, 1..n};
   num c {OBS, 1..n};
   read data indata into OBS=[_N_]
      a {j in 1..n} <b[_N_,j]=col('b'||j)> {j in 1..n} <c[_N_,j]=col('c'||j)>;

   var X >= 0 <= 1;
   min Z = X;
   num obsThis, aThis, bThis {1..n}, cThis {1..n};
   con IRR:
      aThis = cThis[1] + sum {j in 1..n-1} (bThis[j]-X)*cThis[j]/(1+X)^j + (bThis[n]-X)*cThis[n]/X/(1+X)^(n-1);

   num Xsol {OBS};
   do obsThis = OBS;
      put obsThis=;
      aThis = a[obsThis];
      for {j in 1..n} do;
         bThis[j] = b[obsThis,j];
         cThis[j] = c[obsThis,j];
      end;
      solve;
      Xsol[obsThis] = X.sol;
   end;
   create data outdata from [i] X=Xsol;
quit;

You can also use a COFOR loop to solve these independent problems concurrently. Just replace the DO statement with these two statements:

   cofor {i in OBS} do;
      obsThis = i;

In SAS Viya 3.5, it is even simpler (no explicit looping) by using the groupBy option in the runOptmodel action:

data sascas1.indata;
   set indata;
   i = _N_;
run;

proc cas noqueue;
   source pgm;
      put _BY_LINE_;

      num n = &n;
      num a;
      num b {1..n};
      num c {1..n};
      read data indata into
         a {j in 1..n} <b[j]=col('b'||j)> {j in 1..n} <c[j]=col('c'||j)>;

      var X >= 0 <= 1;
      min Z = X;
      con IRR:
         a = c[1] + sum {j in 1..n-1} (b[j]-X)*c[j]/(1+X)^j + (b[n]-X)*c[n]/X/(1+X)^(n-1);

      solve;
      create data outdata from X;
   endsource;
   
   action optimization.runOptmodel / code=pgm printlevel=0 groupBy={"i"};
   run;
quit;

Re: Solving an equation , non linear

Thomas_mp — Sun, 21 Jun 2020 18:24:17 GMT

Thank you !!!

Re: Solving an equation , non linear

Thomas_mp — Tue, 23 Jun 2020 04:27:34 GMT

Robert, Your code using proc optmodel works. Thank you. I am trying to merge the solutions for X (in outdata) with my original sample (data - indata). I can create a counting variable n_ +1 in indata.
How can I keep this variable n_ in outdata so I can merge outdata and indata by n_ ?

Thank you again
Tomas

Re: Solving an equation , non linear

RobPratt — Tue, 23 Jun 2020 14:14:02 GMT

No need to merge. Just replace the READ DATA statement with this:

read data indata into OBS=[n_]
   a {j in 1..n} <b[n_,j]=col('b'||j)> {j in 1..n} <c[n_,j]=col('c'||j)>;

And then replace the CREATE DATA statement with this:

create data outdata from [n_] X=Xsol
   a {j in 1..n} <col('b'||j)=b[n_,j]> {j in 1..n} <col('c'||j)=c[n_,j]>;

Re: Solving an equation , non linear

Thomas_mp — Fri, 26 Jun 2020 02:40:03 GMT

Thank you again !!!