Finden der besten Kombination (pulp, scipy.optimize?)

[qojote]

Hallo,

Ich möchte eine Anzahl von (numpy-) Arrays gewichtet miteinander addieren und dadurch möglichst "nahe" an gegebene Zielwerte kommen.
Also mal als kleines Beispiel:

A1 * [0, 1, 1, 1] + A2 * [0, 0, 0, 4] + A3 * [0, 0, 2, 5] === [0, 1, 4, 4]

Die Faktoren A1, A2 und A3 dürfen nur Werte zwischen 0 und 1 annehmen. Ferner wird die Güte der Anpassung durch die Summe der Beträge der Differenzen an allen Punkten bestimmt.
An diesem Problem sitze ich nun schon wirklich lange und würde mich wirklich über Anregungen freuen-


Hier meine bisherigen Ansätze:

ILP-Problem mit pulp:

Code: Alles auswählen

ySoll = np.array([0., 1., 4., 4.])
yValues = np.array([[0, 1, 1, 1],
                    [0, 0, 0, 4],
                    [0, 0, 2, 5]
                    ], float)

num = yValues.shape[0]
numYValues = yValues.shape[1]

A = []
for i in range(num):
    A.append(LpVariable("A%s"%i, lowBound = 0, upBound = 1, cat= LpContinuous)) # amplitude variables


prob = LpProblem("combination", LpMinimize)            
prob += lpSum([ySoll[j] - lpSum(A[i]*yValues[i][j] for i in range(num))  for j in range(numYValues)]  )

for j in range(numYValues):
    prob += lpSum(A[i]*yValues[i][j] for i in range(num)) == ySoll[j], "Equal_%s"%j
   

prob.writeLP("combinate.lp")
prob.solve()

print "Anpassung ", value(prob.objective)
resY = sum(A[i].varValue * yValues[i] for i in range(num))
print "Ziel:     ", ySoll
print "Erreicht: ", resY

Dummerweise kann man in der Zielfunktion (condition) zwar Summen, aber keine Beträge verwenden...

Dann habe ich scipy.optimize.minimize versucht.

Code: Alles auswählen

def peval(p, X):
    return np.dot(X, p)

def objective(p, y, X):
    return np.sum(   np.fabs(y - peval(p, X))      )
    #return np.sum(   (y - peval(p, X)) ** 2      )

yValues = np.array([[0, 1, 4, 4],
                    [1, 0, 0, 4],
                    [4, 1, 1, 4],
                    [4, 1, 1, 4]
                    ], float)                  
ySoll = np.array([0, 1, 4, 4])

p0 = np.zeros(4)
cons = {
    'type': 'eq',
    'fun': np.sum,
}

bnds = ((0, None), (0, None), (0, None), (0, None))

res = optimize.minimize(
    objective, p0, args=(ySoll, yValues),
    method='SLSQP', constraints=cons)

Ohne "Bounds" bekomme ich Lösungen für die Faktoren, jedoch sind die auch negativ. Wenn ich die "Bounds" angebe, wird keine Lösung mehr gefunden.

Herzlichen Dank!

[MagBen]

Du musst dem Optimierer mehr Hilfestellung geben.
Zuerst minimierst Du die Funktion ohne Grenzen, bei mir kommt da
[[0.9537387161503309, 0.20044877439233902, -0.061677180391656775, 0.011564479711942107], 3.663842887215857e-05, 17, 0, 'Optimization terminated successfully.']
raus. Als nächstes suchst Du eine Lösung in der Nähe dieser ersten Lösung aber mit Grenzen, d.h. Du setzt

Code: Alles auswählen

p0 = [0.9537387161503309, 0.20044877439233902, 0, 0.011564479711942107]
bnds = ((0, 2), (0, 2), (0, 2), (0, 2))

Dabei kommt bei mir
[[1.0000000127997202, 4.1984186964315607e-13, 2.4936470880040453e-14, 1.015996387292698e-13], 1.1520051740800775e-07, 9, 0, 'Optimization terminated successfully.']
raus, also ca. [1, 0, 0, 0].
Was ja auch Plausibel ist, wenn ySoll == yValues[0] gilt.

[qojote]

Hallo,

Also wenn ich das Skript ausführe erhalte ich als Ergebnis ohne Grenzen:
status: 0
success: True
njev: 25
nfev: 188
fun: 2.1458641917382337e-05
x: array([ 1.44445079, -1.9259337 , 0.5926002 , -0.11111729])
message: 'Optimization terminated successfully.'
jac: array([ -9., -3., -6., -16., 0.])
nit: 25

Und wenn ich dann p0 = res.x setze und nochmal mit bounds starte:

Code: Alles auswählen

p0 = res.x
res = optimize.minimize(
    objective, p0, args=(ySoll, yValues),
    method='SLSQP', bounds=bnds, constraints=cons)

kommt das gleiche Ergebnis raus... Mache ich da was falsch?

[MagBen]

Der neue Startwert muss innerhalb der Grenzen sein, deshalb solltest Du die negativen Werte der Lösung auf Null setzen. Da in meiner Scipy Version minimize noch nicht drin enthalten ist, habe ich fmin_slsqp direkt aufgerufen:

Code: Alles auswählen

import numpy as np
from scipy import optimize

def peval(p, X):
    return np.dot(X, p)
 
def objective(p, y, X):
    return np.sum(   np.fabs(y - peval(p, X))      )
 
yValues = np.array([[0, 1, 4, 4],
                    [1, 0, 0, 4],
                    [4, 1, 1, 4],
                    [4, 1, 1, 4]
                    ], float)                  
ySoll = np.array([0, 1, 4, 4])
 
p1 = np.zeros(4)
res1 = optimize.fmin_slsqp(func=objective,
                           x0=p1, 
                           full_output=True, args=(ySoll, yValues))
print( res1)
#[[0.9537387161503309, 0.20044877439233902, -0.061677180391656775, 0.011564479711942107], 3.663842887215857e-05, 17, 0, 'Optimization terminated successfully.']

p2 =res1[0]
p2[p2<0]=0 # <======== negative Werte auf Null setzen
bnds = ((0, 2), (0, 2), (0, 2), (0, 2)) # sinnvolle Grenzen
res2 = optimize.fmin_slsqp(func=objective,
                           x0=p2, 
                           bounds=bnds, 
                           full_output=True, args=(ySoll, yValues))
print( res2)
#[[0.99999992612748334, 6.1553294301239063e-16, -3.4682336035881025e-16, 9.9776029780057904e-17], 6.648526489745999e-07, 8, 0, 'Optimization terminated successfully.']

[qojote]

Vielen Dank, das klappt! Super.
Es funktioniert auch gut, wenn man gleich beim ersten Aufruf die Bounds mit angibt.

Ergebnis:

Optimization terminated successfully. (Exit mode 0)
Current function value: 1.79926309743e-06
Iterations: 16
Function evaluations: 124
Gradient evaluations: 16
(array([ 9.99999743e-01, -7.70690451e-15, -1.15153596e-14,
1.37575597e-07]), 1.7992630974339915e-06, 16, 0, 'Optimization terminated successfully.')

Rein interessenshalber: Ist scipy.optimize.fmin_slsqp das gleiche wie scipy.optimize.minimize(..., method='SLSQP', ...)?

[MagBen]

Rein interessenshalber: Ist scipy.optimize.fmin_slsqp das gleiche wie scipy.optimize.minimize(..., method='SLSQP', ...)?

Ja, beides führt zum Aufruf der Funktion _minimize_slsqp.