Correct Matrix Inversion in a Black–Litterman Return Update
Article Quant Q&A · Author: user113156
Summary
The document shows an attempted Black–Litterman calculation in R using a covariance matrix, a view-link matrix, view uncertainty, and equilibrium excess returns. The author’s intermediate figures appear close to a tutorial, but the updated-return vector is implausibly large. The accepted response identifies the problem: elementwise exponentiation with a negative power was used where matrix inversion was required.
The correction uses linear-system solvers to invert the scaled covariance matrix and the view-uncertainty matrix, then applies those inverses in the second part of the Black–Litterman update. The resulting vector is much closer to the author’s expectation. This is a focused implementation lesson; the example does not independently verify the input data, assumptions about market weights, or broader model specification.
Key ideas
- Black–Litterman updates require matrix inverses, not elementwise reciprocals of matrix entries.
- The view uncertainty matrix is formed from the link matrix and the asset covariance matrix.
- Use a matrix-solving operation for both the scaled covariance matrix and view uncertainty matrix.
- A coding error in inversion can produce extreme updated-return values even when earlier calculations look plausible.
- The worked example does not validate the model inputs or assumptions.
Tags
Full text
# Black-Litterman computation in R - where am I going wrong?
# Black-Litterman computation in R - where am I going wrong?
I am trying to compute a small Black Litterman model in R. I am following a Youtube video and translating the excel implementation in R.
I have a var cov S matrix
```
INTC AEP AMZN MRK XOM ^GSPC
INTC 0.0119535151 0.0005721887 0.0072352418 0.0016447926 0.0005925077 0.0024795274
AEP 0.0005721887 0.0042225253 0.0008231236 0.0011854049 0.0010758889 0.0011941026
AMZN 0.0072352418 0.0008231236 0.0191091776 0.0009086193 -0.0002442391 0.0017836173
MRK 0.0016447926 0.0011854049 0.0009086193 0.0063486415 0.0009187387 0.0009943984
XOM 0.0005925077 0.0010758889 -0.0002442391 0.0009187387 0.0027747986 0.0009486789
^GSPC 0.0024795274 0.0011941026 0.0017836173 0.0009943984 0.0009486789 0.0012362303
```
Data:
```
S_cov_var <- structure(c(0.0119535151035911, 0.000572188710022071, 0.00723524182537011,
0.00164479256302833, 0.000592507747499871, 0.00247952741729956,
0.000572188710022071, 0.00422252526205478, 0.000823123610432928,
0.00118540486616208, 0.00107588894445389, 0.00119410264013768,
0.00723524182537011, 0.000823123610432928, 0.0191091775682989,
0.000908619322530227, -0.000244239135715373, 0.00178361731695959,
0.00164479256302833, 0.00118540486616208, 0.000908619322530227,
0.00634864154256473, 0.000918738733973792, 0.00099439837734023,
0.000592507747499871, 0.00107588894445389, -0.000244239135715373,
0.000918738733973792, 0.00277479857981738, 0.000948678870995285,
0.00247952741729956, 0.00119410264013768, 0.00178361731695959,
0.00099439837734023, 0.000948678870995285, 0.00123623026419288
), .Dim = c(6L, 6L), .Dimnames = list(c("INTC", "AEP", "AMZN",
"MRK", "XOM", "^GSPC"), c("INTC", "AEP", "AMZN", "MRK", "XOM",
"^GSPC")))
```
I have a link matrix P:
```
INTC AEP AMZN MRK XOM
View 1 0 1 0 0 -1
View 2 1 0 -1 0 0
```
Data:
```
P <- structure(c(0, 1, 1, 0, 0, -1, 0, 0, -1, 0), .Dim = c(2L, 5L), .Dimnames = list(
c("View 1", "View 2"), c("INTC", "AEP", "AMZN", "MRK", "XOM"
)))
```
I compute Omega as:
```
tau = 1
Omega = tau * P %*% S_cov_var[1:5 ,1:5] %*% t(P)
```
I compute the first part of the formula as:
$$((\tau S)^{-1} + P^{T}\Omega^{-1}P)^{-1}$$
```
first <- ((tau * S_cov_var[1:5 ,1:5])^(-1) + (t(P) %*% Omega^(-1) %*% P))^(-1)
```
Then the second part of the formula as:
$$(\tau S)^{-1}\pi + P^{T}\Omega^{-1}Q$$
Data:
```
Q <- c(0.01, 0.0175) # uncertainty about my views
implied_equilib_excess_rets <- structure(c(0.00933950373355221, 0.0031834850342374, 0.00648459638783838,
0.00560398430525973, 0.00578609504932214), .Dim = c(5L, 1L), .Dimnames = list(
c("INTC", "AEP", "AMZN", "MRK", "XOM"), NULL))
```
Calculation:
```
second <- (tau * S_cov_var[1:5 ,1:5])^(-1) %*% implied_equilib_excess_rets[,1] + (t(P) %*% (Omega^(-1)) %*% Q)
```
Which gives me the result (for the second part):
```
[,1]
INTC 12.274655
AEP 21.034321
AMZN -3.885805
MRK 22.681126
XOM 14.381804
```
Which is completely wrong.
I have checked all my figures up until this point and they almost match the video I am following (he uses adjusted Yahoo prices I use Closing prices since the video is a few years old). I expect the results to not match but they do not match by a long way. For example the expected output should be (for the second part)
```
INTC 1.175
AEP 2.304
AMZN -1.074
MRK 0.448
XOM -0.431
```
Minute 11:27 here shows how the second part of the formula should look like.
Additional:
Here is a dump of the R code I have from the excel video (I get pretty close results based on the video output until the `second` part of the code):
```
library(tsibble)
library(tidyverse)
library(tidyquant)
start_date <- "2002-01-01"
end_date <- "2007-08-01"
symbols <- c("INTC", "AEP", "AMZN", "MRK", "XOM", "^GSPC")
portfolio_prices <- tq_get(
symbols,
from = start_date,
to = end_date,
) %>%
select(symbol, date, close)
portfolio_monthly_prices <- portfolio_prices %>%
group_by(symbol) %>%
tq_transmute(
select = close,
mutate_fun = to.period,
period = "months"
) %>%
pivot_wider(names_from = symbol, values_from = close) %>%
tk_xts(., date_var = date)
portfolio_monthly_returns <- portfolio_prices %>%
group_by(symbol) %>%
tq_transmute(
select = close,
mutate_fun = periodReturn,
period = "monthly",
type = "log",
) %>%
pivot_wider(names_from = symbol, values_from = monthly.returns) %>%
tk_xts(., date_var = date)
portfolio_monthly_returns[,1:5]
Asset_Ave_Rets <- colMeans(portfolio_monthly_returns[, 1:5])
Market_Ave_Rets <- colMeans(portfolio_monthly_returns[, 6])
Market_variance <- var(portfolio_monthly_returns[, 6])
obs <- nrow(portfolio_monthly_returns) - 1
S_cov_var <- as.matrix(cov(portfolio_monthly_returns))
Variance <- diag(S_cov_var)
StandardDev <- sqrt(Variance)
lambda = c(1.5, 1.5, 1.5, 1.5, 1.5)
Market_caps <- data.frame(
stock = c("INTC", "AEP", "AMZN", "MRK", "XOM"),
mkt_cap = c(153.42, 19.2, 36.62, 125.5, 505.49)
) %>%
mutate(
market_weights = mkt_cap / sum(mkt_cap)
)
weights <- as.vector(Market_caps$market_weights)
implied_equilib_excess_rets <- 2*c(lambda) * (S_cov_var[1:5, 1:5] %*% weights[1:5]) # AKA pi
implied_equilib_excess_rets
#VIEW 1: AP outperforms exxon mobile by 1% per month
#VIEW 2: Intel outperforms Amazon by 1.75 % per month
Q <- c(0.01, 0.0175)
VIEWS = matrix(data = 0, nrow = 2, ncol = ncol(S_cov_var[,1:5]))
rownames(VIEWS) = c(paste("View", seq(1:2)))
colnames(VIEWS) = colnames(S_cov_var[, 1:5])
# Fill out the link matrix
VIEWS[1, 2] <- 1
VIEWS[1, 5] <- -1
VIEWS[2, 1] <- 1
VIEWS[2, 3] <- -1
P = as.matrix(VIEWS) # link matrix
tau = 1
Omega = tau * P %*% S_cov_var[1:5 ,1:5] %*% t(P) # uncertainty associated with our views
Omega
# black litterman formula
# part 1:
# expected returns calculation
first <- ((tau * S_cov_var[1:5 ,1:5])^(-1) + (t(P) %*% Omega^(-1) %*% P))^(-1)
first
# part 2:
second <- (tau * S_cov_var[1:5 ,1:5])^(-1) %*% implied_equilib_excess_rets[,1] + (t(P) %*% (Omega^(-1)) %*% Q)
second
```
## Answer by AK88 (score 6, accepted)
https://quant.stackexchange.com/a/49921
I had a quick look at your code and it seems like you are applying the `inverse` function incorrectly in your calculations. For example, in the `second` equation I did the following changes:
```
sub1 = tau * S_cov_var[1:5 ,1:5]
Isub1 = solve(sub1)
IOmega = solve(Omega)
second <- (Isub1 %*% implied_equilib_excess_rets[,1] + (t(P) %*% IOmega %*% Q))
```
Which returns the following:
```
> second
[,1]
INTC 1.7555493
AEP 2.4034124
AMZN -1.0770208
MRK 0.4480916
XOM -0.5300325
```
I hope this is what you are looking for. Otherwise let me know.Shown in full with attribution under the source's licence. Licence: CC BY-SA 4.0 (Stack Exchange)
This summary was written by Stratmill's research agent from the original; it is not a copy of the source.