넘파이 따라하면서 쉽게 배우자

넘파이 따라하면서 쉽게 배우자

2020, Jan 02    
intro_numpy

Introduction to NumPy

by JunPyo Park

Part of the Alpha Square Lecture Series:


넘파이 관련 요약 정리

In [1]:
import numpy as np
import matplotlib.pyplot as plt

Basic NumPy Arrays(어레이)

금융 데이터 분석에서 numpy를 활용한 가장 단순한 예제는 다음과 같이 몇몇 주식들의 기간 수익률이 주어졌을 때 평균을 구하는것 입니다.

In [3]:
stock_list = [3.5, 5, 2, 8, 4.2]

np.array() 함수를 활용하여 리스트를 어레이로 만들 수 있습니다.

In [4]:
returns = np.array(stock_list)
print(returns, type(returns))
[3.5 5.  2.  8.  4.2] <class 'numpy.ndarray'>

타입명을 보시면 그냥 array가 아닌 ndarray 라고 되어있습니다. 이는 NumPy 어레이가 차원(multiple dimensions)을 가질수 있기 때문입니다. np.array() 함수에 리스트 안에 리스트를 넣으면 2차원 어레이가 만들어집니다. 만약 리스트 안에 리스트 안에 리스트를 넣으면 3차원 어레이가 만들어집니다. 이와 같이 차원을 확장시켜 n차원 어레이를 만들 수 있습니다.

In [5]:
A = np.array([[1, 2], [3, 4]])
print(A, type(A))
[[1 2]
 [3 4]] <class 'numpy.ndarray'>

어레이의 차원을 확인하기 위해서는 shape 메서드를 통해 접근이 가능합니다.

In [6]:
print(A.shape)
(2, 2)

위의 returns 어레이는 리스트와 비슷한 방식으로 인덱싱이 되어있습니다. 0에서 시작해서 n-1 번재에서 끝나는 방식입니다. 여기서 n은 어레이의 길이를 의미합니다.

In [7]:
print(returns[0], returns[len(returns)-1])
3.5 4.2

리스트 슬라이싱의 방법을 똑같이 적용하여 원하는 데이터에 접근이 가능합니다.

In [8]:
print(returns[1:3])
[5. 2.]

이차원 이상의 어레이에서도 이런 슬라이싱 기법을 확장하여 적용이 가능합니다.

In [9]:
# 첫번째 열의 데이터를 출력합니다.
print(A[:,0])
[1 3]
In [10]:
# 첫번째 행의 데이터를 출력합니다.
print(A[0,:])
[1 2]

위와 같이 슬라이싱을 통해 얻은 값도 array 타입을 가지게 됩니다.

In [11]:
print(type(A[0,:]))
<class 'numpy.ndarray'>

다음과 같은 방식으로도 행에 접근이 가능합니다.

In [12]:
print(A[0])
[1 2]

특정 위치의 데이터에 접근하고자 한다면 다음과 같이 행과 열의 인덱스 번호를 통해 접근이 가능합니다.

In [13]:
print(A[1,1])
4

어레이 함수(Array Function)

넘파이에서 만들어진 대부분의 함수들은 입력값으로 어레이를 사용할 수 있습니다. 예를 들어 np.log()함수는 다음과 같이 어레이에 적용이 가능합니다.

In [14]:
print(np.log(returns))
[1.25276297 1.60943791 0.69314718 2.07944154 1.43508453]

returns 어레이에 들어있는 각 데이터 값에 자연로그(고등학교 수학시간에 배웠던 $\text{ln}\,x$ 입니다)가 일괄적으로 적용되었습니다.

mean()max() 같이 어레이를 입력받아 하나의 값만 출력하는 함수들도 존재합니다.

In [15]:
# 입력받은 어레이의 평균을 계산합니다.
print(np.mean(returns))
4.54
In [16]:
# 입력받은 어레이에서 가장 큰 값을 계산합니다.
print(np.max(returns))
8.0

더 많은 함수들에 대한 설명은 공식 Documentation을 확인해 주세요.

다시 returns 어레이로 돌아와 봅시다

처음에 만들었던 1차원 어레이인 returns에 숫자를 더하면 또는 곱하면 어떻게 될까요? 다음과 같이 어레이의 모든 데이터에 일괄적으로 더해지거나 곱해집니다.

