Pergunta

Eu sou novo na otimização. Estou tentando resolver um problema de mínimos quadrados lineares usando o fmin_slsqp função em scipy.optimize.

Eu tenho a função objetiva como norma de Frobenius de |q0_T*P-q1_T| quadrado, onde q0_T é a transposição de um nX1 vetor e P é nXn matriz e q1_T é uma transposição de um nX1 vetor. Este é basicamente um processo de Markov com q vetores como a probabilidade de estar em um estado e P sendo a matriz das probabilidades de transição.

A função objetiva seria minimizada wrt P Onde estão as restrições:

1) Todos os elementos em P deve ser não negativo

2) todas as fileiras em P deve resumir para 1

Eu defini essa função objetiva que não tenho certeza se está correta:

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

O segundo argumento no FMIN_SLSQP solicita um 1D NDARRAY X0, que é o palpite inicial para a variável independente. Mas aqui, como minha variável independente seria P, eu precisaria de uma matriz 2D para o meu palpite inicial. Não tenho certeza se estou enquadrando o problema corretamente ou terei que usar outra função. Obrigado.

Foi útil?

Solução

Portanto, há um problema em que ele aceita apenas uma matriz 1D como a variável independente, uma solução hackish seria passá -la como 1D e remodelar na função.

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)
Licenciado em: CC-BY-SA com atribuição
Não afiliado a StackOverflow
scroll top