Tikhonov regularized least-squares with unit simplex constraint

1k Views Asked by At

First, I got this $\ell_2$ norm problem, which is easy to solve

$$w_0 = \arg\min_⁡ \| −_ \|_2^2+ \lambda \|\|_2^2 \tag{1}$$

where $X$ is a matrix and $\lambda$ is a fixed scalar. Then, I added the following constraints

$$∈[0,1], \qquad \thinspace \sum_ _i =1 \tag{2}$$

At first I found this optimization problem similar to fmincon function of matlab (reference here), but I keep getting poor results (significant higher error than results in Eq(1)).

I already try with other matlab build-in functions such as lsqlin, lsqnonneg but all return bad results (err = $\|X*w_o - _ \|_2^2 $) such as

Eq2-fmincon: 75596.826, Eq2-lsqlin: 386.2777, Eq2-lsnonneg:311.5511, Eq1: 0.0019292

Is there any suggestion on how to solve this problem? Thank you for your time.


My matlab code for fmincon are follow

gamma = 1; 
fun = @(w) ((norm(X * w  - xc, 2))^2 + gamma * (norm(w, 2))^2);                    
w0  = rand(size(X, 2), 1); %w0 = w0./sum(w0); 
A   = [];                   b   = []; 
Aeq = ones(1, size(X, 2));  beq = 1; 
lb  = zeros(size(X, 2), 1); ub  = ones(size(X, 2), 1); 
w1  = fmincon(fun,w0,A,b,Aeq,beq,lb,ub) ;
w2  = lsqlin(X, xc, A, b, Aeq, beq, lb, ub); 
w3  = lsqnonneg(X, xc);
w4  = pinv(X' * X + gamma * eye(size(X,2))) *X' * xc;

tmp1 = norm(X * w1 - xc, 2);     tmp2 = norm(X * w2 - xc, 2); 
tmp3 = norm(X * w3 - xc, 2);     tmp4 = norm(X * w4 - xc, 2); 

display(['Eq2-fmincon: ' num2str(tmp1) ', Eq2-lsqlin: ' num2str(tmp2) ', Eq2-lsnonneg:' ...
            num2str(tmp3) ', Eq1: ' num2str(tmp4)]);
2

There are 2 best solutions below

1
On BEST ANSWER

This is a regularized least-squares (RLS) problem subject to the standard $(n-1)$-simplex. We have the following quadratic program (QP)

$$\begin{array}{ll} \text{minimize} & \| \mathrm A \mathrm x - \mathrm b \|_2^2 + \lambda \| \mathrm x \|_2^2\\ \text{subject to} & 1_n^{\top} \mathrm x = 1\\ & \mathrm x \geq \mathrm 0_n\end{array}$$

which can be rewritten as follows

$$\begin{array}{ll} \text{minimize} & \mathrm x^{\top} \left( \mathrm A^{\top} \mathrm A + \lambda \mathrm I_n \right) \mathrm x - 2 \mathrm b^{\top} \mathrm A \mathrm x + \mathrm b^{\top} \mathrm b\\ \text{subject to} & 1_n^{\top} \mathrm x = 1\\ & \mathrm x \geq \mathrm 0_n\end{array}$$

Do note that if $\lambda < 0$, then the QP may be non-convex.

In MATLAB, one can use function quadprog to solve this QP.

0
On

You can also solve this by Projected Gradient Descent since the projection onto the Unit Simplex is known.

The Code:

vX = pinv(mA) * vB;

mAA = mA.' * mA;
mAb = mA.' * vB;

for ii = 1:numIterations
    stepSize = stepSizeBase / sqrt(ii);
    vG = (mAA * vX) - mAb + (2 * paramLambda * vX);
    vX = vX - (stepSize * vG);
    vX = ProjectSimplex(vX, simplexRadius, stopThr);
end

objVal = (0.5 * sum((mA * vX - vB) .^ 2)) + (paramLambda * sum(vX .^ 2));

disp([' ']);
disp(['Projected Gradient Solution Summary']);
disp(['The Optimal Value Is Given By - ', num2str(objVal)]);
disp(['The Optimal Argument Is Given By - [ ', num2str(vX.'), ' ]']);
disp([' ']);

The full code (Validated by CVX) is in my StackExchange Mathematics Q1983030 GitHub Repository.

You can even accelerate this using Nesterov / FISTA to get really fast and efficient algorithm.