In [17]:
returns
Out[17]:
array([3.5, 5. , 2. , 8. , 4.2])
In [18]:
returns + 100
Out[18]:
array([103.5, 105. , 102. , 108. , 104.2])
In [19]:
returns * 10
Out[19]:
array([35., 50., 20., 80., 42.])
In [20]:
returns * 2 + 5
Out[20]:
array([12. , 15. ,  9. , 21. , 13.4])

다음과 같이 평균과 표준편차를 구할 수 있습니다.

In [21]:
print("Mean: ", np.mean(returns), "Std Dev: ", np.std(returns))
Mean:  4.54 Std Dev:  1.9915822855207364

이제 NumPy 함수들을 사용하여 가상의 자산 포트폴리오를 만들어 보도록 하겠습니다. returns 라는 어레이를 만들어 종목 마다 월별 가격 변동률을 저장합니다. 이 때 월별 가격 변동은 평균이 1.01이고 표준편차가 0.03인 정규분포를 따른다고 가정하도록 합니다. 다음과 같이 코드를 작성하여 100개월 간의 가격변동을 나타낼 수 있습니다.

In [22]:
R_base = np.random.normal(1.01, 0.03, 100)
returns = R_base
returns
Out[22]:
array([1.0043993 , 1.02186788, 0.98120378, 1.0168387 , 1.00910664,
       1.01251239, 1.07151284, 0.98287918, 0.98185529, 1.06192513,
       1.01079382, 1.05248844, 1.00783297, 1.02339855, 1.00504771,
       0.99045917, 1.071858  , 1.03403789, 1.05072032, 0.98482139,
       1.01145033, 1.05627951, 0.98383025, 0.97562154, 0.96302575,
       0.98630084, 1.01510072, 1.03679348, 1.01373306, 0.99987467,
       1.02561501, 1.03926775, 1.03513315, 1.03266337, 1.01516756,
       0.97869877, 1.02777437, 1.0023664 , 0.99050348, 0.99198994,
       0.9923167 , 1.01239267, 1.05277046, 1.0876452 , 0.98027066,
       0.99700062, 1.00553525, 1.02606698, 1.02964694, 1.03582226,
       0.98266918, 0.97667876, 1.0369319 , 1.01546044, 1.03973852,
       0.99802166, 1.01295953, 1.01387421, 1.01648131, 0.98881554,
       1.06585671, 1.04646783, 1.01919849, 1.02391643, 0.92899265,
       1.02838783, 1.06542481, 1.03070301, 1.03910764, 0.98480794,
       1.02369659, 1.02073538, 1.06196596, 0.99965963, 0.99446137,
       1.03578263, 1.01423252, 0.98328582, 1.02711772, 1.0125088 ,
       0.99553255, 0.95413701, 0.97384931, 1.01703512, 0.98955365,
       1.07523977, 0.9919207 , 1.02766834, 1.00623367, 1.00244756,
       0.96629972, 0.98504908, 1.02144493, 1.00960517, 1.02978203,
       0.97729322, 1.02584354, 1.02613512, 0.96691902, 1.0104309 ])

NumPy의 random(랜덤)모듈은 이런 난수 값을 발생시킬때 유용합니다. 랜덤 모듈에는 다양한 샘플링을 위한 확률 분포(Probability Distribution)가 내장되어있습니다. 자세한 내용은 확률 변수(Random Variables)를 다루는 강의에서 설명하도록 하겠습니다. 여기에서 $R_{base}$는 평균이 $1.01$이고 표준편차가 $0.03$인 정규분포에서 100개의 샘플을 저장한 어레이(벡터)라고 이해하시면 됩니다.

위 출력 결과를 해석 하면 다음과 같습니다. 첫째달은 전월대비 가격이 0.43993% 상승하였고 둘째달은 전월대비 2.186788% 상승 ... 마지막인 100개월 이후에는 전월 대비 가격이 1.04309% 상승한 것 입니다.

누적 가격 변동은 월별 가격 변동을 누적하여 곱해주면 됩니다. np.cumprod()함수를 사용하여 이런 누적 곱셈에 대한 계산이 손쉽게 가능합니다.

