Domanda

Sono nuovo per l'ottimizzazione. Sto cercando di risolvere un problema lineare dei minimi quadrati usando il fmin_slsqp funzione in scipy.optimize.

Ho la funzione oggettiva come Frobenius Norm di |q0_T*P-q1_T| quadrato, dove q0_T è la trasposizione di a nX1 vettore e P è nXn matrice e q1_T è una trasposizione di a nX1 vettore. Questo è fondamentalmente un processo di Markov con q vettori come probabilità di essere in uno stato e P essere la matrice delle probabilità di transizione.

La funzione obiettivo sarebbe ridotta al minimo WRT P Dove sono i vincoli:

1) Tutti gli elementi in P deve essere non negativo

2) Tutte le righe P deve sommarsi a 1

Ho definito questa funzione oggettiva che non sono sicuro sia corretto:

def func(P, q):
  return (np.linalg.norm(np.dot(transpose(q[0,]),P)-transpose(q[1,])))**2

Il secondo argomento in FMIN_SLSQP chiede un 1D NdArray X0 che è l'ipotesi iniziale per la variabile indipendente. Ma qui, poiché la mia variabile indipendente sarebbe P, avrei bisogno di un array 2D per la mia ipotesi iniziale. Non sono sicuro di inquadrare correttamente il problema o dovrò usare un'altra funzione. Grazie.

È stato utile?

Soluzione

Quindi c'è un problema in quanto accetta solo un array 1D come variabile indipendente, una soluzione hackish sarebbe quella di passarla come 1D e rimodellare nella funzione.

import numpy as np
from scipy.optimize import fmin_slsqp

# some fake data
d = 3 # dimensionality of the problem
P0 = np.arange(d*d).reshape(d, d)
P0 /= P0.sum(1, keepdims=True) # so that each row sums to 1
q = np.random.rand(2, d)       # assuming this is the structure of your q
q /= q.sum(1, keepdims=True)

# the function to minimize
def func(P, q):
    n = q.shape[-1] # or n = np.sqrt(P.size)
    P = P.reshape(n, n)
    return np.linalg.norm(np.dot(q[0], P) - q[1])**2  # no changes here, just simplified syntax

def row_sum(P0, q):
    """ row sums - 1 are zero """
    n = np.sqrt(P0.size)
    return P0.reshape(n,n).sum(1) - 1.

def non_neg(P0, q):
    """ all elements >= 0 """
    return P0

P_opt = fmin_slsqp(func, P0.ravel(), args=(q,), f_eqcons=row_sum, f_ieqcons=non_neg).reshape(d, d)

assert np.allclose(P_opt.sum(1), 1)
assert np.all(P_opt >= 0)
Autorizzato sotto: CC-BY-SA insieme a attribuzione
Non affiliato a StackOverflow
scroll top