0

我试图找到一个在乘以矩阵时最小化残差平方和的向量。

我知道 scipy 的优化包(它具有最小化功能)。但是,我的代码有一个额外的限制。w 的所有条目的总和(见下面的函数)必须等于 1,并且 w 的任何条目都不能小于 0。有没有一个包可以为我做这个?如果没有,我该怎么做?

试图最小化 w:

def w_rss(w,x0,x1):
    predictions = np.dot(x0,w)
    errors = x1 - predictions
    rss = np.dot(errors.transpose(),errors).item(0)

    return rss

X0 = np.array([[3,4,5,3],
               [1,2,2,4],
               [6,5,3,7],
               [1,0,5,2]])  

X1 = np.array([[4],
               [2],
               [4],
               [2]]) 

W = np.array([[.0],
              [.5],
              [.5],
              [.0]])

print w_rss(W,X0,X1)

到目前为止,这是我循环遍历 w 的可能值的最佳尝试,但它无法正常工作。

def get_w(x0,x1):

J = x0.shape[1]
W0 = np.matrix([[1.0/J]*J]).transpose()
rss0 = w_rss(W0,x0,x1)
loop = range(J)
for i in loop:
    W1 = W0
    rss1 = rss0
    while rss0 == rss1:
        den = len(loop)-1
        W1[i][0] += 0.01
        for j in loop:
            if i == j:
                continue
            W1[j][0] -= 0.01/den
            if W1[j][0] <= 0:
                loop.remove(j)
        rss1 = w_rss(W1,x0,x1)
        if rss1 < rss0:
            #print W1
            W0 = W1
            rss0 = rss1
        print '--'
        print rss0
        print W0

return W0,rss0
4

1 回答 1

3

scipy 中的 SLSQP 代码可以做到这一点。可以使用scipy.optimize.minimizewith method='SLSQP,也可以fmin_slsqp直接使用函数。在下文中,我使用fmin_slsqp.

scipy 求解器通常将一维数组传递给目标函数,因此为了保持一致,我将更改为一维数组,W并且X1我将编写目标函数(现在称为w_rss1)以期望一维论据w

使用参数指定所有元素w必须介于 0 和 1 之间bounds的条件,使用参数指定总和必须为 1 的条件f_eqcons。约束函数返回np.sum(w) - 1,所以元素和为 1 时为 0。

这是代码:

import numpy as np
from scipy.optimize import fmin_slsqp


def w_rss1(w, x0, x1):
    predictions = np.dot(x0, w)
    errors = x1 - predictions
    rss = (errors**2).sum()
    return rss


def sum1constraint(w, x0, x1):
    return np.sum(w) - 1


X0 = np.array([[3,4,5,3],
               [1,2,2,4],
               [6,5,3,7],
               [1,0,5,2]])  

X1 = np.array([4, 2, 4, 2]) 

W = np.array([.0, .5, .5, .0])

result = fmin_slsqp(w_rss1, W, f_eqcons=sum1constraint, bounds=[(0.0, 1.0)]*len(W),
                    args=(X0, X1), disp=False, full_output=True)
Wopt, fW, its, imode, smode = result

if imode != 0:
    print("Optimization failed: " + smode)
else:
    print(Wopt)

当我运行它时,输出是

[ 0.05172414  0.55172414  0.39655172  0.        ]
于 2016-12-14T16:52:42.163 回答