In [23]:
np.cumprod([1,2,3,4,5])
Out[23]:
array([  1,   2,   6,  24, 120], dtype=int32)
In [24]:
assets = np.cumprod(R_base)
assets
Out[24]:
array([1.0043993 , 1.02636338, 1.00707163, 1.0240294 , 1.03335487,
       1.04628461, 1.12110739, 1.10191311, 1.08191921, 1.1489172 ,
       1.1613184 , 1.22227419, 1.23184823, 1.2606717 , 1.26703521,
       1.25494664, 1.3451246 , 1.39090979, 1.46145718, 1.43927429,
       1.45575447, 1.53768362, 1.51281966, 1.47593945, 1.4213677 ,
       1.40189616, 1.4230658 , 1.47542535, 1.49568745, 1.4955    ,
       1.53380724, 1.5940364 , 1.65003992, 1.70393578, 1.72978032,
       1.69293388, 1.73995405, 1.74407147, 1.72750887, 1.71367143,
       1.70050478, 1.72157858, 1.81242706, 1.9712776 , 1.9323856 ,
       1.92658964, 1.9372538 , 1.98775216, 2.04668294, 2.11999976,
       2.08325842, 2.03467426, 2.10981865, 2.14243737, 2.22757467,
       2.22316777, 2.25197898, 2.28322341, 2.32085392, 2.29489643,
       2.44603076, 2.55969249, 2.60883472, 2.67122872, 2.48155184,
       2.5519977 , 2.71896167, 2.80244197, 2.91203888, 2.86779902,
       2.93575606, 2.99663008, 3.18231915, 3.18123598, 3.1636163 ,
       3.27681881, 3.3234562 , 3.26790734, 3.35652554, 3.39851164,
       3.38332897, 3.22815938, 3.14374079, 3.19729479, 3.16389474,
       3.40194544, 3.3744601 , 3.46782581, 3.4894431 , 3.49798371,
       3.38010067, 3.32956504, 3.40096733, 3.43363418, 3.53589477,
       3.45560597, 3.54491105, 3.63755771, 3.51722373, 3.55391155])

위 결과를 통해 100개월 뒤 주식의 가격은 처음 시점 대비 255.391155% 상승했음을 알 수 있습니다.

이 아이디어를 확장하여 주식 $10$종목과 $100$개월 동안의 기간에 해당하는 데이터를 가상으로 만들어 보도록 하겠습니다. 다음과 같이 종목 수 인 $N$ 값에 10을 지정하고 100개월간 데이터를 저장하기 위해 2차원 어레이를 생성합니다. 다음과 같이 np.zeros(shape)를 사용하여 모든 값이 0인 어레이를 원하는 shape에 맞게 만들 수 있습니다. 여기서 $N$ 값이 10이기 때문에 모든 값이 0으로 초기화 된 10행 100열 짜리 2차원 어레이($N\times100$)가 생성됩니다.

In [27]:
N = 10
returns = np.zeros((N,100))
assets = np.zeros((N,100))

첫 행의 데이터(첫번째 종목)를 위에서 계산하였던 값으로 채워줍니다.

In [28]:
returns[0] = R_base
assets[0] = np.cumprod(R_base)

나머지 9종목 각각에 대한 100개월치 데이터를 채우기 위해 for 루프를 사용하여 값을 채워줍니다. 이 때 나머지 9개의 자산이 첫 번째 자산과 가격움직임이 유사하도록(correlated) 설정하기 위해 다음과 같이 처음에 만든 $R_{base}$ 벡터에 노이즈 움직임을 나타내는 $R_i$(Random Noise) 벡터를 더해 다음과 같이 값을 채워줍니다.

In [29]:
# 처음 자산에 대한 데이터는 R_base를 활용하여 채웠기 때문에
# 1부터 시작합니다.
for i in range(1, N):
    
    # 자산가격 움직임이 유사하도록 만들기 위해 다음과 같이 R_i 벡터를 설정합니다.
    R_i = R_base + np.random.normal(0.001, 0.02, 100)
    # i번째 행(i번째 자산)을 채워줍니다.   
    returns[i] = R_i
    assets[i] = np.cumprod(R_i)

