Minimum Volatility Portfolio (Python Code: scipy.optimize & cvxopt)

Minimum Volatility Portfolio (MVP: 최소분산 포트폴리오)는

이론적으로는 간단합니다.

정해진 종목들을 바탕으로

포트폴리오의 변동성을 최대한 낮게 만듭니다.

포트폴리오 변동성은 σ(p)^2로 나타나며

행렬의 형태로 표현하면 w'Ωw 로 나타납니다.

(w는 종목당 비중, Ω는 variance covariance matrix 입니다.)

이를 계산하는 패키지는 두가지가 있습니다.

scipy.optimize의 minimize 패키지를 이용하는 방법과,

cvxopt 패키지를 이용하는 방법

두가지를 모두 보도록 합니다.

먼저 scipy.optimize 패키지를 이용합니다.

사실 지난번에 Risk Parity [Link] 에서 봤던 방법과 매우 비슷합니다.

먼저 필요한 패키지들을 쭉쭉쭉 불러옵니다.
핵심이 되는건
from scipy.optimize import minimize
요 녀석 입니다.

from pandas import Series, DataFrame
import pandas as pd
from pandas_datareader import data
from pandas.tseries.offsets import Day, MonthEnd

import numpy as np
import sys

import matplotlib.pyplot as plt

from scipy.stats import rankdata
from scipy.stats import stats
from scipy.optimize import minimize

먼저 목적함수 입니다.

포트폴리오의 variance는 위에서 말했던 것 처럼

w'Ωw로 나타낼 수 있습니다.

포트폴리오의 sigma(standard deviation)은 이 값의 루트 값이겠죠

우리의 목적은 이 포트폴리오의 sigma (혹은 variance 자체)를 

가장 낮게 만드는 겁니다.

def MinVol_objective(x) :
        
        variance = x.T @ covmat @ x
        sigma = variance ** 0.5
        return (sigma)

이번에는 제약함수 입니다.

x값들의 합이 1이 되게 만들어 줍니다.

sum - 1 = 0 이니

sum = 1이 되는 거겠죠?

def weight_sum_constraint(x) :
        return(x.sum() - 1.0 )

최적화를 구하는 함수입니다.

입력값은 covmat(분산-공분산행렬), lb(최소비중), ub(최대비중)

3가지 값이 필요합니다.

x0는 초기값 셋팅으로 일단 동일비중 만큼 줍니다.

lbound는 최소비중이며,

lb에 입력된 값들을 종목수 만큼 생성합니다.

만일 최소비중을 0을 입력하였고, 종목이 5개 라면

repeat 함수에 의해 lbound는 [0, 0, 0, 0, 0] 이 됩니다.

uboud는 최대비중이며, 만일 1을 입력했다면

[1, 1, 1, 1, 1] 이 되겠죠?

bnds는 위에서 입력한 최소비중과 최대비중을 통해

각 종목별 최소, 최대 비중 짝을 만들어 줍니다.

[0, 0, 0, 0, 0] 과 [1, 1, 1, 1, 1] 에 zip 함수를 이용하면

[0, 1], [0, 1], [0, 1], [0, 1], [0, 1]

이러한 형태로 짝이 지어집니다.

이를 튜플형태로 변형해 줍니다.

constraints는 제약조건입니다.

종목별 비중이 0이상 인 것은 이미 bnds를 통해 셋팅했으니

종목의 합이 1인 제약조건만 필요합니다.

sum weight = 1 이므로

'eq' 로 지정해 줘야 겠죠?

option은 tolerance level과 max iteration 횟수입니다.

위에서 입력한 변수들을 쓰까쓰까 입력하여

minimize 명령어를 통해 값을 구합니다.

Risk Parity 구하는 법을 완벽히 이해했다면,

이정도 목적함수는 너무 쉽죠? :D

def MinVol(covmat, lb, ub) :
    x0 = np.repeat(1/covmat.shape[1], covmat.shape[1]) 
    lbound  = np.repeat(lb, covmat.shape[1])
    ubound  = np.repeat(ub, covmat.shape[1])
    bnds = tuple(zip(lbound, ubound))
    
    constraints = ({'type': 'eq', 'fun': weight_sum_constraint})
    options = {'ftol': 1e-20, 'maxiter': 800}
    
    result = minimize(fun = MinVol_objective,
                      x0 = x0,
                      method = 'SLSQP',
                      constraints = constraints,
                      options = options,
                      bounds = bnds)
    return(result.x)

두번째는 CVXOPT 패키지를 이용합니다.

해당 패키지는 윈도우 / 파이썬 3.6 버젼에는 지원이 안됩니다.

3.4 버젼에서 실행하든지, 2.7 버젼에서 실행하든지 해야합니다.

cvxopt의 solvers.qp 설명서는

강의교재 [Link] 를 보면 매우 쉽게 이해됩니다.

R의 solve.qp를 사용해보신분은 사실 거의 동일하여

매우 쉽게 이해가 되실 겁니다.

equality 와 inequality condition이 R에서는 같은 매트릭스에 쓰는데

python에서는 따로 쓰는 차이 밖에 없습니다.

먼저 패키지를 불러 옵니다.

cvxopt 에서 matrix와 solvers가 필요합니다.

from cvxopt import matrix
from cvxopt import solvers

아래식이 해당 패키지의 quadratic programming 식입니다.

x가 결과값이므로 우리가 입력해야 할 값은

P, q, G, h, A, b 입니다.

먼저 목적함수인 min 1/2 x'Px 를 봅니다.

1/2는 생략해도 되며, x는 결과값이므로

