주식( tidyquant, tqk), 시계열 (forecast) 패키지 통계분석 예제

  가장 얻기는 쉽지만,  활용하기 어려운 데이터 셋 주식 데이터 이다.   데이터의 시계열의 흐름을 보고 미래 사항을 예측 하지만, 아무리 큰 시장에도 심리적인 사항이 있어 정확하게 예측 하기는 어렵다.  하지만,  포토 폴리오를 예측을 할때,  해당 종목이 계절성을 띄고 있는지 선형성을 띄고 있는지 분석 하여, 전략을 수립 할 수 있다. 

우리는 아래 2가지 파생상품을 가지고 이야기 할 것이다.

• 미국 나스닥(코드 : 133690)
• 골드 선물 (코드 : 319640)

둘다 어떻게 하겠다고,  예측은 할 수 없지만,  기다려서 버는 시점을 알 수 있다. 

지금은  예제를 사용하였기 때문에,  코드 자체가 간소화 되어 있다.  상세한 것은 별도의 코너를 만들겠다.  

라이브러리를 불러온다.  아래 코드에서 패키지를 로딩 할때,  시계열 예측(forecast), 데이터를 핸들링 할( dplyr),  그래픽을 시각화 할 ggplot2이 없다. 이유는    위의 3개의 기능을 하는 패키지가 이미 tidyquant에 내장 되어 있기 때문이다. 

## 불필요 객체 지우기  
  rm(list = ls())

# library loading 
  library(tidyquant)  
  library(tqk) 

 일반적으로 반복되는 연산을 할 때 함수를 사용한다.  아래는 주식을 가져오는 함수 인데,  인자(argument)를  시작일,  종목코드 2개를 받아서 계산을 한다. 여기에서 데이터를 가져올때는 tqk에 내장된 함수와 tidyquant에 내장된 dplyr 함수를 이용한다.   시계열 데이터는 일반적으로 월단위로 집계한다. 

# 주식 가져오는 함수 
  data_get_fun <- function(days, stock_code){  

    # Data get from tqk 
      stock <- tqk_get(stock_code, from= days)  
      stock$종목코드 <- stock_code
  
    # R은 데이터를 가져와 전저리 할때 단순함      
      stock <-   stock %>% 
        arrange(date) %>% 
        mutate(전일비 = c(NA,  diff(adjusted))) %>% 
        filter(!is.na(전일비)) %>% 
        mutate(전일시세 = adjusted - 전일비) %>% 
        mutate(등락율 =  round((adjusted - 전일시세)/전일시세 * 100, 2))
      
    return(stock)
  }
  
# 133690(나스닥),  319640 (골드선물)
# 오늘 보다 2019년 4월 부터 데이터를 가져온다. 
  days = as.Date("2019-04-01")
  
  nasdaq <- data_get_fun(days, "133690")
  gold <- data_get_fun(days, "319640")  

# 데이터 프레임으로 취합 한다.  
  df <- rbind(nasdaq, gold)

# 월 단위로 묶는다.  
df_sum <-    df %>% 
         mutate(월 = paste0(substr(date, 1,4), substr(date, 6,7))) %>%           group_by(종목코드, 월) %>% 
          summarise(가격 = mean(close))

  

  시계열 데이터를 만들려면,  월별로 집계된 값을 전처리 하여, 벡터로 만들어야 한다.  "$" 서브세팅은 데이터 프레임에 있는 특정 필드 값을 벡터화하여 가져오는 것이다. 

  미국 나스닥을 선형 추세 기법으로 하였을 때 RSME(평균 제곱근 오차)가 가장 Fitting화 되기 때문에 더블지수 평할로 분석 하고 그린 것이다.  여기에서 h=5 라는 파라미터가 있는데, 앞으로 5개월을 예측 한다는 것이다.  그리고 더블지수 감쇠 라는 파라미터가 있는데 damped = TRUE, phi = 0.9 로 분석 한 것이다.  즉 예측량을 0.9로 10% 감쇠한다는 의미이다.

# 미국 나스닥    
  x <-ts( (df_sum %>% filter(종목코드 == "133690"))$가격 ,
          start = c(2019, 04), 
          end = c(substr(today(),1,4),substr(today(),6,7)),
           frequency= 12)
  
  # 더블지수 평활
  fd <- holt(x, h=5)
  fd1 <- holt(x, damped = TRUE, phi = 0.9, h=5)

  # 더블지수 평활  
  autoplot(x) +
    autolayer(fd, series="홀트 기법", PI=FALSE) +
    autolayer(fd1, series="감쇠 홀트 기법", PI=FALSE) +
    ggtitle("더블지수 평활로  얻은 예측값") + xlab("연도") +
    ylab("미국나스닥 가격") +
    guides(colour=guide_legend(title="예측값"))

더블지수 평활

 

   골드선물의 케이스는 계절성 패턴을 가지고 있어,  가장 낮은 시점이나, 떨어지려고 하는 시점에 구매 해서  올라 갈 때 까지 계속 가지고 있다고,  어느 정도 올랐다고 하면, 파는 전략으로 할 것이기에  매수 지점을 예측 하고, 매도 지점을 예측 하는 계절적인 패턴을 가지고 있어,  투기적인 성향을 띄고 있어 보인다.   여기에서 덧셈기법은 최근의 실적을 더 많이 반영 하고,  곱셈기법은 과거의 실적을 균등하게 반영한다. 

 # 골드 선물
  x <-ts( (df_sum %>% filter(종목코드 == "319640"))$가격 ,
          start = c(2019, 04), 
          end = c(substr(today(),1,4),substr(today(),6,7)),
          frequency= 12)

  
  # 홀트 -윈터스의 덧셈 
  fit1 <- hw(x, seasonal = "additive" , h =5)
  fit2 <- hw(x,seasonal="multiplicative", h =5)
  
  autoplot(x) +
    autolayer(fit1, series="HW 덧셈 예측", PI=FALSE) +
    autolayer(fit2, series="HW 곱셈 예측",
              PI=FALSE) +
    ggtitle("홀트 윈터스로  얻은 예측값") + xlab("연도") +
    ylab("골드 선물 가격") +
    guides(colour=guide_legend(title="예측"))

홀트 윈터스 평활