# 使用 Black Litterman 方法生成用于优化投资组合的代码

### 关于 hlblacklitterman 函数

hlblacklitterman.m 函数读入关于投资组合的财务信息，并使用 Black Litterman 方法执行投资组合优化。

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
%

%#codegen 指令指示 MATLAB 代码用于代码生成。

### 生成用于测试的 MEX 函数

codegen hlblacklitterman -args {0, zeros(1, 7), zeros(7,7), 0, zeros(1, 7), 0, 0}
Code generation successful.

### 运行 MEX 函数

testMex();
View 1
Country        P       mu      w*
Australia	     0	 4.328	 1.524
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
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.050145 seconds
Execution Time - MEX function   : 0.012961 seconds

### 生成 C 代码

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

codegen 与指定的 -config cfg 选项结合使用生成独立 C 库。

### 检查生成的代码

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

### 检查 hlblacklitterman.c 函数的 C 代码

type codegen/lib/hlblacklitterman/hlblacklitterman.c
/*
* File: hlblacklitterman.c
*
* MATLAB Coder version            : 23.2
* C/C++ source code generated on  : 03-Aug-2023 18:53:10
*/

/* Include Files */
#include "hlblacklitterman.h"
#include "inv.h"
#include "rt_nonfinite.h"
#include "rt_nonfinite.h"
#include <math.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])
{
double b_er_tmp[49];
double dv[49];
double posteriorSigma[49];
double ts[49];
double er_tmp[7];
double pi[7];
double y_tmp[7];
double absxk;
double b_P;
double b_y_tmp;
double nrm;
double scale;
int br;
int i;
int ib;
int ic;
int ps_tmp;
boolean_T p;
/*  Reverse optimize and back out the equilibrium returns */
/*  This is formula (12) page 6. */
for (i = 0; i < 7; i++) {
nrm = 0.0;
for (ic = 0; ic < 7; ic++) {
nrm += weq[ic] * sigma[ic + 7 * i];
}
pi[i] = nrm * 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. */
b_y_tmp = 0.0;
b_P = 0.0;
for (i = 0; i < 7; i++) {
nrm = 0.0;
scale = 0.0;
for (ic = 0; ic < 7; ic++) {
absxk = P[ic];
nrm += ts[i + 7 * ic] * absxk;
scale += absxk * ts[ic + 7 * i];
}
y_tmp[i] = scale;
er_tmp[i] = nrm;
nrm = P[i];
b_y_tmp += scale * nrm;
b_P += nrm * pi[i];
}
absxk = 1.0 / (b_y_tmp + Omega);
scale = Q - b_P;
/*  We can also do it the long way to illustrate that d1 + d2 = I */
b_y_tmp = 1.0 / Omega;
/*  Compute posterior estimate of the uncertainty in the mean */
/*  This is a simplified and combined version of formulas (9) and (15) */
nrm = 0.0;
for (i = 0; i < 7; i++) {
er[i] = pi[i] + er_tmp[i] * absxk * scale;
nrm += y_tmp[i] * P[i];
}
absxk = 1.0 / (nrm + Omega);
for (i = 0; i < 7; i++) {
for (ic = 0; ic < 7; ic++) {
b_er_tmp[ic + 7 * i] = er_tmp[ic] * absxk * P[i];
}
}
for (i = 0; i < 7; i++) {
for (ic = 0; ic < 7; ic++) {
nrm = 0.0;
for (ps_tmp = 0; ps_tmp < 7; ps_tmp++) {
nrm += b_er_tmp[i + 7 * ps_tmp] * ts[ps_tmp + 7 * ic];
}
ps_tmp = i + 7 * ic;
ps[ps_tmp] = ts[ps_tmp] - nrm;
}
}
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 (ic = 0; ic < 7; ic++) {
nrm = 0.0;
for (ps_tmp = 0; ps_tmp < 7; ps_tmp++) {
nrm += dv[i + 7 * ps_tmp] * ps[ps_tmp + 7 * ic];
}
ts[i + 7 * ic] = nrm;
}
}
b_P = 0.0;
for (ps_tmp = 0; ps_tmp < 7; ps_tmp++) {
b_P += ts[ps_tmp + 7 * ps_tmp];
}
theta[0] = b_P / 7.0;
for (i = 0; i < 7; i++) {
for (ic = 0; ic < 7; ic++) {
b_er_tmp[ic + 7 * i] = P[ic] * b_y_tmp * P[i];
}
}
for (i = 0; i < 7; i++) {
for (ic = 0; ic < 7; ic++) {
nrm = 0.0;
for (ps_tmp = 0; ps_tmp < 7; ps_tmp++) {
nrm += b_er_tmp[i + 7 * ps_tmp] * ps[ps_tmp + 7 * ic];
}
ts[i + 7 * ic] = nrm;
}
}
b_P = 0.0;
for (ps_tmp = 0; ps_tmp < 7; ps_tmp++) {
b_P += ts[ps_tmp + 7 * ps_tmp];
}
theta[1] = b_P / 7.0;
for (i = 0; i < 7; i++) {
for (ic = 0; ic < 7; ic++) {
b_er_tmp[ic + 7 * i] = P[ic] * b_y_tmp * P[i];
}
}
for (i = 0; i < 7; i++) {
for (ic = 0; ic < 7; ic++) {
nrm = 0.0;
for (ps_tmp = 0; ps_tmp < 7; ps_tmp++) {
nrm += b_er_tmp[i + 7 * ps_tmp] * ps[ps_tmp + 7 * ic];
}
ts[i + 7 * ic] = nrm;
}
}
b_P = 0.0;
for (ps_tmp = 0; ps_tmp < 7; ps_tmp++) {
b_P += ts[ps_tmp + 7 * ps_tmp];
}
theta[2] = b_P / 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++) {
nrm = 0.0;
for (ic = 0; ic < 7; ic++) {
nrm += er[ic] * dv[ic + 7 * i];
}
w[i] = nrm;
}
/*  Compute posterior weights based on uncertainty in mean and covariance */
for (i = 0; i < 49; i++) {
posteriorSigma[i] *= delta;
}
inv(posteriorSigma, dv);
/*  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*. */
for (i = 0; i < 7; i++) {
nrm = 0.0;
for (ic = 0; ic < 7; ic++) {
nrm += pi[ic] * dv[ic + 7 * i];
}
pw[i] = nrm;
er_tmp[i] = P[i];
}
p = true;
for (ps_tmp = 0; ps_tmp < 7; ps_tmp++) {
pi[ps_tmp] = 0.0;
if (p) {
nrm = P[ps_tmp];
if (rtIsInf(nrm) || rtIsNaN(nrm)) {
p = false;
}
} else {
p = false;
}
}
if (!p) {
for (i = 0; i < 7; i++) {
pi[i] = rtNaN;
}
} else {
nrm = 0.0;
scale = 3.3121686421112381E-170;
for (ps_tmp = 0; ps_tmp < 7; ps_tmp++) {
absxk = fabs(P[ps_tmp]);
if (absxk > scale) {
b_P = scale / absxk;
nrm = nrm * b_P * b_P + 1.0;
scale = absxk;
} else {
b_P = absxk / scale;
nrm += b_P * b_P;
}
}
nrm = scale * sqrt(nrm);
if (nrm > 0.0) {
if (P[0] < 0.0) {
absxk = -nrm;
} else {
absxk = nrm;
}
if (fabs(absxk) >= 1.0020841800044864E-292) {
scale = 1.0 / absxk;
for (ps_tmp = 0; ps_tmp < 7; ps_tmp++) {
er_tmp[ps_tmp] *= scale;
}
} else {
for (ps_tmp = 0; ps_tmp < 7; ps_tmp++) {
er_tmp[ps_tmp] /= absxk;
}
}
er_tmp[0]++;
absxk = -absxk;
} else {
absxk = 0.0;
}
for (ps_tmp = 0; ps_tmp < 7; ps_tmp++) {
y_tmp[ps_tmp] = er_tmp[ps_tmp];
}
if (absxk != 0.0) {
for (ps_tmp = 0; ps_tmp < 7; ps_tmp++) {
y_tmp[ps_tmp] = -y_tmp[ps_tmp];
}
y_tmp[0]++;
nrm = fabs(absxk);
scale = absxk / nrm;
absxk = nrm;
for (ps_tmp = 0; ps_tmp < 7; ps_tmp++) {
y_tmp[ps_tmp] *= scale;
}
} else {
for (ps_tmp = 0; ps_tmp < 7; ps_tmp++) {
y_tmp[ps_tmp] = 0.0;
}
y_tmp[0] = 1.0;
}
if (rtIsInf(absxk)) {
scale = rtNaN;
} else if (absxk < 4.4501477170144028E-308) {
scale = 4.94065645841247E-324;
} else {
frexp(absxk, &br);
scale = ldexp(1.0, br - 53);
}
if (absxk > 7.0 * scale) {
scale = 1.0 / absxk;
for (ps_tmp = 0; ps_tmp < 7; ps_tmp++) {
i = ps_tmp + 1;
for (ic = i; ic <= i; ic++) {
pi[ic - 1] = 0.0;
}
}
br = 0;
for (ps_tmp = 0; ps_tmp < 7; ps_tmp++) {
br++;
for (ib = br; ib <= br; ib += 7) {
i = ps_tmp + 1;
for (ic = i; ic <= i; ic++) {
pi[ic - 1] += y_tmp[ib - 1] * scale;
}
}
}
}
}
*lambda = 0.0;
for (ps_tmp = 0; ps_tmp < 7; ps_tmp++) {
*lambda += pi[ps_tmp] * (w[ps_tmp] * (tau + 1.0) - weq[ps_tmp]);
}
}

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