Source code for qf_lib.portfolio_construction.black_litterman.black_litterman

#     Copyright 2016-present CERN – European Organization for Nuclear Research
#
#     Licensed under the Apache License, Version 2.0 (the "License");
#     you may not use this file except in compliance with the License.
#     You may obtain a copy of the License at
#
#         http://www.apache.org/licenses/LICENSE-2.0
#
#     Unless required by applicable law or agreed to in writing, software
#     distributed under the License is distributed on an "AS IS" BASIS,
#     WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied.
#     See the License for the specific language governing permissions and
#     limitations under the License.

from typing import Tuple

import numpy as np
from numpy import power, sqrt
from numpy.linalg import inv

from qf_lib.containers.dataframe.qf_dataframe import QFDataFrame
from qf_lib.containers.series.qf_series import QFSeries


[docs]class BlackLitterman: """ Creates an object allowing calculation of the distribution of the returns using Black Litterman model Parameters ----------- hist_cov: QFDataFrame covariance matrix of historical excess returns of assets; should be annualized weights: QFSeries weights of assets in the market cap index sample_size: int number of simple returns used to estimate the historical covariance sharpe: float average sharpe ratio of the market (by default value of 0.5 is used) """ def __init__(self, hist_cov: QFDataFrame, weights: QFSeries, sample_size: int, sharpe: float = 0.5): self.hist_cov = hist_cov self.weights = weights self.tau = self._caltulate_tau(sample_size) self.sharpe = sharpe # properties used in view implementation self._p = None self._q = None self._omega = None @property def P(self) -> np.array: """ Matrix of views P. nr of rows = number of different views nr of columns = number of assets Example for 4 assets and 2 views: - (0, 1, 0, 0) means go long second asset - (0, 1, 0, -1) means go long second asset while shorting first asset """ return self._p @property def Q(self) -> np.array: """ Vector of expected performance Q. Also sometimes refered as vector v nr of rows = number of different views nr of columns = 1 For example: - [0.10] - [0.03] means that we have two views and we expect to realise 10% above risk free rate on first view and 3% on second """ return self._q @property def Omega(self) -> np.array: """ Matrix of variance of the views Matrix is square and diagonal. Each value Omega[i,i] corresponds to the variance of a specific view nr of rows = number of different views nr of columns = number of different views For example: - (0.0001, 0.0000) means we expect volatility of the first view to be sqrt(0.0001) = 0.01 = 1% - (0.0000, 0.0025) means we expect volatility of the first view to be sqrt(0.0025) = 0.05 = 5% NOTE: self._omega is stored as horizontal vector. Values are copied on diagonal """ if self._omega is not None: N = self._omega.size omega = np.identity(N) for i in range(N): omega[i, i] = self._omega[i] return omega return None
[docs] def add_absolute_view(self, asset_index: int, outperformance: float, view_vol: float): """ Adds a new view that will be taken into account in calculating the posterior Parameters ------------- asset_index: int index of the asset that we expect to outperform the market (indexing starts at 0) outperformance: float how much are we expecting one asset to outperform the market view_vol: float volatility of the view. This is the measure of a standard deviation of the outperformance value. For example: 0.02 means that the real outperformance will be provided outperformance +- 2% within one standard deviation confidence interval. """ new_p_row = np.zeros((1, self.weights.size)) new_p_row[0, asset_index] = 1 new_q_elem = np.ones((1, 1)) * outperformance new_omega_elem = np.ones(1) * view_vol * view_vol # Omega corresponds to variance of view so use the vol^2 self._update_views(new_p_row, new_q_elem, new_omega_elem)
[docs] def add_relative_view(self, outperforming_asset_index: int, underperforming_asset_index: int, outperformance: float, view_vol: float): """ Adds a new view that will be taken into account in calculating the posterior The investor believs that outperforming_asset will outperform the underperforming_asset by outperformance with the volatility of the view of view_vol Parameters ---------- outperforming_asset_index: int index of the asset that you believe will outperform (indexing starts at 0) underperforming_asset_index: int index of the asset that you believe will underperform (indexing starts at 0) outperformance: float how much are we expecting one asset to outperform the other view_vol: float volatility of the outperformance. This is the measure of a standard deviation of the outperformance value. For example: 0.02 means that the real outperformance will be provided outperformance +- 2% within one standard deviation confidence interval. """ new_p_row = np.zeros((1, self.weights.size)) new_p_row[0, outperforming_asset_index] = 1 new_p_row[0, underperforming_asset_index] = -1 new_q_elem = np.ones((1, 1)) * outperformance new_omega_elem = np.ones(1) * view_vol * view_vol # Omega corresponds to variance of view so use the vol^2 self._update_views(new_p_row, new_q_elem, new_omega_elem)
[docs] def calculate_prior(self) -> Tuple[QFSeries, QFDataFrame]: """ Function calculates the prior for the BL model: - prior_mean = lambda * hist_cov * market_cap_weights - prior_cov = tau * hist_cov Returns -------- Tuple[QFSeries, QFDataFrame] Tuple(prior_mean, prior_cov) """ lambda_ = self.calculate_lambda() prior_mean = lambda_ * self.hist_cov.dot(self.weights) prior_cov = self.tau * self.hist_cov return prior_mean, prior_cov
[docs] def calculate_posterior(self) -> Tuple[QFSeries, QFDataFrame]: """ Calculate shape of posterior distribution including the views. It is calculated using a numerically stable formula posterior_mean = prior_mean + tau * hist_cov * P' * ( tau * P * hist_cov * P' + Omega)^(-1) * (Q - P * prior_mean) posterior_cov = (1 + tau) * hist_cov - tau^2 * hist_cov * P' * (tau * P * hist_cov * P' + Omega)^(-1) * P * hist_cov Returns -------- Tuple[QFSeries, QFDataFrame] """ prior_mean, prior_sigma = self.calculate_prior() if self.P is not None: # calculate the posterior according to the formula above inverse_mat = self.tau * self.P.dot(self.hist_cov).dot(self.P.T) + self.Omega inverse_mat = inv(inverse_mat) P_dot_mean = self.P.dot(prior_mean).reshape(self.Q.shape) # P * prior_mean with a correct shape. # for some reason numpy calculates [3x9] * [1x9] -> [1x3], but it should be [3x1] instead !!! posterior_mean_update = self.tau * self.hist_cov.dot(self.P.T).dot(inverse_mat).dot(self.Q - P_dot_mean) posterior_mean_update = posterior_mean_update.squeeze() posterior_cov_update = power(self.tau, 2) * self.hist_cov.dot( self.P.T.dot(inverse_mat).dot(self.P).dot(self.hist_cov)) # set the column names that disappeared in the calculation posterior_cov_update.columns = posterior_cov_update.index else: # There are no views. posterior_mean_update = 0 posterior_cov_update = 0 posterior_mean = prior_mean + posterior_mean_update posterior_cov = (1 + self.tau) * self.hist_cov - posterior_cov_update return posterior_mean, posterior_cov
def _update_views(self, new_p_row, new_q_elem, new_omega_elem): if self.P is None: self._p = new_p_row self._q = new_q_elem self._omega = new_omega_elem else: self._p = np.concatenate((self._p, new_p_row), axis=0) self._q = np.concatenate((self._q, new_q_elem), axis=0) self._omega = np.concatenate((self._omega, new_omega_elem), axis=0)
[docs] def calculate_lambda(self): """ Calculates lambda corresponding to risk aversion in the formula: max(w): w'R - lambda/2 * w' * Cov * w: lambda is calculated in the followign way: lambda = sharpe_ratio_of_market / vol_of_market. """ vol = self._calculate_volatility() lambda_ = self.sharpe / vol return lambda_
def _calculate_volatility(self): # calculate historical volatility of the portfolio with market weights vol = sqrt(self.weights.dot(self.hist_cov.dot(self.weights.T))) return vol def _caltulate_tau(self, len_of_returns: int): """ In many implementations it approximated by 1 / nr of observations, which corresponds to standard error That is the error of estimation of the covariance matrix Parameters ----------- len_of_returns: int number of returns used to calculate the historical covariance matrix Returns -------- float coefficient of uncertainty in the prior estimate of the mean """ return 1 / len_of_returns def print(self): print("=====> VOL") print(self._calculate_volatility()) print("=====> TAU") print(self.tau) print("=====> P") print(self.P) print("=====> Q") print(self.Q) print("=====> OMEGA") print(self.Omega)