How can I add constraints in a solution of a function?

Hi everyone,

I would like to solve the following function by adding two constraints.

S=randn(6,6); y=randn(6,1); ONE = ones(6,1); rft=randn(1,1); K = randn(1,1);     
    x = inv(S)*(y- ONE*rft)*0.1/sqrt(K);

I would like to include 2 additional constraints, -1<sum(x)<2 . I am not sure how I should use fmincon.

[Walter Roberson]

What variable are you minimizing with respect to?

[gsourop]

Apologies, I missed that part. So, I want to maximize ReTurn

ReTurn = (1-x'*ONE)*rft+x'*y;

[Bruno Luong]

What S has anything to do then in your command

x = inv(S)*(y- ONE*rft)*0.1/sqrt(K)

[gsourop]

The first constraint of

ReTurn = (1-x'*ONE)*rft+x'*y;

is

 (0.1)^2=x'Sx,

the algebric solution of which is

x=(0.1/sqrt(K)) * inv(S) * (y -  ONE*rft)

[Bruno Luong]

You still haven't told us everything relevant: the cost function (return) and the constraint do not have K anywhere, yet you claim the algebraic solution has some expression which is ~ 1/sqrt(K).

Why?

[gsourop]

I missed it because I thought we will go far away from the practical issue and dive into theory. However, the theoretical problem is the following: an investor is investing into 6 foreign risk free assets and a domestic risk free asset. Hence, the expected return for the next period should be

y = ONE * rf_domestic + forecast_currency_y  % y is a 6x1 vector, rf_domestic is a scalar.

assume that S the conditional variance-covariance matrix of y. At the end of the period the investor allocates the weights on each asset in order to maximize the returns of the portfolio

ReTurn = (1-x'*ONE)*rft+x'*y;

s.t.

(0.1)^2=x'Sx

the algebric solution, namely the weight x on each foreign asset, is given by:

x=(0.1/sqrt(K)) * inv(S) * (y -  ONE*rft)  % x is a 6x1 vector

where

K=(y-ONE*rft)' * inv(S) * (y-ONE*rft); %y-ONE*rft is a 6x1 vector

Now I want to impose 2 additional constraits in the existing one. The new constraint should not allow the sum of weights to exceed

-1<sum(x)<2

[Bruno Luong]

And S is symmetric definite positive matrix?

I think it is much better if you give an complete and accurate description of your problem instead of hiding relevant information that you think it's not.

[gsourop]

OK, I understand. It must be. Otherwise, I would get an orange message about the singularity or bad scaling of the matrix, right?

[Bruno Luong]

Yes, even the sqrt() might return complex result.

[gsourop]

I see. But is there a way to impose the additional constraints without solving the relationship arithmetically?

Answer 1

[Bruno Luong]

n = 6;
L = randn(n);
S = L'*L;
y = randn(n,1);
rft = randn(1,1);
su = 2;
sl = -1;

% T = V*V'
T = S/(0.1^2);
T = 0.5*(T+T');
[V,D] = eig(T);
V = V.*sqrt(diag(D)');
A = inv(V');
% xx = V'*x
% x = V'\xx = W'*xx
% |xx| = 1
yy = V \ (-rft+y);
SX = sum(A,1);

SXu = SX/su;
SXl = SX/sl;

% ReTurn = (1-x'*ONE)*rft+x'*y
% return = rtf - x'*(rtf + y)
% = rft + xx'*yy
% with yy = W*(-rft+y)

% maximize (xx'*yy)
% such that 
%   |xx| = 1
%   SXu*xx <= 2
%   SXl*xx <= 1

xx = yy/norm(yy);
if SXu*xx > 1
    n2 = (SXu*SXu');
    cs2 = 1/n2;
    sn = sqrt(1-cs2);
    Tu = null(SXu);
    yyu = Tu*(Tu'*yy);
    xx = cs2*SXu' + (sn/norm(yyu))*yyu;
elseif SXl*xx > 1
    n2 = (SXl*SXl');
    cs2 = 1/n2;
    sn = sqrt(1-cs2);
    Tl = null(SXl);
    yyl = Tl*(Tl'*yy);
    xx = cs2*SXl' + (sn/norm(yyl))*yyl;
end

x = V' \ xx

% Check constraints
sum(x)
x'*S*x

[Bruno Luong]

Impose the additional constraints but still able to solve the relationship "arithmetically".

[gsourop]

Thank you for answering the question. I think there is a bug in this line,

V = V.*sqrt(diag(D)');

Are you trying to multiply each row of V with the row vector sqrt(diag(D)') ?

[Bruno Luong]

No I want to scale each column of V by sqrt(diag(D)'), so the code does what I intend.

If you think there is a bug you have to explain why.

But you can check, it gives the same result than your formula when the sum does not hit the limits -1 or 2

[gsourop]

I am getting an error on the dimensions of the matrix. Hence, I replaces this line with

Scaling = sqrt(diag(D)');
for i = 1:6
    for j =1 : 6
        V(i,j) = V(i,j).* Scaling(j);
    end
end

but the validations at the end do not give me the expected result.

[Bruno Luong]

"I am getting an error on the dimensions of the matrix"

Then apply my code wrongly. I provide the code working with some fake data

S: sym-def-pos matrix (6 x 6)
y: vector (6 x 1)