# 각 자산별 평균 수익률과 표준편차를 계산합니다.
mean_returns = [(np.mean(R) - 1)*100 for R in returns]
return_volatilities = [np.std(R) for R in returns]

자산별 평균 수익률을 plotting 해보면 다음과 같습니다.

In [30]:
plt.bar(np.arange(len(mean_returns)), mean_returns)
plt.xlabel('Stock')
plt.ylabel('Returns')
plt.title('Returns for {0} Random Assets'.format(N));

기대 수익률 계산하기

지금 까지 여러 자산의 가격 변동 및 수익률에 대한 데이터를 생성해 보았습니다. 이제 이 자산들을 포트폴리오에 담아 기대 수익률을 계산해 보도록 하겠습니다. 먼저 $N$개의 자산들 각각의 가중치(weight)를 생성하여 봅시다.

In [31]:
weights = np.random.uniform(0,1,N)
weights = weights/np.sum(weights)

위의 weights 벡터는 각 자산별 보유 비중을 의미합니다. 0과 1사이 에서 난수를 N개 발생 시켜 벡터를 만든 뒤, 총합이 1이 되도록 만들어 주기 위해 np.sum(weights)값으로 각 원소들을 나누어 줍니다.

In [32]:
weights
Out[32]:
array([0.04811237, 0.03987759, 0.0628113 , 0.06317026, 0.18560434,
       0.09190413, 0.0807558 , 0.11033189, 0.17573278, 0.14169953])
In [33]:
weights.sum()
Out[33]:
1.0

위와 같이 10개의 자산에 대한 보유 비중이 랜덤으로 설정되었으며 총합이 1로 맞추어 졌습니다.

포트폴리오의 기대 수익률을 계산하려면 위에서 만든 수익률 데이터에 가중치 값을 곱하여 더해 주어야 합니다. 이 계산을 하기위해 반복문인 for문을 사용할 수도 있지만 numpydot()함수를 사용하여 손쉽게 계산이 가능합니다. dot()함수는 같은 크기를 같는 두 벡터의 스칼라 곱(dot product)를 반환 합니다.

만약 $v = \left[ 1, 2, 3 \right]$ 이고 $w = \left[4, 5, 6 \right]$ 으로 주어져 있다면

$$ v \cdot w = 1 \times 4 + 2 \times 5 + 3 \times 6 $$

와 같이 계산되어 집니다.

일차원 벡터의 경우, 스칼라 곱은 각각의 원소들을 끼리끼리 곱하여 더한 값을 갖게 됩니다. 위에서 만든 가중치 벡터 $\omega = \left[ \omega_1, \omega_2, \dots \omega_N\right]$ 와 자산별 평균 수익률 벡터 $\mu = \left[ \mu_1, \mu_2, \dots, \mu_N\right]$를 dot()함수를 통해 계산하면 다음과 같이 포트폴리오 전체의 평균 수익률을 얻을 수 있습니다.

$$ \omega \cdot \mu = \omega_1\mu_1 + \omega_2\mu_2 + \dots + \omega_N\mu_N = \mu_P $$
In [34]:
p_returns = np.dot(weights, mean_returns)
print("Expected return of the portfolio: ", p_returns)
Expected return of the portfolio:  1.2995199099515848

다음으로는 포트폴리오의 분산(variance)를 계산해 보도록 하겠습니다. 이를 위해서는 선형대수 관련 지식이 필요합니다. 간단하게 정리해 보도록 하겠습니다.

NaN 값들 다루기

NumPy를 활용하여 실수(Real Number) 값들을 다루다 보면 종종 어레이 상의 nan(Not a Number)값들을 마주치게 됩니다. 이런 nan값들은 데이터가 비어있거나 존재하지 않는 경우를 나타냅니다. 이런 nan값들을 그대로 두고 함수 계산을 적용했다가는 에러가 발생하거나 nan값이 그대로 출력될 가능성이 높습니다. 따라서 에러를 발생시키지 않기 위해 nan값들을 처리하는 방법들을 알아보도록 하겠습니다.

In [35]:
v = np.array([1, 2, np.nan, 4, 5])
print(v)
[ 1.  2. nan  4.  5.]

