Portfolio Optimization (Black Litterman Approach)

This example shows how to generate a MEX function and C source code from MATLAB code that performs portfolio optimization using the Black Litterman approach.

Prerequisites

There are no prerequisites for this example.

About the hlblacklitterman Function

The hlblacklitterman.m function reads in financial information regarding a portfolio and performs portfolio optimization using the Black Litterman approach.

type hlblacklitterman
function [er, ps, w, pw, lambda, theta] = hlblacklitterman(delta, weq, sigma, tau, P, Q, Omega)%#codegen
% hlblacklitterman
%   This function performs the Black-Litterman blending of the prior
%   and the views into a new posterior estimate of the returns as
%   described in the paper by He and Litterman.
% Inputs
%   delta  - Risk tolerance from the equilibrium portfolio
%   weq    - Weights of the assets in the equilibrium portfolio
%   sigma  - Prior covariance matrix
%   tau    - Coefficiet of uncertainty in the prior estimate of the mean (pi)
%   P      - Pick matrix for the view(s)
%   Q      - Vector of view returns
%   Omega  - Matrix of variance of the views (diagonal)
% Outputs
%   Er     - Posterior estimate of the mean returns
%   w      - Unconstrained weights computed given the Posterior estimates
%            of the mean and covariance of returns.
%   lambda - A measure of the impact of each view on the posterior estimates.
%   theta  - A measure of the share of the prior and sample information in the
%            posterior precision.

