Post List

2016년 7월 9일 토요일

How to backtest? : Low volatility example




install.packages("zoo")
install.packages("PerformanceAnalytics")
install.packages("quantmod")

library(zoo) 
library(PerformanceAnalytics)
library(quantmod)

먼저, install.packages("~~~") 패키지를 받아 설치하는 명령어 입니다.
PC 설치해 패키지는 계속 저장되어 있습니다.

다음으로 library(~~~) 명령어를 통해 설치한 패키지를 불러옵니다.
요건 매번 실행해주어야 하므로, 자주 쓰는 패키지들은
메모장 같은곳에 정리해 두었다가, 모든 코딩 최상단에 박아두는 것이 편합니다.



price_m = read.csv("price_m.csv", row.names=1, header=TRUE)
price_m[is.na(price_m)] = 0
price_kospi = read.csv("kospi.csv", row.names = 1, header = TRUE)
ro = dim(price_m)[1]
co = dim(price_m)[2] 

read.csv 명령어는 csv 파일을 읽어오는 명령어 입니다.
(
이상하게 R 에서는 xls 파일 땡겨오면 뻑이 잘납니다.)
뒤의 row.names n 번째 행까지를 이름으로 설정할 것인지,
header
TRUE 헤드 테이블을 만듭니다.

처음에 csv 파일에는 다음과 같은 형태의 raw data 있습니다.

A00000
A00001
A00002
1994-12-31
57261
26622
24669
1995-01-31
57369
25050
22004
1995-02-28
61409
25148
23429
1995-03-31
74429
26425
26156
1995-04-30
77853
27309
28822
1995-05-31
71393
27899
27954
1995-06-30
78823
27899
29007
1995-07-31
90800
27702
26776

위의 read.csv 명령어를 사용하면, 날짜들은 행이름에 입력됩니다.
종목명들은 헤드 테이블에 입력됩니다.

is.na
문장으로 쓰면if XX is "NA" data 입니다.
주가가 없는 데이터의 경우 NA 저장되므로,
위의 명령어를 처주면 price_m 매트릭스에서
NA
데이터는 모두 0 으로 바꿔줍니다.

dim(xx) 
행과 열의 갯수를 반환합니다.
dim(xx)[1]
갯수이며, 이는 nrow(xx) 결과가 같습니다.
dim(xx)[2]
갯수이며, 이는 ncol(xx) 결과가 같습니다.


ret_kospi = Delt(price_kospi)
ret_1m = matrix(0,ro,co)

for (i in 1:dim(ret_1m)[2]) {
  ret_1m[, i] = Delt((price_m)[, i])
}  

Delt(XX) quantmod 패키지의 명령어로써, % 변화값을 반환합니다.
raw data
주가 이므로, 여기서는 수익률이 반환됩니다.
(
패키지를 설치한 , ?Delt 입력하면, 해당 명령어에 대한 자세한 설명이 나옵니다)

matrix(k, x1, x2) 
k 구성된 x1 , x2 매트릭스를 만들어 줍니다.
ret_1m
이라는 수익률이 저장될 zero 매트릭스를 미리 설정해 줍니다.

Delt
명령어의 경우, 1 밖에 실행이 되지 않으므로
for
구문을 통해, 열을 하나씩 땡깁니다.

for (i in 1:dim(ret_1m)[2]) 
1부터 ret_1m 매트릭스의 열갯수까지
명령어를 실행하겠다는 겁니다.
price_m[, i]
모든 행의 i 이라는 의미 입니다.
price_m 1, 2, 3, ... , n 까지 Delt 함수를 계속 적용해줘서
ret_m
1, 2, 3,,,, n 열에 수익률을 저장시켜 줍니다.

Matlab
R 차이점은, Matlab 경우 매트릭스 결과값이 저장될 매트릭스를 미리 만들어 주지 않아도 되지만, R 미리 만들어 주어야 된다는 점입니다.



std = rbind( matrix(0,59,dim(ret_1m)[2]),rollapply(ret_1m,60,sd) ) 
std[std == 0] = NA