위와 같이 np.nan값이 들어있는 어레이의 평균을 계산해 보도록 하겠습니다.

In [36]:
print(np.mean(v))
nan

어레이 상에 nan값이 들어 있기 때문에 전체 어레이의 평균이 nan으로 출력됨을 볼 수 있습니다. 특정 어레이에 nan값이 들어 있는지 다음과 같이 np.isnan()함수를 사용하여 확인 할 수 있습니다.

In [37]:
np.isnan(v)
Out[37]:
array([False, False,  True, False, False])

isnan()함수를 사용하면 위와 같이 True/False가 담긴 어레이가 출력됩니다. nan값이 들어 있는 위치에는 True가 아닌곳에는 False값을 저장하여 반환합니다. 그렇다면 어레이에 들어있는 nan값들을 어떻게 제거 할 수 있을까요?

저번 시간에 파이썬 기초 강의에서 ~기호가 어떤 statement의 부정형을 반환한다는 것을 배웠습니다. 따라서 다음과 같이 ~np.isnan()을 사용한다면 nan이 아닌곳에 해당하는 인덱스에 True 값을 얻을 수 있고 이 결과를 어레이의 인덱스에 집어넣어주면 True값에 해당하는 값들만 골라서 선택이 가능합니다.

In [38]:
ix = ~np.isnan(v) # ix에 nan값이 없는 곳에 True 값들을 지정합니다.
print(v[ix]) # 인덱스 어레이가 True 인곳만 선택이 됩니다.
[1. 2. 4. 5.]
In [39]:
print(np.mean(v[ix]))
3.0

NumPy는 위 과정을 하나로 묶어놓은 nanmean()함수를 제공합니다.

In [40]:
print(np.nanmean(v))
3.0

위의 np.nanmean()처럼 NumPy에는 nan값들을 처리하고 특정 연산들을 수행해 주는 함수들이 많이 존재합니다. 자세한 사항은 NumPy documentation을 통해 확인이 가능합니다.

선형대수(Linear Algebra)에 대한 간단한 소개

선형대수는 금융 분야뿐만 아니라 다른 분야에서 널리 사용되고 있습니다. 포트폴리오 이론(Mordern Portfolio Theory)에서 최적의 가중치(Optimal Weight)를 계산할 때 선형대수의 기법이 사용됩니다. NumPy는 이런 선형대수 계산에 필요한 많은 기능들을 제공합니다. 선형대수에 대한 간단한 소개와 NumPy로 어떻게 선형대수 관련 연산이 가능한지 알아보도록 하겠습니다.

먼저 스칼라(Scalar)와 행렬의 곱셈, 합성에 대해 알아보겠습니다. 스칼라는 수학적으로 엄밀히 정의하면 벡터공간(Vector Space)을 이루는 체(Field)의 원소를 의미합니다. 여기서는 단순하게 행렬에 곱해지는 실수(Real Number)라고 생각하시면 됩니다. 스칼라 값을 행렬에 곱함으로써 행렬을 스케일링(scaling)을 할 수 있습니다. 스칼라를 행렬에 곱한다는 것은 행렬 각각의 원소에 스칼라 값을 곱해주는 것 입니다.

행렬은 값들을 모아놓은 형태를 가집니다. $m \times n$ 행렬은 $m$행 $n$열로 구성되어있음을 의미합니다. 어떤 $m \times n$ 행렬에서 $m=n$이라면 정사각 행렬(square matrix)이라 부릅니다. 또 $m=1$이거나 $n=1$인 특수한 경우 해당 행렬을 벡터라고 부릅니다. numpy에는 행렬을 따로 표현하는 matrix객체가 존재하지만 array를 통해서도 행렬의 모든것이 표현 가능하기 때문에 여기서는 array만 사용하기로 합니다. 따라서 이 장에서는 행렬과 어레이를 같은것으로 생각하셔도 됩니다.

일반적으로 행렬이 들어간 방정식을 다음과 같이 표현 가능합니다.

$$ y = A\cdot x $$

