2

我正在尝试以八度解决以下优化问题

公式

第一个约束是 A 是半正定的。S 是一组数据点，如果 (xi,xj) 在 S 中，则 xi 与 xj 相似，D 是一组数据点，如果 (xi,xj) 在 D 中，则 xi 和 xj 不相似。请注意，上面的公式是 2 个单独的总和，第二个总和没有嵌套。还假设 xi 和 xj 是长度为 N 的列向量。

因为这是一个非线性优化，我正在尝试使用 octave 的非线性程序求解器 sqp。 问题是，如果我只是为它提供优化的功能，在一些小玩具测试中，找到Hessian的 BFGS 方法会 失败。因此，我尝试提供自己的 Hessian 函数，但现在出现此问题

error: __qp__: operator *: nonconformant arguments (op1 is 2x2, op2 is 3x1)
error: called from:
error:   /usr/share/octave/3.6.3/m/optimization/qp.m at line 393, column 26
error:   /usr/share/octave/3.6.3/m/optimization/sqp.m at line 414, column 32

当我对 sqp 进行以下调用时

[A, ~, Info] = sqp(initial_guess, {@toOpt, @CalculateGradient,@CalculateHessian},
[],[],0,[],maxiter);

我通过仅求解对角线条目并将所有对角线条目限制为 >=0 来简化 A 为半正定和对角线的约束。initial_guess 是 N 长的向量。

这是我计算我认为是 Hessian 矩阵的代码

%Hessian = CalculateHessian(A)
%calculates the Hessian of the function we are optimizing as follows
%H(i,j) = (sumsq(D(:,i),1) * sumsq(D(:,j),1)) / (sum(A.*sumsq(D,1))^2)
%where D is a matrix of of differences between observations that are dissimilar, with one difference on each row
%and sumsq is the sum of the squares
%input A: the current guess for A
%output Hessian: The hessian of the function we are optimizing
function Hessian = CalculateHessian(A)
    global HessianNumerator; %this is a matrix with the numerator of H(i,j)
    global Dsum_of_squares; %the sum of the squares of the differences of each dimensions of the dissimilar observations

    if(iscolumn(A)) %if A is a column vector
        A = A'; %make it a row vector. necessary to prevent broadcasting
    endif

    if(~isempty(Dsum_of_squares)) %if disimilar constraints were provided   
        Hessian = HessianNumerator / (sum(A.*Dsum_of_squares)^2)
    else
        Hessian = HessianNumerator; %the hessian is a matrix of 0s
    endif

endfunction

和 Dsum_of_squares 和 HessianNumertor 是

[dissimilarRow,dissimilarColumn] = find(D); %find which observations are dissimilar to each other
DissimilarDiffs = X(dissimilarRow,:) - X(dissimilarColumn,:); %take the difference between the dissimilar observations
Dsum_of_squares = sumsq(DissimilarDiffs,1);
HessianNumerator = Dsum_of_squares .* Dsum_of_squares'; %calculate the numerator of the Hessian. it is a constant value

X 是一个 M x N 矩阵，每行有一个观测值。

D 是 M x M 相异矩阵。如果 D(i,j) 为 1，则 X 的第 i 行与第 j 行不同。否则为 0。

我相信我的错误出现在以下领域之一（从最不可能到最有可能）

  1. 我用来推导 Hessian 函数的数学是错误的。我使用的公式在我对该函数的评论中。
  2. 我的数学实现。
  3. sqp 想要的 Hessian 矩阵与 Hessian Matrix Wikipedia 页面上描述的不同。

任何帮助将不胜感激。如果您需要我发布更多代码，我很乐意这样做。目前，尝试解决优化问题的代码量约为 160 行。

这是我正在运行的导致代码失败的测试用例。如果我只将梯度函数传递给它，它就可以工作。

X = [1 2 3; 
     4 5 6; 
     7 8 9; 
     10 11 12];
S = [0 1 1 0; 
     1 0 0 0; 
     1 0 0 0; 
     0 0 0 0]; %this means row 1 of X is similar to rows 2 and 3
D = [0 0 0 0; 
     0 0 0 0;
     0 0 0 1;
     0 0 1 0]; %this means row 3 of X is dissimilar to row 4
gml(X,S,D, 200); %200 is the maximum number of iterations for sqp to run
4

0 回答 0