% Reverse optimize and back out the equilibrium returns
% This is formula (12) page 6.
pi = weq * sigma * delta;
% We use tau * sigma many places so just compute it once
ts = tau * sigma;
% Compute posterior estimate of the mean
% This is a simplified version of formula (8) on page 4.
er = pi' + ts * P' * inv(P * ts * P' + Omega) * (Q - P * pi');
% We can also do it the long way to illustrate that d1 + d2 = I
d = inv(inv(ts) + P' * inv(Omega) * P);
d1 = d * inv(ts);
d2 = d * P' * inv(Omega) * P;
er2 = d1 * pi' + d2 * pinv(P) * Q;
% Compute posterior estimate of the uncertainty in the mean
% This is a simplified and combined version of formulas (9) and (15)
ps = ts - ts * P' * inv(P * ts * P' + Omega) * P * ts;
posteriorSigma = sigma + ps;
% Compute the share of the posterior precision from prior and views,
% then for each individual view so we can compare it with lambda
theta=zeros(1,2+size(P,1));
theta(1,1) = (trace(inv(ts) * ps) / size(ts,1));
theta(1,2) = (trace(P'*inv(Omega)*P* ps) / size(ts,1));
for i=1:size(P,1)
    theta(1,2+i) = (trace(P(i,:)'*inv(Omega(i,i))*P(i,:)* ps) / size(ts,1));
end
% Compute posterior weights based solely on changed covariance
w = (er' * inv(delta * posteriorSigma))';
% Compute posterior weights based on uncertainty in mean and covariance
pw = (pi * inv(delta * posteriorSigma))';
% Compute lambda value
% We solve for lambda from formula (17) page 7, rather than formula (18)
% just because it is less to type, and we've already computed w*.
lambda = pinv(P)' * (w'*(1+tau) - weq)';
end

% Black-Litterman example code for MatLab (hlblacklitterman.m)
% Copyright (c) Jay Walters, blacklitterman.org, 2008.
%
% Redistribution and use in source and binary forms, 
% with or without modification, are permitted provided 
% that the following conditions are met:
%
% Redistributions of source code must retain the above 
% copyright notice, this list of conditions and the following 
% disclaimer.
% 
% Redistributions in binary form must reproduce the above 
% copyright notice, this list of conditions and the following 
% disclaimer in the documentation and/or other materials 
% provided with the distribution.
%  
% Neither the name of blacklitterman.org nor the names of its
% contributors may be used to endorse or promote products 
% derived from this software without specific prior written
% permission.
%  
% THIS SOFTWARE IS PROVIDED BY THE COPYRIGHT HOLDERS AND 
% CONTRIBUTORS "AS IS" AND ANY EXPRESS OR IMPLIED WARRANTIES, 
% INCLUDING, BUT NOT LIMITED TO, THE IMPLIED WARRANTIES OF 
% MERCHANTABILITY AND FITNESS FOR A PARTICULAR PURPOSE ARE 
% DISCLAIMED. IN NO EVENT SHALL THE COPYRIGHT OWNER OR 
% CONTRIBUTORS BE LIABLE FOR ANY DIRECT, INDIRECT, INCIDENTAL, 
% SPECIAL, EXEMPLARY, OR CONSEQUENTIAL DAMAGES (INCLUDING, 
% BUT NOT LIMITED TO, PROCUREMENT OF SUBSTITUTE GOODS OR 
% SERVICES; LOSS OF USE, DATA, OR PROFITS; OR BUSINESS 
% INTERRUPTION) HOWEVER CAUSED AND ON ANY THEORY OF LIABILITY, 
% WHETHER IN CONTRACT, STRICT LIABILITY, OR TORT (INCLUDING 
% NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY OUT OF THE USE 
% OF THIS SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF SUCH 
% DAMAGE.
%
% This program uses the examples from the paper "The Intuition 
% Behind Black-Litterman Model  Portfolios", by He and Litterman,
% 1999.  You can find a copy of this  paper at the following url.
%     http:%papers.ssrn.com/sol3/papers.cfm?abstract_id=334304
%
% For more details on the Black-Litterman model you can also view
% "The BlackLitterman Model: A Detailed Exploration", by this author
% at the following url.
%     http:%www.blacklitterman.org/Black-Litterman.pdf
%

The %#codegen directive indicates that the MATLAB code is intended for code generation.

Generate the MEX Function for Testing

Generate a MEX function using the codegen command.

codegen hlblacklitterman -args {0, zeros(1, 7), zeros(7,7), 0, zeros(1, 7), 0, 0}

Before generating C code, you should first test the MEX function in MATLAB to ensure that it is functionally equivalent to the original MATLAB code and that no run-time errors occur. By default, codegen generates a MEX function named hlblacklitterman_mex in the current folder. This allows you to test the MATLAB code and MEX function and compare the results.

Run the MEX Function

Call the generated MEX function

testMex();
View 1
Country        P       mu      w*
Australia	     0	 4.328	 1.524
Canada   	     0	 7.576	 2.095
France   	 -29.5	 9.288	-3.948
Germany  	   100	 11.04	 35.41
Japan    	     0	 4.506	 11.05
UK       	 -70.5	 6.953	-9.462
USA      	     0	 8.069	 58.57
q        	     5
omega/tau	     0.0213
lambda   	     0.317
theta   	     0.0714
pr theta  	     0.929


View 1
Country        P       mu      w*
Australia	     0	 4.328	 1.524
Canada   	     0	 7.576	 2.095
France   	 -29.5	 9.288	-3.948
Germany  	   100	 11.04	 35.41
Japan    	     0	 4.506	 11.05
UK       	 -70.5	 6.953	-9.462
USA      	     0	 8.069	 58.57
q        	     5
omega/tau	     0.0213
lambda   	     0.317
theta   	     0.0714
pr theta  	     0.929

Execution Time - MATLAB function: 0.12721 seconds
Execution Time - MEX function   : 0.017795 seconds

Generate C Code

cfg = coder.config('lib');
codegen -config cfg hlblacklitterman  -args {0, zeros(1, 7), zeros(7,7), 0, zeros(1, 7), 0, 0}

Using codegen with the specified -config cfg option produces a standalone C library.

Inspect the Generated Code

By default, the code generated for the library is in the folder codegen/lib/hbblacklitterman/.

The files are:

dir codegen/lib/hlblacklitterman/
.                              hlblacklitterman_terminate.o   
..                             hlblacklitterman_types.h       
buildInfo.mat                  interface                      
codeInfo.mat                   inv.c                          
codedescriptor.dmr             inv.h                          
compileInfo.mat                inv.o                          
examples                       pinv.c                         
hlblacklitterman.a             pinv.h                         
hlblacklitterman.c             pinv.o                         
hlblacklitterman.h             rtGetInf.c                     
hlblacklitterman.o             rtGetInf.h                     
hlblacklitterman_data.c        rtGetInf.o                     
hlblacklitterman_data.h        rtGetNaN.c                     
hlblacklitterman_data.o        rtGetNaN.h                     
hlblacklitterman_initialize.c  rtGetNaN.o                     
hlblacklitterman_initialize.h  rt_nonfinite.c                 
hlblacklitterman_initialize.o  rt_nonfinite.h                 
hlblacklitterman_ref.rsp       rt_nonfinite.o                 
hlblacklitterman_rtw.mk        rtw_proj.tmw                   
hlblacklitterman_terminate.c   rtwtypes.h                     
hlblacklitterman_terminate.h   

Inspect the C Code for the hlblacklitterman.c Function

type codegen/lib/hlblacklitterman/hlblacklitterman.c
/*
 * File: hlblacklitterman.c
 *
 * MATLAB Coder version            : 5.0
 * C/C++ source code generated on  : 28-Feb-2020 22:57:44
 */

/* Include Files */
#include "hlblacklitterman.h"
#include "hlblacklitterman_data.h"
#include "hlblacklitterman_initialize.h"
#include "inv.h"
#include "pinv.h"
#include "rt_nonfinite.h"

/* Function Definitions */

/*
 * hlblacklitterman
 *    This function performs the Black-Litterman blending of the prior
 *    and the views into a new posterior estimate of the returns as
 *    described in the paper by He and Litterman.
 *  Inputs
 *    delta  - Risk tolerance from the equilibrium portfolio
 *    weq    - Weights of the assets in the equilibrium portfolio
 *    sigma  - Prior covariance matrix
 *    tau    - Coefficiet of uncertainty in the prior estimate of the mean (pi)
 *    P      - Pick matrix for the view(s)
 *    Q      - Vector of view returns
 *    Omega  - Matrix of variance of the views (diagonal)
 *  Outputs
 *    Er     - Posterior estimate of the mean returns
 *    w      - Unconstrained weights computed given the Posterior estimates
 *             of the mean and covariance of returns.
 *    lambda - A measure of the impact of each view on the posterior estimates.
 *    theta  - A measure of the share of the prior and sample information in the
 *             posterior precision.
 * Arguments    : double delta
 *                const double weq[7]
 *                const double sigma[49]
 *                double tau
 *                const double P[7]
 *                double Q
 *                double Omega
 *                double er[7]
 *                double ps[49]
 *                double w[7]
 *                double pw[7]
 *                double *lambda
 *                double theta[3]
 * Return Type  : void
 */
void hlblacklitterman(double delta, const double weq[7], const double sigma[49],
                      double tau, const double P[7], double Q, double Omega,
                      double er[7], double ps[49], double w[7], double pw[7],
                      double *lambda, double theta[3])
{
  int i;
  double b;
  int i1;
  double ts[49];
  double pi[7];
  double y_tmp;
  double b_P;
  double b_b;
  double er_tmp[7];
  double b_y_tmp[7];
  double unusedExpr[7];
  double b_er_tmp[49];
  double posteriorSigma[49];
  int ps_tmp;
  double dv[49];
  if (!isInitialized_hlblacklitterman) {
    hlblacklitterman_initialize();
  }

  /*  Reverse optimize and back out the equilibrium returns */
  /*  This is formula (12) page 6. */
  for (i = 0; i < 7; i++) {
    b = 0.0;
    for (i1 = 0; i1 < 7; i1++) {
      b += weq[i1] * sigma[i1 + 7 * i];
    }

    pi[i] = b * delta;
  }

  /*  We use tau * sigma many places so just compute it once */
  for (i = 0; i < 49; i++) {
    ts[i] = tau * sigma[i];
  }

  /*  Compute posterior estimate of the mean */
  /*  This is a simplified version of formula (8) on page 4. */
  y_tmp = 0.0;
  b_P = 0.0;
  for (i = 0; i < 7; i++) {
    b = 0.0;
    b_b = 0.0;
    for (i1 = 0; i1 < 7; i1++) {
      b += ts[i + 7 * i1] * P[i1];
      b_b += P[i1] * ts[i1 + 7 * i];
    }

    b_y_tmp[i] = b_b;
    er_tmp[i] = b;
    y_tmp += b_b * P[i];
    b_P += P[i] * pi[i];
  }

  b_b = 1.0 / (y_tmp + Omega);
  b = Q - b_P;
  for (i = 0; i < 7; i++) {
    er[i] = pi[i] + er_tmp[i] * b_b * b;
  }

  /*  We can also do it the long way to illustrate that d1 + d2 = I */
  y_tmp = 1.0 / Omega;
  pinv(P, unusedExpr);

  /*  Compute posterior estimate of the uncertainty in the mean */
  /*  This is a simplified and combined version of formulas (9) and (15) */
  b = 0.0;
  for (i = 0; i < 7; i++) {
    b += b_y_tmp[i] * P[i];
  }

  b_b = 1.0 / (b + Omega);
  for (i = 0; i < 7; i++) {
    for (i1 = 0; i1 < 7; i1++) {
      b_er_tmp[i1 + 7 * i] = er_tmp[i1] * b_b * P[i];
    }
  }

  for (i = 0; i < 7; i++) {
    for (i1 = 0; i1 < 7; i1++) {
      b = 0.0;
      for (ps_tmp = 0; ps_tmp < 7; ps_tmp++) {
        b += b_er_tmp[i + 7 * ps_tmp] * ts[ps_tmp + 7 * i1];
      }

      ps_tmp = i + 7 * i1;
      ps[ps_tmp] = ts[ps_tmp] - b;
    }
  }

  for (i = 0; i < 49; i++) {
    posteriorSigma[i] = sigma[i] + ps[i];
  }

  /*  Compute the share of the posterior precision from prior and views, */
  /*  then for each individual view so we can compare it with lambda */
  inv(ts, dv);
  for (i = 0; i < 7; i++) {
    for (i1 = 0; i1 < 7; i1++) {
      b = 0.0;
      for (ps_tmp = 0; ps_tmp < 7; ps_tmp++) {
        b += dv[i + 7 * ps_tmp] * ps[ps_tmp + 7 * i1];
      }

      ts[i + 7 * i1] = b;
    }
  }

  b = 0.0;
  for (ps_tmp = 0; ps_tmp < 7; ps_tmp++) {
    b += ts[ps_tmp + 7 * ps_tmp];
  }

  theta[0] = b / 7.0;
  for (i = 0; i < 7; i++) {
    for (i1 = 0; i1 < 7; i1++) {
      b_er_tmp[i1 + 7 * i] = P[i1] * y_tmp * P[i];
    }
  }

  for (i = 0; i < 7; i++) {
    for (i1 = 0; i1 < 7; i1++) {
      b = 0.0;
      for (ps_tmp = 0; ps_tmp < 7; ps_tmp++) {
        b += b_er_tmp[i + 7 * ps_tmp] * ps[ps_tmp + 7 * i1];
      }

      ts[i + 7 * i1] = b;
    }
  }

  b = 0.0;
  for (ps_tmp = 0; ps_tmp < 7; ps_tmp++) {
    b += ts[ps_tmp + 7 * ps_tmp];
  }

  theta[1] = b / 7.0;
  for (i = 0; i < 7; i++) {
    for (i1 = 0; i1 < 7; i1++) {
      b_er_tmp[i1 + 7 * i] = P[i1] * y_tmp * P[i];
    }
  }

  for (i = 0; i < 7; i++) {
    for (i1 = 0; i1 < 7; i1++) {
      b = 0.0;
      for (ps_tmp = 0; ps_tmp < 7; ps_tmp++) {
        b += b_er_tmp[i + 7 * ps_tmp] * ps[ps_tmp + 7 * i1];
      }

      ts[i + 7 * i1] = b;
    }
  }

  b = 0.0;
  for (ps_tmp = 0; ps_tmp < 7; ps_tmp++) {
    b += ts[ps_tmp + 7 * ps_tmp];
  }

  theta[2] = b / 7.0;

  /*  Compute posterior weights based solely on changed covariance */
  for (i = 0; i < 49; i++) {
    b_er_tmp[i] = delta * posteriorSigma[i];
  }

  inv(b_er_tmp, dv);
  for (i = 0; i < 7; i++) {
    b = 0.0;
    for (i1 = 0; i1 < 7; i1++) {
      b += er[i1] * dv[i1 + 7 * i];
    }

    w[i] = b;
  }

  /*  Compute posterior weights based on uncertainty in mean and covariance */
  for (i = 0; i < 49; i++) {
    posteriorSigma[i] *= delta;
  }

  inv(posteriorSigma, dv);
  for (i = 0; i < 7; i++) {
    b = 0.0;
    for (i1 = 0; i1 < 7; i1++) {
      b += pi[i1] * dv[i1 + 7 * i];
    }

    pw[i] = b;
  }

  /*  Compute lambda value */
  /*  We solve for lambda from formula (17) page 7, rather than formula (18) */
  /*  just because it is less to type, and we've already computed w*. */
  pinv(P, er_tmp);
  *lambda = 0.0;
  for (ps_tmp = 0; ps_tmp < 7; ps_tmp++) {
    *lambda += er_tmp[ps_tmp] * (w[ps_tmp] * (tau + 1.0) - weq[ps_tmp]);
  }

  /*  Black-Litterman example code for MatLab (hlblacklitterman.m) */
  /*  Copyright (c) Jay Walters, blacklitterman.org, 2008. */
  /*  */
  /*  Redistribution and use in source and binary forms,  */
  /*  with or without modification, are permitted provided  */
  /*  that the following conditions are met: */
  /*  */
  /*  Redistributions of source code must retain the above  */
  /*  copyright notice, this list of conditions and the following  */
  /*  disclaimer. */
  /*   */
  /*  Redistributions in binary form must reproduce the above  */
  /*  copyright notice, this list of conditions and the following  */
  /*  disclaimer in the documentation and/or other materials  */
  /*  provided with the distribution. */
  /*    */
  /*  Neither the name of blacklitterman.org nor the names of its */
  /*  contributors may be used to endorse or promote products  */
  /*  derived from this software without specific prior written */
  /*  permission. */
  /*    */
  /*  THIS SOFTWARE IS PROVIDED BY THE COPYRIGHT HOLDERS AND  */
  /*  CONTRIBUTORS "AS IS" AND ANY EXPRESS OR IMPLIED WARRANTIES,  */
  /*  INCLUDING, BUT NOT LIMITED TO, THE IMPLIED WARRANTIES OF  */
  /*  MERCHANTABILITY AND FITNESS FOR A PARTICULAR PURPOSE ARE  */
  /*  DISCLAIMED. IN NO EVENT SHALL THE COPYRIGHT OWNER OR  */
  /*  CONTRIBUTORS BE LIABLE FOR ANY DIRECT, INDIRECT, INCIDENTAL,  */
  /*  SPECIAL, EXEMPLARY, OR CONSEQUENTIAL DAMAGES (INCLUDING,  */
  /*  BUT NOT LIMITED TO, PROCUREMENT OF SUBSTITUTE GOODS OR  */
  /*  SERVICES; LOSS OF USE, DATA, OR PROFITS; OR BUSINESS  */
  /*  INTERRUPTION) HOWEVER CAUSED AND ON ANY THEORY OF LIABILITY,  */
  /*  WHETHER IN CONTRACT, STRICT LIABILITY, OR TORT (INCLUDING  */
  /*  NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY OUT OF THE USE  */
  /*  OF THIS SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF SUCH  */
  /*  DAMAGE. */
  /*  */
  /*  This program uses the examples from the paper "The Intuition  */
  /*  Behind Black-Litterman Model  Portfolios", by He and Litterman, */
  /*  1999.  You can find a copy of this  paper at the following url. */
  /*      http:%papers.ssrn.com/sol3/papers.cfm?abstract_id=334304 */
  /*  */
  /*  For more details on the Black-Litterman model you can also view */
  /*  "The BlackLitterman Model: A Detailed Exploration", by this author */
  /*  at the following url. */
  /*      http:%www.blacklitterman.org/Black-Litterman.pdf */
  /*  */
}

/*
 * File trailer for hlblacklitterman.c
 *
 * [EOF]
 */