$A$는 $m \times n$ 행렬, $y$는 $m \times 1$ 벡터, 그리고 $x$는 $n \times 1$ 벡터를 나타냅니다. 위 식의 우측에서 행렬에 벡터를 곱하는 연산을 사용하였습니다. 어떤 연산을 사용하려면 먼저 그 연산을 정의해야 합니다. 스칼라와 행렬을 곱하는 연산은 위에서 각 원소마다 스칼라를 곱하는 방식으로 정의하였습니다. 그렇다면 행렬과 행렬(또는 벡터)의 곱은 어떻게 정의되는것 일까요?

행렬 곱셈

실수 끼리의 곱셈에서는 $2 \times 3 = 3 \times 2 = 6$ 과 같이 곱셈의 순서를 거꾸로 해도 같은 값이 계산됩니다. 이를 "교환법칙이 성립한다(commutative)"고 합니다. 하지만 행렬의 곱셈에서는 교환법칙이 성립하지 않습니다. 즉 곱해지는 순서에 따라 곱셈 결과가 달라지거나 심지어 정의가 되지 않기도 합니다.

In [41]:
A = np.array([
        [1, 2, 3, 12, 6],
        [4, 5, 6, 15, 20],
        [7, 8, 9, 10, 10]        
    ])
B = np.array([
        [4, 4, 2],
        [2, 3, 1],
        [6, 5, 8],
        [9, 9, 9]
    ])

위에서 생성된 행렬 A는 $3\times 5$, B는 $4 \times 3$의 shape를 가지게 됩니다. 즉 서로 다른 shape를 가지게 되는것입니다. 행렬의 곱셈이 정의되기 위해서는 이 shape가 매우 중요합니다. 곱하기 연산자($\cdot$) 왼쪽 행렬의 열(column)수가 오른쪽 행렬의 행(row)수와 일치해야지만 곱셈이 정의 됩니다. 만약 $m \times n$행렬에 $p \times q$행렬을 우측에 곱한다면 $n=p$인 경우에만 곱셈이 정의 되며 결과로 계산되는 행렬의 사이즈는 다음과 같이 계산되어 집니다.

$$ (m \times n) \cdot (p \times q) = (m \times q) $$
In [42]:
print(np.dot(A,B))
---------------------------------------------------------------------------
ValueError                                Traceback (most recent call last)
<ipython-input-42-88a46c9c9dcc> in <module>()
----> 1 print(np.dot(A,B))

ValueError: shapes (3,5) and (4,3) not aligned: 5 (dim 1) != 4 (dim 0)

위에서 살펴보았듯 A는 $3\times 5$, B는 $4 \times 3$의 shape를 가집니다. 5와 4가 일치하지 않기 때문에 이 곱셈은 정의되지 않습니다. 하지만 A와 B의 순서를 바꾸게 되면 $ (4 \times 3) \cdot (3 \times 5) = (4 \times 5) $ 이기 때문에 곱셈이 정의되며 다음과 같이 4행 5열의 행렬이 계산됩니다.

In [43]:
print(np.dot(B,A))
[[ 34  44  54 128 124]
 [ 21  27  33  79  82]
 [ 82 101 120 227 216]
 [108 135 162 333 324]]

포트폴리오의 분산(Portfolio Variance)

다시 포트폴리오의 예제로 돌아와 봅시다. NumPy의 난수 생성 기능을 통해 랜덤으로 가중치를 두어 기대 수익률을 계산했습니다. 그렇다면 포트폴리오의 분산은 어떻게 계산할까요? 포트폴리오가 각 종목과 해당되는 가중치만큼 자산이 배분되어 구성되었기에 다음과 같이 식을 쓸 수 있습니다.

$$ VAR[P] = VAR[\omega_1 S_1 + \omega_2 S_2 + \cdots + \omega_N S_N] $$

$S_0, \cdots, S_N$는 포트폴리오에 들어있는 자산들을 의미합니다. 만약 모든 자산들이 독립(independent)이라면 자산들 끼리의 공분산(covariance)이 0이 되기 때문에 위 값은 단순하게 다음 식으로 계산이 가능합니다.

$$ VAR[P] = VAR[\omega_1 S_1] + VAR[\omega_2 S_2] + \cdots + VAR[\omega_N S_N] = \omega_1^2\sigma_1^2 + \omega_2^2\sigma_2^2 + \cdots + \omega_N^2\sigma_N^2 $$

