eine multivariate, differenzierbare Funktion minimiert wird mit scipy.optimize

Question

Ich versuche, die folgende Funktion mit scipy.optimize zu minimieren:

, deren Steigung ist dies:

(für diejenigen, die daran interessiert sind, dies ist die Wahrscheinlichkeitsfunktion eines Bradley-Terry-Luce-Modells für paarweise Vergleiche. Sehr eng mit der logistischen Regression verbunden.)

Es ist ziemlich klar, dass das Hinzufügen einer Konstante zu allen Parametern den Wert der Funktion nicht ändert. Daher lasse ich \ theta_1 = 0. Hier sind die Umsetzung der Zielfunktionen und die Steigung in Python (Theta wird x hier):

def objective(x): 
    x = np.insert(x, 0, 0.0) 
    tiles = np.tile(x, (len(x), 1)) 
    combs = tiles.T - tiles 
    exps = np.dstack((zeros, combs)) 
    return np.sum(cijs * scipy.misc.logsumexp(exps, axis=2)) 

def gradient(x): 
    zeros = np.zeros(cijs.shape) 
    x = np.insert(x, 0, 0.0) 
    tiles = np.tile(x, (len(x), 1)) 
    combs = tiles - tiles.T 
    one = 1.0/(np.exp(combs) + 1) 
    two = 1.0/(np.exp(combs.T) + 1) 
    mat = (cijs * one) + (cijs.T * two) 
    grad = np.sum(mat, axis=0) 
    return grad[1:] # Don't return the first element 

Hier ist ein Beispiel dafür, was cijs aussehen könnte:

[[ 0 5 1 4 6] 
[ 4 0 2 2 0] 
[ 6 4 0 9 3] 
[ 6 8 3 0 5] 
[10 7 11 4 0]] 

Dies ist der Code, den ich ausführen, um die Minimierung auszuführen:

x0 = numpy.random.random(nb_items - 1) 
# Let's try one algorithm... 
xopt1 = scipy.optimize.fmin_bfgs(objective, x0, fprime=gradient, disp=True) 
# And another one... 
xopt2 = scipy.optimize.fmin_cg(objective, x0, fprime=gradient, disp=True) 

es ist jedoch nicht immer in der ersten Iteration:

Warning: Desired error not necessarily achieved due to precision loss. 
     Current function value: 73.290610 
     Iterations: 0 
     Function evaluations: 38 
     Gradient evaluations: 27 

Ich kann nicht herausfinden, warum es fehlschlägt. Der Fehler wird wegen dieser Zeile angezeigt:

So scheint diese "Wolfe Zeile Suche" nicht gelingen, aber ich habe keine Ahnung, wie Sie von hier aus fortfahren ... Jede Hilfe wird geschätzt!

Answer 1

Als @pv. Als Kommentar darauf hingewiesen, habe ich bei der Berechnung des Gradienten einen Fehler gemacht. Zunächst einmal ist der korrekte (mathematische) Ausdruck für die Steigung meiner Zielfunktion:

(das Minuszeichen bemerken.) Darüber hinaus mein Python-Implementierung, über die Zeichen Fehler völlig falsch war.Hier ist mein aktualisiert Gradient:

def gradient(x): 
    nb_comparisons = cijs + cijs.T 
    x = np.insert(x, 0, 0.0) 
    tiles = np.tile(x, (len(x), 1)) 
    combs = tiles - tiles.T 
    probs = 1.0/(np.exp(combs) + 1) 
    mat = (nb_comparisons * probs) - cijs 
    grad = np.sum(mat, axis=1) 
    return grad[1:] # Don't return the first element. 

Um es zu debuggen, habe ich:

• scipy.optimize.check_grad: zeigten, dass meine Gradientenfunktion wurde Ergebnisse sehr weit weg von einer angenähert (Finite-Differenzen) Gradienten zu erzeugen.
• scipy.optimize.approx_fprime, um eine Vorstellung von den Werten zu erhalten sollte aussehen.
• ein paar handverlesene einfache Beispiele, die bei Bedarf von Hand analysiert werden könnten, und einige Wolfram Alpha-Abfragen zur Überprüfung der Funktionsfähigkeit.

Answer 2

Es scheint, dass Sie es in ein (nicht-lineares) Least-Square Problem verwandeln können. Auf diese Weise müssen Sie Intervalle für jede der Variablen n und die Anzahl der Abtastpunkte für jede Variable definieren, um die Matrix der Koeffizienten zu erstellen.

In diesem Beispiel habe ich die gleiche Anzahl von Punkten verwende und das gleiche Intervall für alle Variablen:

from scipy.optimize import leastsq 
from numpy import exp, linspace, zeros, ones 

n = 4 
npts = 1000 
xs = [linspace(0, 1, npts) for _ in range(n)] 

c = ones(n**2) 

a = zeros((n*npts, n**2)) 
def residual(c): 
    a.fill(0) 
    for i in range(n): 
     for j in range(n): 
      for k in range(npts): 
       a[i+k*n, i*n+j] = 1/(exp(xs[i][k] - xs[j][k]) + 1) 
       a[i+k*n, j*n+i] = 1/(exp(xs[j][k] - xs[i][k]) + 1) 

    return a.dot(c) 

popt, pconv = leastsq(residual, x0=c) 
print(popt.reshape(n, n)) 
#[[ -1.24886411 1.07854552 -2.67212118 1.86334625] 
# [ -7.43330057 2.0935734 37.85989442 1.37005925] 
# [ -3.51761322 -37.49627917 24.90538136 -4.23103535] 
# [ 11.93000731 2.52750715 -14.84822686 1.38834225]] 

---

EDIT: mehr Details über die Koeffizienten oben gebaut Matrix: