Find k of n assets that "minimize" the correlation matrix

Question

I'm trying to find an efficient way to select $k$ from $n$ risky assets that are the least correlated with each other. I know that I can perform a brute-force search of all $k$-sized combinations of the $n$ assets but this doesn't scale as $n$ grows so I'm wondering if there's an optimization for this problem.

For example, you have the following candidate 10 assets SHW, GOOG, AMZN, WMT, XOM, JNJ, UPS, AMT, AAPL, and NEE. What python, matlab, or R code (you pick) would you run to collect the daily returns for each asset and find the subset of 5 assets that minimizes the square root of the sum of the squares of the entries of the correlation matrix of those 5 assets' returns.

I think this is a binary integer programming problem (or maybe convex optimization since we're searching for a minimum?), please correct the following formulation if it's wrong.

Given

1. a returns matrix $R$ where $\left(r_{ij}\right) \in \mathbb{R}^{m \times n}$ is the return for the $i$-th day and the $j$-th asset and
2. the correlation matrix $C = \text{corr}\left( R \right)$ where $\left(c_{ij}\right) \in \mathbb{R}^{n \times n}$ is the correlation coefficient between the $i$-th and $j$-th assets

We want to find $\vec{x} \in \left\{0,1\right\}^n $ s.t.

1. $\sum \limits_i^n x_i = k$ and
2. $\vec{x}$ minimizes $\sqrt{\sum \limits_{i,j}^n {c'}_{ij}^{2}}$, where $\left({c'}_{ij}\right) \in \mathbb{R}^{n \times n}$ is the entry from the modified correlation matrix $C'$ for assets selected by $\vec{x}$, given by $C'= \left(\vec{x} \otimes \vec{x}\right) \odot C$.

That is to say, $C'$ is $C$ with the rows and columns of the rejected assets "zeroed-out", given by the outer product of $\vec{x}$ with itself (to get a matrix $X$ with $x_{i,j} \in \left\{0,1\right\} = 0$ when $\vec{x}_i = 0$ and $x_{i,j} = 1$ when $\vec{x}_i = 1$) Hadamard multiplied by $C$. Finally I chose the square root of the sum of the squares of the entries of $C'$ but I think any distance metric would do (like the average of the absolute values of the entries).

Thanks.

Answer 1

(I take it that 5 out of 10 assets is just an example, because in this case all combinations could easily be checked.)

Here would be an example how to do it in R with an algorithm called Threshold Accepting.

library("neighbours")  ## https://github.com/enricoschumann/neighbours
library("NMOF")        ## https://github.com/enricoschumann/NMOF

I create a correlation matrix C from random returns. (Note that I use few observations, so that the correlations differ substantially.)

C <- cor(randomReturns(na = 10, ns = 20, rho = 0.5, sd = 0.01))
##        [,1]  [,2]  [,3]  [,4]  [,5]  [,6]  [,7]  [,8]  [,9] [,10]
##  [1,] 1.000 0.390 0.715 0.645 0.621 0.587 0.701 0.437 0.394 0.219
##  [2,] 0.390 1.000 0.359 0.447 0.496 0.725 0.499 0.218 0.469 0.244
##  [3,] 0.715 0.359 1.000 0.570 0.592 0.506 0.427 0.602 0.495 0.326
##  [4,] 0.645 0.447 0.570 1.000 0.785 0.798 0.646 0.419 0.442 0.655
##  [5,] 0.621 0.496 0.592 0.785 1.000 0.677 0.474 0.319 0.414 0.471
##  [6,] 0.587 0.725 0.506 0.798 0.677 1.000 0.539 0.370 0.378 0.559
##  [7,] 0.701 0.499 0.427 0.646 0.474 0.539 1.000 0.324 0.417 0.433
##  [8,] 0.437 0.218 0.602 0.419 0.319 0.370 0.324 1.000 0.450 0.295
##  [9,] 0.394 0.469 0.495 0.442 0.414 0.378 0.417 0.450 1.000 0.295
## [10,] 0.219 0.244 0.326 0.655 0.471 0.559 0.433 0.295 0.295 1.000

The objective function: Given a logical vector x and a correlation matrix, compute the square root of the sum of the squared correlations.

mean_cor <- function(x, C) {
    tmp <- C[x, x][lower.tri(C[x, x])]
    sqrt(sum(tmp*tmp))
}

The so-called neighbourhood function. It will change a given solution by randomly changing one asset.

nb <- neighbourfun(type = "logical", kmin = 5, kmax = 5)

Running the algorithm.

x <- TAopt(mean_cor,
           list(x0 = rep(c(TRUE, FALSE), 5),
                neighbour = nb,
                nI = 200),
           C = C)$xbest

The solution.

which(x)
## [1]  1  2  8  9 10

C[x, x]
##        [,1]  [,2]  [,3]  [,4]  [,5]
## [1,] 1.000 0.390 0.437 0.394 0.219
## [2,] 0.390 1.000 0.218 0.469 0.244
## [3,] 0.437 0.218 1.000 0.450 0.295
## [4,] 0.394 0.469 0.450 1.000 0.295
## [5,] 0.219 0.244 0.295 0.295 1.000

More about the algorithm is in this tutorial at SSRN. (Disclosure: I am the maintainer of packages NMOF and neighbours that I used in the example.)

Answer 2

This is indeed a convex, quadratic, integer programming problem. Most of those are NP-hard, so don't beat your brains to find the optimum efficiently. That said, nowadays these problems can be solved efficiently, approximately but with high accuracy with convex solvers (e.g. Mosek). From my own experience, this could work up to sizes of thousands of assets, say, with selected subportfolios of sizes in the hundreds or less on standard retail hardware.

But there are really simple heuristic solutions as well. Those are very easy to implement, work for almost any size of $n$ or $k$, and may give you a "good enough" solution.

Some ideas for heuristics:

• Trial and error: Just sample randomly $k$ out of the $n$ assets. Calculate total variance, repeat until you loose patience, retain the best.
• Greedy bottom-up: Start with the least variance asset, find the next among the remaining $n-1$ which gives you smallest variance for the pair, select the third from the remaining $n-2$ and so on.
• Greedy top-down: same as bottom-up, but now you drop an asset in each step.
• Relaxation: Drop the constraint of the weights being zero-one and just find the minimum variance pottfolio. Retain the $k$ assets with largest weight. I think this is what @Kermittfrog had in mind.

I am certain there are more possibilities (the algorithm by @Enrico Schumann, genetic/evolutionary algorithms, simulated annealing, ...).

See also this related question of mine on MSE.