하지만 자산들의 데이터를 생성할 때 나머지 9개의 자산이 첫 번째 자산과 가격움직임이 유사하도록(correlated) 설정되었기 때문에 각 자산들은 서로 독립적이지 않습니다. 그렇기 때문에 포트폴리오의 분산을 계산할 때 각 자산들 끼리의 공분산들을 더해주어야 합니다. 포트폴리오 분산에 대한 계산식에 다음과 같이 공분산들의 값들이 추가됩니다.

$$ VAR[P] = \sigma_P^2 = \sum_i \omega_i^2\sigma_i^2 + \sum_i\sum_{i\neq j} \omega_i\omega_j\sigma_i\sigma_j\rho_{i, j}, \ i, j \in \lbrace 1, 2, \cdots, N \rbrace $$

$\rho_{i,j}$는 $S_i$ 와 $S_j$의 상관 계수(correlation)을 의미합니다, $\rho_{i, j} = \frac{COV[S_i, S_j]}{\sigma_i\sigma_j}$. 위 공식은 계산하기 매우 복잡해 보입니다. 만약 반복문을 써서 계산해야 한다면 번거로운 과정이 될것입니다. 하지만 NumPy 어레이와 행렬계산을 통해 위 값을 손쉽게 계산할 수 있습니다. 먼저 모든 종목들의 관련성을 살펴볼 수 있는 공분산 행렬(covariance matrix)을 계산합니다.

In [44]:
cov_mat = np.cov(returns)
print(cov_mat)
[[0.00082889 0.00078472 0.00096386 0.00085766 0.00089229 0.00081201
  0.00081845 0.00080856 0.00080851 0.00078648]
 [0.00078472 0.00116103 0.00093683 0.00077023 0.00081475 0.00076332
  0.00073996 0.0007523  0.00071866 0.00076813]
 [0.00096386 0.00093683 0.00151296 0.00100758 0.00102075 0.00090619
  0.000942   0.00096616 0.00095783 0.00088538]
 [0.00085766 0.00077023 0.00100758 0.0012806  0.00093174 0.00087377
  0.00081194 0.00083563 0.00084589 0.00082239]
 [0.00089229 0.00081475 0.00102075 0.00093174 0.00137311 0.0009276
  0.00081326 0.00088601 0.0008449  0.00096766]
 [0.00081201 0.00076332 0.00090619 0.00087377 0.0009276  0.00117648
  0.00077608 0.00081147 0.00080017 0.00079353]
 [0.00081845 0.00073996 0.000942   0.00081194 0.00081326 0.00077608
  0.00114571 0.00080778 0.00083488 0.00077121]
 [0.00080856 0.0007523  0.00096616 0.00083563 0.00088601 0.00081147
  0.00080778 0.00103371 0.00078191 0.0007657 ]
 [0.00080851 0.00071866 0.00095783 0.00084589 0.0008449  0.00080017
  0.00083488 0.00078191 0.00125131 0.00075222]
 [0.00078648 0.00076813 0.00088538 0.00082239 0.00096766 0.00079353
  0.00077121 0.0007657  0.00075222 0.00109969]]

출력 결과가 포맷팅이 되어있지 않아 약간 지저분하게 느껴질 수 있지만 위의 공분산 행렬은 매우 중요한 정보들을 담고 있습니다. 위에 출력된 공분산 행렬은 아래의 값을 계산하여 출력된 것입니다.

$$ \left[\begin{matrix} VAR[S_1] & COV[S_1, S_2] & \cdots & COV[S_1, S_N] \\ COV[S_2, S_1] & VAR[S_2] & \cdots & COV[S_2, S_N] \\ \vdots & \vdots & \ddots & \vdots \\ COV[S_N, S_1] & COV[S_N, S_2] & \cdots & VAR[S_N] \end{matrix}\right] $$

위 행렬에서 대각 성분에 위치한 값들은 그 순서에 해당하는 자산의 분산을 나타내며 대각 성분이 아닌곳의 값들은 해당 행과 열 번호에 해당하는 자산들의 공분산을 나타냅니다. 여기서 중요한 점은 위의 행렬을 계산을 통해 얻었다면 아래의 식을 통해 포트폴리오 분산을 손쉽게 할 수 있다는 것입니다.

