Incorrect negative portfolio weights for constrained minimum-variance portfolios
Hi Community,
I implemented the Jagannathan & Ma (Journal of Finance, 2003) algorithm, which reduces the error in the sample covariance matrix to improve the out-of-sample performance of portfolios (at the end of this post is the excerpt from the paper with the theoretical derivation). The resulting portfolio weights are supposed to be nonnegative and in particular it should work for large portfolios where the number of assets exceed the number of observations (N>T).
However, the results from the Gurobi optimizer show negative weights. And in particular if I increase N to N=100, the results get even worse, i.e. more negative.
Could you please check the implementation as it appears to be an issue with the Lagrangian multipliers, especially if N>T?
Below is the code I am using as well as the set up of the Gurobi function.
Thank you very much in advance.
% This is a simple code to reproduce the issue
% Generate test data
N=30;
data = randn(60,N).*repmat(rand(1,N)/10,60,1);
% Construct sample covariance matrix
S = cov(data);
N = size(data, 2);
ConstrL = 0; % non-negativity constraint
% Run optimization using Gurobi function
[x,fval,exitflag,output,lambda]=gurobiquadprog(S,zeros(1,N),[],[],ones(1,N),1, repmat(ConstrL,1,N));
% Jagannathan and Ma correction (New JF-Article, no impact on treatment of shortsale constraints)
J = -(lambda.upper*ones(1,N)+ones(N,1)*lambda.upper')...
+(lambda.lower*ones(1,N)+ones(N,1)*lambda.lower');
% Shrink the covariance matrix
SS = S-J;
% Compute portfolio weights
w = (SS\ones(N,1))/(ones(N,1)'/SS*ones(N,1))
% This is the gurobi optimization function
function [x,fval,exitflag,output,lambda] = gurobiquadprog(H,f,A,b,Aeq,beq,lb,ub,x0,options)
% Initialize missing arguments
if nargin < 10
options = struct();
if nargin == 9
if isa(x0,'struct') || isa(x0,'optim.options.SolverOptions')
options = x0; % x0 was omitted and options were passed instead
x0 = [];
end
else
x0 = [];
if nargin < 8
ub = [];
if nargin < 7
lb = [];
if nargin < 6
beq = [];
if nargin < 5
Aeq = [];
end
end
end
end
end
end
% Warn user if x0 argument is ignored
if ~isempty(x0)
warning('LINPROG will ignore non-empty starting point X0');
end
% Build Gurobi model
model.obj = f;
model.Q = sparse(H);
model.A = [sparse(A); sparse(Aeq)]; % A must be sparse
model.sense = [repmat('<',size(A,1),1); repmat('=',size(Aeq,1),1)];
model.rhs = full([b(:); beq(:)]); % rhs must be dense
if ~isempty(lb)
model.lb = lb;
else
model.lb = -inf(size(model.A,2),1); % default lb for MATLAB is -inf
end
if ~isempty(ub)
model.ub = ub;
end
% Extract relevant Gurobi parameters from (subset of) options
params = struct();
if isfield(options,'Display') || isa(options,'optim.options.SolverOptions')
if any(strcmp(options.Display,{'off','none'}))
params.OutputFlag = 0;
end
end
if isfield(options,'MaxTime') || isa(options,'optim.options.SolverOptions')
params.TimeLimit = options.MaxTime;
end
params.outputflag = 0;
params.NonConvex = 2; % avoid issue for non-convex cases but leads to
% matlab freeze
params.PSDTol = inf;
% Solve model with Gurobi
result = gurobi(model,params);
% Resolve model if status is INF_OR_UNBD
if strcmp(result.status,'INF_OR_UNBD')
params.DualReductions = 0;
warning('Infeasible or unbounded, resolve without dual reductions to determine...');
result = gurobi(model,params);
end
% Collect results
x = [];
output.message = result.status;
output.constrviolation = [];
if isfield(result,'x')
x = result.x;
if nargout > 3
slack = model.A*x-model.rhs;
violA = slack(1:size(A,1));
violAeq = norm(slack((size(A,1)+1):end),inf);
viollb = model.lb(:)-x;
violub = 0;
if isfield(model,'ub')
violub = x-model.ub(:);
end
output.constrviolation = max([0; violA; violAeq; viollb; violub]);
end
end
fval = [];
if isfield(result,'objval')
fval = result.objval;
end
if strcmp(result.status,'OPTIMAL')
exitflag = 1; % converged to a solution
elseif strcmp(result.status,'UNBOUNDED')
exitflag = -3; % problem is unbounded
elseif strcmp(result.status,'ITERATION_LIMIT')
exitflag = 0; % maximum number of iterations reached
else
exitflag = -2; % no feasible point found
end
lambda.lower = [];
lambda.upper = [];
lambda.ineqlin = [];
lambda.eqlin = [];
if nargout > 4
if isfield(result,'rc')
lambda.lower = max(0,result.rc);
lambda.upper = -min(0,result.rc);
end
if isfield(result,'pi')
if ~isempty(A)
lambda.ineqlin = -result.pi(1:size(A,1));
end
if ~isempty(Aeq)
lambda.eqlin = -result.pi((size(A,1)+1):end);
end
end
end
% Local Functions =========================================================
function [f,A,b,Aeq,beq,lb,ub,x0,options] = probstruct2args(s)
%PROBSTRUCT2ARGS Get problem structure fields ([] is returned when missing)
f = getstructfield(s,'f');
A = getstructfield(s,'Aineq');
b = getstructfield(s,'bineq');
Aeq = getstructfield(s,'Aeq');
beq = getstructfield(s,'beq');
lb = getstructfield(s,'lb');
ub = getstructfield(s,'ub');
x0 = getstructfield(s,'x0');
options = getstructfield(s,'options');
function f = getstructfield(s,field)
%GETSTRUCTFIELD Get structure field ([] is returned when missing)
if isfield(s,field)
f = getfield(s,field);
else
f = [];
end
-
Official comment
This post is more than three years old. Some information may not be up to date. For current information, please check the Gurobi Documentation or Knowledge Base. If you need more help, please create a new post in the community forum. Or why not try our AI Gurobot?. -
Hi Christian,
Thank you for participating in our community forum! Could you please edit your post to make it more readable? Especially the code section should be put into a code block to preserve indentation. Please see our forum instructions for additional info.
Thanks,
Matthias0 -
HI Matthias,
Thanks a lot for your reply and tips. I have adjusted the post accordingly and also added a screenshot from the paper for more clarity.
Thanks,
Christian
0
Post is closed for comments.
Comments
3 comments