Search code examples
optimizationpyomo

pyomo matrix product

I would like to use pyomo to solve a multiple linear regression under constraint in pyomo. to do so I have 3 matrices :

  • X (noted tour1 in the following code) My inputs (600x13)(bureaux*t1 in pyomo)
  • Y (noted tour2 in the following code) the matrix I want to predict (6003)(bureauxt2 inpyomo)
  • T (noted transfer In the code) (13x3)(t1*t2 in pyomo)

I would like to do the following

ypred = XT
minimize (ypred-y)**2
subject to 
0<T<1 
and Sum_i(Tij)=1

To that effect, I started the following code

from pyomo.environ import * 

tour1=pd.DataFrame(np.random.random(size=(60,13)),columns=["X"+str(i) for i in range(13)],index=["B"+str(i) for i in range(60)])
tour2=pd.DataFrame(np.random.random(size=(60,3)),columns=["Y"+str(i) for i in range(3)],index=["B"+str(i) for i in range(60)])

def gettour1(model,i,j):
  return tour1.loc[i,j]

def gettour2(model,i,j):
  return tour2.loc[i,j]

def cost(model):
  return sum((sum(model.tour1[i,k] * model.transfer[k,j] for k in model.t1) - model.tour2[i,j] )**2 for i in model.bureaux for j in model.tour2)


model = ConcreteModel()
model.bureaux = Set(initialize=tour1.index.tolist())
model.t1 = Set(initialize=tour1.columns)
model.t2 = Set(initialize=tour2.columns)
model.tour1 = Param(model.bureaux, model.t1,initialize=gettour1)
model.tour2 = Param(model.bureaux, model.t2,initialize=gettour2)
model.transfer = Var(model.t1,model.t2,bounds=[0,1])
model.obj=Objective(rule=cost, sense=minimize)

I unfortunately get an error at this stage :

KeyError: "Index '('X0', 'B0', 'Y0')' is not valid for indexed component 'transfer'"

anyone knows how I can calculate the objective ?

furthermore any help for the constrains would be appreciated :-)

Solution

  • A couple things...

    First, the error you are getting. There is information in that error statement that should help identify the problem. The construction appears to be trying to index transfer with a 3-part index (x, b, y). That clearly is outside of t1 x t2. If you look at the sum equation you have, you are mistakenly using model.tour2 instead of model.t2.

    Also, your bounds parameter needs to be a tuple.

    While building the model, you should be pprint()-ing the model very frequently to look for these types of issues. That only works well if you have "small" data. 60 x 13 may be the normal problem size, but it is a total pain to troubleshoot. So, start with something tiny, maybe 3 x 4. Make a Set, pprint(). Make a Constraint, pprint()... Once the model computes/solves with "tiny" data, just pop in the real stuff.