rollapply 함수는 zoo 패키지의 함수로써, 롤링 값을 계산하는 명령어 입니다.
rollapply(x,n,f)
x 매트릭스에서, n구간의, sd 구하겠다는 뜻입니다.
여기서는, 월간 수익률의 60개월 변동성을 구합니다.

matrix(0,59,dim(ret_1m)[2] 
0행렬을 만들어 주는 이유는,
60
개월 롤링 변동성을 구해주면, 1~59행은 계산되지 않기 때문에.
수익률 매트릭스의 갯수보다 59개가 부족해 지게 됩니다.
나중에 간의 Time Period 같도록 만들어 주기 위해,
부분은 0으로 잡아 줍니다.

rbind 
여러가지 매트릭스들을 행으로 쌓아주는 역할을 합니다.
59
행의 zero matrix, 그리고 계산해준 롤링 변동성을 합쳐주면
price
행렬이 같은 std 매트릭스가 만들어 집니다.

std==0 
std 0 값입니다.
60
개월 변동성이 0이라는 말은 주가가 0 이라는 뜻입니다.
이는 dummy data 이므로 NA 바꿔줍니다.


ret_lowvol_quan = matrix(0,ro,5)
ret_lowvol = matrix(0,ro,1)

rownames(ret_lowvol_quan) = rownames(price_m)
rownames(ret_lowvol) = rownames(price_m)

먼저 ret_lowvol_quan ret_lowvol 이라는 zero matrix 만들어 줍니다.
여기는 나중에 계산된 분위수별 수익률, 상위종목 수익률이 저장됩니다.

rownames 
행이르메 관련된 함수입니다.
각각 만들어준 매트릭스의 행이름을
price_m
행이름, 날짜로 저장해주겠다는 것입니다.
(Date
행이름을 정해주지 않으면 PerformanceAnalytics 명령어가 실행되지 않습니다.)



for (i in (which(rownames(price_m) == "1999-12-31")) : (ro-1)) {

con = matrix(0,2,co)
con[1,] = std[i, ]

con[, which(price_m[i,] == 0)] = NA

qt = matrix(0,1,5)
for (q in 1 : 5) {
  qt[, q] = quantile(con[1,], (q-1)*2/10, na.rm=TRUE) 
  
for (w in 1 : co) {
  if ( ( con[1, w] > qt[, q] ) & (!is.na(con[1,w])) ) { 
    con[2, w] = q }
  }
}

for (q in 1 : 5) { 
  ret_lowvol_quan[i+1,q] = mean(ret_1m[i+1, which ((con[2, ] == q) ) ])
}

  ret_lowvol[i+1] = mean(ret_1m[i+1, which(rank(con[1, ]) <= 20 ) ])
}

사실 위는 전부터 input 데이터를 만드는 단계였고,
여기부터가 백테스트 단계 입니다.

가장 먼저 which 괄호 안의 데이터가 몇번째에 있냐는 숫자를 반환하는 명령어 입니다.
여기서는 price_m 행이름이 "1999-12-31" 곳이 몇번째 인가를 찾습니다.
데이터 기준으로는 61 번째 부터 입니다.
백테스트 기간의 시작시점을 셋팅할 날짜를 바꿔주면 됩니다.

다음, for 구문 안에 2 X 열갯수 만큼의 zero 행렬을 만들어 줍니다.
if
안에 있으므로, i = 61 번째 계산이 다되고,
i = 62
번째 계산이 시작되면, 매번 zero 행렬이 새로 만들어 집니다.

con[1,] = std[i,] 
std i번째 데이터 전부를 con 첫번째 행에 가져오는 겁니다.
con[, which(price_m[i,] == 0)] = NA 
price_m i번째 행에서 데이터가 0 값들의 위치를 찾고, con (현재는 std 값들이 저장되 있음) 해당 데이터들을 모두 NA 만들어 주는 겁니다

이제 데이터 들에서 Low Vol 뽑아냅니다.

먼저 1 X 5 짜리 qt 행렬을 만들어 줍니다.
quantile
분위수를 뽑는 명령어 입니다.
(q-1)2*/10
해주면, std 기준 각각 0%, 20%, 40%, 60%, 80% 분위 (1분위, 2분위, 3분위, 4분위, 5분위 criteria ) qt 행렬에 저장됩니다.

na.rm = TRUE 
경우, NA 데이터가 있으면, 이를 무시하고 계산하겠 다는 뜻입니다. 컨디션을 입력해주지 않을 경우, NA 데이터들로 인해 값이 계산되지 않습니다. (이는 NA 있는 데이터들을 처리해줄 굉장히 유용한 컨디션 입니다.)

다음 줄로 내려와서 !is.na  NA 데이터가 아닐 경우 라는 뜻입니다.
알에서 ! if not  의미입니다.

, 모든 열에서 con[1,w] (여기서는 std) qt q열의 분위수 보다 크고, con[1,w] NA 아닐 경우 con 2번째 행에 q (분위수) 입력합니다.

기본적으로 변동성이 가장 작은 종목들의 경우, 0% (1분위) 기준 보다는 크고 20% (2분위) 기준 보다는 작을 것이므로 1 이라는 값이 저장됩니다. 변동성이 가장 종목은 반대로 80% (5분위) 기준 보다 크므로, 5 라는 값이 저장됩니다.

다음으로는 con 2번째 행에서 q값에 해당되는 열들을 찾아, ret_1m i+1 (다음달 수익률) 값들의 mean 구합니다.

이를 ret_lowvol_quan 함수의 i+1번째 , q번째 열에 저장합니다.
, 분위수 별로 다음달 수익률의 평균값을 계산하는 겁니다.

마지막 줄은, con 첫번째 행에서, 순위가 20 이하인 종목들의 열들을 찾아, ret_1m i+1 값들의 mean 구합니다.

사실, 분위수별로 분석할 것이 아니라 단순히 상위 종목을 찾을 거면 이거 한줄만 있어도 상관없습니다.

참고로 R에서 랭크 함수는 기본적으로 오름차순 기준입니다. Low vol 경우 변동성이 낮은 종목을 찾는거라 오름차순이 맞지만, momentum 처럼 내림차순으로 찾아야 경우 rank(1/con[1,]) 처럼 역수의 형태로
해주거나, rank(-con[1,]) 처럼 음수를 붙여주면 됩니다.

ret_lowvol_quan 매트릭스에는 분위수별로 저장된 수익률 값이,
ret_lowvol 에는 상위 20 종목의 수익률 값이 저장되게 됩니다.



ret_lowvol_quan = ret_lowvol_quan[ which(rownames(ret_lowvol_quan) == "2000-01-31") : ro, ]
ret_lowvol = as.matrix(ret_lowvol[ which(rownames(ret_lowvol) == "2000-01-31") : ro, ])
ret_kospi = as.matrix(ret_kospi[ which(rownames(ret_kospi) == "2000-01-31") : ro, ])

우리는 2000 1 부터 수익률을 보고 싶으므로, 앞의 데이터들은 날려 버립니다. which 함수를 통해 행이름이 "2000-01-31" 곳을 찾고, 곳부터 끝까지를 다시 저장하는 겁니다. (사실 안해도 상관 없습니다.)



chart.CumReturns(ret_lowvol_quan)   
legend('topleft', c("1st","2nd","3rd","4th","5th"),col=1:5, lty=1, horiz = TRUE, cex = 0.5, bty="n")

chart.CumReturns(ret_lowvol)   

legend('topleft',c("Lowvol Index"),col=1:2, lty=1, horiz = TRUE, bty="n")

chart.Drawdown(ret_lowvol_quan)

legend('bottom', c("1st","2nd","3rd","4th","5th"),col=1:5, lty=1, horiz = TRUE, cex = 0.5, bty="n")

table.AnnualizedReturns(ret_lowvol_quan)

chart.CumReturns PerformanceAnalytics 패키지가 제공하는 굉장히 강력한 명령어 입니다. 누적 수익률 그래프를 그려줍니다.

legend
범례를 만들어 줍니다. topleft 범례가 위치할 , c("1st","2nd","3rd","4th","5th") 범례의 내용, col=1:5 범례의 색깔, lty 범례의 모양, horiz 행으로 나타낼 것인지 열로 나타낼 것인지, cex 범례의 크기, bty 투명하게 것인지의 여부입니다.

분위 수별 그래프와, 상위 20종목 그래프를 각각 그립니다.

chart.Drawdown 
역시 PerformanceAnalytics  제공하는 함수로써, Drawdown 그래프를 그려줍니다.

table.AnnualizedReturns
  PerformanceAnalytics  제공하는 함수로써, 연율화 수익률, 변동성, 샤프지수 (Rf=0% 가정) 계산해 줍니다.


앞에 써놓은 코드들을 펼치면 아래와 같습니다.
라인 by 라인으로 해보시면 쉽게 느껴질 겁니다.



install.packages("zoo")
install.packages("PerformanceAnalyrics")
install.packages("quantmod")

library(zoo) 
library(PerformanceAnalytics)
library(quantmod)

price_m = read.csv("price_m.csv", row.names=1, header=TRUE)
price_m[is.na(price_m)] = 0
price_kospi = read.csv("kospi.csv", row.names = 1, header = TRUE)

ro = dim(price_m)[1]
co = dim(price_m)[2] 

ret_kospi = Delt(price_kospi)
ret_1m = matrix(0,ro,co)

for (i in 1:dim(ret_1m)[2]) {
  ret_1m[, i] = Delt((price_m)[, i])
}  

std = rbind( matrix(0,59,dim(ret_1m)[2]),rollapply(ret_1m,60,sd) ) 
std[std == 0] = NA

ret_lowvol_quan = matrix(0,ro,5)
ret_lowvol = matrix(0,ro,1)

rownames(ret_lowvol_quan) = rownames(price_m)
rownames(ret_lowvol) = rownames(price_m)

for (i in (which(rownames(price_m) == "1999-12-31")) : (ro-1)) {

con = matrix(0,2,co)
con[1,] = std[i, ]

con[, which(price_m[i,] == 0)] = NA

qt = matrix(0,1,5)
for (q in 1 : 5) {
  qt[, q] = quantile(con[1,], (q-1)*2/10, na.rm=TRUE) 
  
for (w in 1 : co) {
  if ( ( con[1, w] > qt[, q] ) & (!is.na(con[1,w])) ) { 
    con[2, w] = q }
  }
}

for (q in 1 : 5) { 
  ret_lowvol_quan[i+1,q] = mean(ret_1m[i+1, which ((con[2, ] == q) ) ])
}

  ret_lowvol[i+1] = mean(ret_1m[i+1, which(rank(con[1, ]) <= 20 ) ])
}

ret_lowvol_quan = ret_lowvol_quan[ which(rownames(ret_lowvol_quan) == "2000-01-31") : ro, ]
ret_lowvol = as.matrix(ret_lowvol[ which(rownames(ret_lowvol) == "2000-01-31") : ro, ])
ret_kospi = as.matrix(ret_kospi[ which(rownames(ret_kospi) == "2000-01-31") : ro, ])

chart.CumReturns(ret_lowvol_quan)   
legend('topleft', c("1st","2nd","3rd","4th","5th"),col=1:5, lty=1, horiz = TRUE, cex = 0.5, bty="n")

chart.CumReturns(ret_lowvol)   
legend('topleft',c("Lowvol Index"),col=1:2, lty=1, horiz = TRUE, bty="n")

chart.Drawdown(ret_lowvol_quan)
legend('bottom', c("1st","2nd","3rd","4th","5th"),col=1:5, lty=1, horiz = TRUE, cex = 0.5, bty="n")


table.AnnualizedReturns(ret_lowvol_quan)

댓글 2개:

  1. 위 데이터는 어디에서 구할 수 있나요?

    답글삭제
  2. 주식 수익률을 구하는 함수 Delt 명령어가 price_kospi 데이터프레임을 인수로 받아서 오류가 발생하는 것 같습니다.
    인터넷에 찾아보니 Delt 명령어는 (m X 1), 즉 벡터를 인수로 받는다고 나와있는데
    데이터 프레임 price_kospi에서 종목별로 수익률을 구하려면 어떻게 해야 할까요?

    답글삭제