$$ \sigma_p^2 = \omega \ C \ \omega^\intercal $$

$C$는 위에서 계산한 공분산 행렬을 의미하며 $\omega$는 자산들의 가중치를 담은 $1 \times N$ 벡터를 의미합니다. 두 번째 $\omega$뒤의 윗첨자로 달려있는 $\intercal$ 기호는 전치행렬(transpose)를 의미합니다. 포트폴리오의 분산 행렬식에 관한 자세한 정보는 Modern Portfolio Theory 항목을 참고하시면 되겠습니다.

전치행렬(transpose)는 행과 열을 뒤집은 새로운 행렬을 의미합니다. 대각 성분에 대칭 시킨 행렬이라고 생각하시면 이해하기 쉽습니다. 앞에서 정의했던 행렬 A를 통해 이를 살펴보겠습니다.

In [45]:
print(A)
[[ 1  2  3 12  6]
 [ 4  5  6 15 20]
 [ 7  8  9 10 10]]

전치행렬은 np.transpose()를 통해 구할 수 있습니다.

In [46]:
print(np.transpose(A))
[[ 1  4  7]
 [ 2  5  8]
 [ 3  6  9]
 [12 15 10]
 [ 6 20 10]]

A는 $3 \times 5$의 어레이 이며 이를 transpose한 $A^\intercal$은 $3 \times 5$의 어레이가 됩니다. 포트폴리오 분산 행렬식에서 $\omega$는 $1 \times n$ 모양을 가진 행벡터(row vector)이 입니다. 따라서 이를 transpose하면 $n \times 1$ 모양을 갖는 열벡터(column vector)가 됩니다.

이야기가 길어졌는데 위의 행렬식에서 포트폴리오 분산의 shape를 계산해 보도록 하겠습니다. 위의 크기 정보들을 종합하면 다음과 같이 하나의 스칼라값($1 \times 1$의 모양)이 얻어짐을 알 수 있습니다.

$$ \text{Dimensions}(\sigma_p^2) = \text{Dimensions}(\omega C \omega^\intercal) = (1 \times N)\cdot (N \times N)\cdot (N \times 1) = (1 \times 1)$$

정리하자면, 공분산 행렬의 왼쪽에 가중치 행벡터를 곱하고 오른쪽에는 가중치 열벡터를 곱함으로써 포트폴리오 분산에 해당하는 스칼라 값을 얻을 수 있습니다. NumPy 함수들을 사용하여 직접 계산해 보도록 하겠습니다.

In [47]:
# 포트폴리오 분산을 계산합니다.
var_p = np.dot(np.dot(weights, cov_mat), weights.T)
vol_p = np.sqrt(var_p)
print("Portfolio volatility: ", vol_p)
Portfolio volatility:  0.02989249721024448

위 값을 검증하기 위해 아래와 같이 다른 방법을 통해 포트폴리오의 변동성(volatility)을 구해보도록 하겠습니다.

In [48]:
# Confirming calculation
vol_p_alt = np.sqrt(np.var(np.dot(weights, returns), ddof=1))
print("Portfolio volatility: ", vol_p_alt)
Portfolio volatility:  0.029892497210244475

위의 ddof 인자는 Delta Degrees of Freedom 으로 정수값을 가집니다. 이는 조금더 통계적 지식이 필요한 내용이라서 추후 통계 부분을 다루는 콘텐츠에서 설명하도록 하겠습니다. 어찌되었든 위의 행렬로 계산한 값과 일치하는 값이 나옴을 확인할 수 있습니다.

NumPy의 함수들을 사용해서 포트폴리오 이론과 그 계산과정을 코드로 옮겨보았습니다. 선형대수 지식과 행렬을 다루는 기술은 금융분야 뿐만 아니라 데이터를 다루는 다른 분야에서도 필수적이기 때문에 위 과정들을 하나씩 잘 소화하시고 행렬을 코드로 구현하는 법을 숙지해 두신다면 큰 도움이 될것입니다. 읽어주셔서 감사합니다. 다음 콘텐츠에서는 NumPy와 함께 파이썬에서 가장 많이 쓰이는 Pandas 모듈에 대해 다뤄보도록 하겠습니다.

</html>