우리가 입력해야 하는 값은 P 입니다.

우리가 최소화 시키고자 하는 값은 포트폴리오의 변동성이고

포트폴리오의 변동성은 w'Ωw 로 나타낼 수 있습니다.

따라서 P에는 오메가 값인 variance-covariance 매트릭스를 입력하면 됩니다.

여기서 q는 목적식에서 필요가 없으므로, 0을 넣습니다.

P = matrix(np.array(covmat), tc = 'd')
q = matrix(np.zeros(len(covmat)), tc = 'd')

제약함수 중에 첫번째인 Gx <= h 를 보도록 합니다.
(보통 최적화 프로그램과 부등호가 반대네요...)

여기에는 lower bound와 upper bound를 각각 입력해줄 수 있습니다.

아래처럼 입력하면, G에는 왼쪽 매트릭스를, h는 오른쪽 매트릭스를 입력합니다.

위에는 lower bound 제약식이며, 아래는 upper bound 제약식이죠

행렬을 계산하면 아래와 같이 정리됩니다.
lb 보다는 크고, ub 바운드 보다는 작게 정리가 되지요?

이를 코드로 표현하면 아래와 같습니다.

lb_diag는 -1의 diagonal matrix이며,

ub_diag는 1의 diagonal matrix 입니다.

concatenate 명령어를 통해 위의 두 함수를 포개주며,

매트릭스 형태로 바꿔줍니다.

h역시 -lb와 ub를 종목수 만큼 반복해주며

concatenate를 통해 포개줍니다.

lb_diag = np.diag(np.repeat(-1, len(covmat)))
ub_diag = np.diag(np.repeat(1, len(covmat)))
    
G = matrix(np.concatenate((lb_diag, ub_diag)) , tc = 'd')
h = matrix(np.concatenate((np.repeat(-lb, len(covmat)), np.repeat(ub, len(covmat)))), tc = 'd')

마지막으로 Ax = b 는 등호가 있는 제약조건입니다.

여기는 비중이 합이 1인 제약조건을 입력하면 됩니다.

아래 붉게 표시된 매트릭스가 A이며, [1] 이 b 입니다.

이를 정리하면 모든 비중의 합 = 1의 형태로 표현이 됩니다.

이를 코드로 표현하면 아래와 같습니다.
A는 1을 종목수만큼 반복한 후 Transpose 해주며,

b에는 그냥 1을 입력해 줍니다.

둘다 역시 매트릭스 형태로 바꿔주고요

A = matrix(np.repeat(1, len(covmat)), tc = 'd').T
b = matrix(1, tc = 'd')

solvers.qp 함수에 위에서 만든 변수들을 죄다 집어넣으면

결과값이 나옵니다.

solution = solvers.qp(P, q, G, h, A, b)
solution['x']))

이를 함수의 형태로 만들면 다음과 같습니다.
분산-공분산 행렬, 최소값, 최대값을 입력하면
결과값으로 비중을 반환하는 형태입니다.

def MinVol(covmat, lb, ub) :
       
    P = matrix(np.array(covmat), tc = 'd')
    q = matrix(np.zeros(len(covmat)), tc = 'd')
    
    lb_diag = np.diag(np.repeat(-1, len(covmat)))
    ub_diag = np.diag(np.repeat(1, len(covmat)))
    
    G = matrix(np.concatenate((lb_diag, ub_diag)) , tc = 'd')
    h = matrix(np.concatenate((np.repeat(-lb, len(covmat)), np.repeat(ub, len(covmat)))), tc = 'd')
    A = matrix(np.repeat(1, len(covmat)), tc = 'd').T
    b = matrix(1, tc = 'd')
    
    solution = solvers.qp(P, q, G, h, A, b)
    return(np.array(solution['x']))

그럼 위에서 구한 두가지 방법의 결과값을 비교해 봅니다.

	cvxopt	minimize	Diff_abs
1	15.27%	15.40%	0.13%
2	0.04%	0.00%	0.04%
3	0.10%	0.00%	0.10%
4	0.01%	0.00%	0.01%
5	0.00%	0.00%	0.00%
6	77.92%	78.10%	0.18%
7	0.00%	0.00%	0.00%
8	0.07%	0.00%	0.07%
9	0.25%	0.00%	0.25%
10	6.34%	6.50%	0.16%
sum	100%	100%	0.93%

cvxopt(solvers.qp)가 매우 미세한 값까지 찾기는 하지만

크게 차이는 없습니다..

다만 문제는, 0에 가까운 종목이 많이 나오는 코너해 문제가 발생합니다.

이번에는 lower bound를 5%, upper bound를 15%를 준 후 값을 구합니다.

	cvxopt	minimize	Diff_abs
1	13.67%	14.80%	1.13%
2	5.00%	5.00%	0.00%
3	6.33%	5.20%	1.13%
4	5.00%	5.00%	0.00%
5	15.00%	15.00%	0.00%
6	15.00%	15.00%	0.00%
7	5.00%	5.00%	0.00%
8	5.01%	5.00%	0.01%
9	15.00%	15.00%	0.00%
10	14.99%	15.00%	0.01%
sum	100%	100%	2.29%

둘 간의 차이가 약간 벌어지기는 하지만

역시나 크게 차이는 없습니다.

비중이 가장 작은 종목은 5%, 가장 큰 종목은 15%가 나와

제약식 역시 잘 작동하는 모습이 확인됩니다.

두패키지 중 어느걸 사용할 지는 개취인듯 합니다.

저는 3.6에서도 잘 작동하는 scipy.optimize 가 더 좋겠네요 :D