2014-09-23 28 views
1

我是scipy.optimize模塊的新手,需要幫助嘗試使用公式V的最小化函數與矩陣一起工作,並且有2個約束,但我不確定是否我我正在處理功能的形成或其中一個約束。我如何使用Scipy最大限度地減少約束和動態功能

M是N×N的矩陣,但比如我會給4x4矩陣

[[ 0.0.00011974 0.00067403 0.00060845] 
[ 0.00011974 0.00736665 0.00165674 0.00053581] 
[ 0.00067403 0.00165674 0.00281547 0.00129646] 
[ 0.00060845 0.00053581 0.00129646 0.00153195]] 

X是1×N個矩陣,其中所有數字都是正的,必須加起來爲1,就是我解決了

約束:

X[0] + X[1] + X[2] + … X[n] = 1 

E = X[0]*R[0] + X[1]*R[1] + X[2]*R[2] + … X[n]*R[n] (E is user input) 

我試圖最小化的功能是:

V = X[0]**2 * M[0][0] + 2*X[0]*X[1]*M[0][1] + 2*X[0]*X[2]*M[0][2] + 2 * X[0] * X[3] * M[0][3] * X[1]**2 * M[1][1] + ….etc 

我到目前爲止的代碼是

cons = ({'type': 'eq', 'fun': lambda x: sum(w) - 1}) 

    args = n, m, w 

    answer = scipy.optimize.minimize(func(*args), w, args=(), method='SLSQP', 
            bounds=((0, None), (0, None)), 
            constraints=cons) 
def func(n, m, x): 
    row = 0 
    col = 0 
    formula = [] 
    while row < n: 
     while col < n: 
      if col == row: 
       item = (x[row]**2) * m[row][col] 
       formula.append(item) 
       col += 1 
      else: 
       item2 = 2 * x[row] * x[col] * m[row][col] 
       formula.append(item2) 
       col += 1 

     row += 1 
     col = 0 

    return sum(formula) 

什麼我無法理解的是如何表達的第二個約束。我也不確定我是否正確處理了公式的創建。任何幫助非常感謝

+0

你的任務看起來像一個二次規劃問題(https://en.wikipedia.org/wiki/Quadratic_programming)。 CVXopt在這裏可能很有用(http://cvxopt.org/userguide/coneprog.html#quadratic-programming) – Dietrich 2014-09-23 23:15:15

回答

3

正如在評論中提到的,這可以通過專門的凸解算器來解決。


from __future__ import print_function 

import numpy as np 
import cvxpy 
from scipy.optimize import minimize 

# Problem data. 
n = 4 
M = [ 
     [0.0, 0.00011974, 0.00067403, 0.00060845], 
     [0.00011974, 0.00736665, 0.00165674, 0.00053581], 
     [0.00067403, 0.00165674, 0.00281547, 0.00129646], 
     [0.00060845, 0.00053581, 0.00129646, 0.00153195]] 
np.random.seed(1) 
R = np.random.randn(n) 
E = np.random.randn() 

print('R:', R) 
print('E:', E) 
print() 

A = np.matrix([np.ones(n), R]) 
b = np.matrix([[1], [E]]) 


def solve_cvxpy(): 

    # Construct the problem. 
    x = cvxpy.Variable(n) 
    objective = cvxpy.Minimize(cvxpy.quad_form(x, M)) 
    constraints = [0 <= x, A*x == b] 
    problem = cvxpy.Problem(objective, constraints) 

    # Solve the problem. 
    result = problem.solve() 
    print('cvxpy solution:') 
    print(x.value) 
    print() 


def solve_scipy(): 
    def simplex_constraint(x): 
     return x.sum() - 1 
    def dot_constraint(x): 
     return x.dot(R) - E 
    def objective(x): 
     return x.dot(M).dot(x) 
    x0 = np.ones(n)/n 
    bounds = [(0, np.inf)]*n 
    constraints = (
      dict(type='eq', fun=simplex_constraint), 
      dict(type='eq', fun=dot_constraint)) 
    result = minimize(
      objective, x0, method='SLSQP', tol=1e-8, 
      bounds=bounds, constraints=constraints) 
    print('scipy solution:') 
    print(result.x) 
    print() 

solve_cvxpy() 
solve_scipy() 

R: [ 1.62434536 -0.61175641 -0.52817175 -1.07296862] 
E: 0.865407629325 

cvxpy solution: 
[[ 7.03877116e-01] 
[ 8.62848827e-02] 
[ 5.54383296e-06] 
[ 2.09832457e-01]] 

scipy solution: 
[ 7.03877583e-01 8.62886986e-02 -1.24818775e-17 2.09833718e-01] 
相關問題