はじめに
かねてからYahoo!ファイナンスなどの情報に、株価について期間を指定したリターンや標準偏差や回帰分析のあてはまりなどが無いので、そのようなデータを簡単に見れるようにしたいという願望がありました。加えて、できるだけ多くの銘柄の情報をみたいという願望もありました。Rを使えばその全てが簡単に実現できるのでチャレンジしてみました。
・データ収集
・データテーブルの整形
・成長率の計算
・標準偏差の計算
・自己回帰の計算
と非常にシンプルな工程ですが、3600銘柄分のこれらの情報を一気に手に入れれるというのは魅力的だなぁと思います。
データの収集
RFinanceYJというパッケージをインストールして、以下のrfinance_module.Rのスクリプトを置いておきます。
|
1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60 61 62 63 64 65 66 67 68 69 70 71 72 73 74 75 76 77 78 79 80 81 82 83 84 85 86 87 88 89 90 91 92 93 94 95 96 97 98 99 100 101 102 103 104 105 106 107 108 109 110 111 112 113 114 115 116 117 118 119 120 121 122 123 124 125 |
library(RFinanceYJ) #API quoteStockTsData <- function(x, since=NULL,start.num=0,date.end=NULL,time.interval='daily') { time.interval <- substr(time.interval,1,1) function.stock <- function(quote.table.item){ if( xmlSize(quote.table.item) < 5) return(NULL) d <- convertToDate(xmlValue(quote.table.item[[1]]),time.interval) o <- as.number(xmlValue(quote.table.item[[2]])) h <- as.number(xmlValue(quote.table.item[[3]])) l <- as.number(xmlValue(quote.table.item[[4]])) c <- as.number(xmlValue(quote.table.item[[5]])) v <- ifelse(xmlSize(quote.table.item) >= 6,as.number(xmlValue(quote.table.item[[6]])),0) a <- ifelse(xmlSize(quote.table.item) >= 7,as.number(xmlValue(quote.table.item[[7]])),0) return(data.frame(date=d,open=o,high=h,low=l,close=c,volume=v, adj_close=a)) } return(quoteTsData(x,function.stock,since,start.num,date.end,time.interval,type="stock")) } quoteFundTsData <- function(x, since=NULL,start.num=0,date.end=NULL,time.interval='daily') { time.interval <- substr(time.interval,1,1) function.fund <- function(quote.table.item){ d <- convertToDate(xmlValue(quote.table.item[[1]]),time.interval) if(time.interval=='monthly'){ d <- endOfMonth(d) } c <- as.number(xmlValue(quote.table.item[[2]])) v <- as.number(xmlValue(quote.table.item[[3]])) return(data.frame(date=d,constant.value=c,NAV=v)) } return(quoteTsData(x,function.fund,since,start.num,date.end,time.interval,type="fund")) } quoteFXTsData <- function(x, since=NULL,start.num=0,date.end=NULL,time.interval='daily') { time.interval <- substr(time.interval,1,1) function.fx <- function(quote.table.item){ d <- convertToDate(xmlValue(quote.table.item[[1]]),time.interval) o <- as.number(xmlValue(quote.table.item[[2]])) h <- as.number(xmlValue(quote.table.item[[3]])) l <- as.number(xmlValue(quote.table.item[[4]])) c <- as.number(xmlValue(quote.table.item[[5]])) return(data.frame(date=d,open=o,high=h,low=l,close=c)) } return(quoteTsData(x,function.fx,since,start.num,date.end,time.interval,type="fx")) } ###### private functions ##### #get time series data from Yahoo! Finance. quoteTsData <- function(x,function.financialproduct,since,start.num,date.end,time.interval,type="stock"){ r <- NULL result.num <- 51 financial.data <- data.frame(NULL) start <- (gsub("([0-9]{4,4})-([0-9]{2,2})-([0-9]{2,2})","&sy=\\1&sm=\\2&sd=\\3",since)) end <- (gsub("([0-9]{4,4})-([0-9]{2,2})-([0-9]{2,2})","&ey=\\1&em=\\2&ed=\\3",date.end)) if(!any(time.interval==c('d','w','m'))) stop("Invalid time.interval value") extractQuoteTable <- function(r,type){ if(type %in% c("fund","fx")){ tbl <- r[[2]][[2]][[7]][[3]][[3]][[9]][[2]] } else{ tbl <- r[[2]][[2]][[7]][[3]][[3]][[10]][[2]] } return(tbl) } #while( result.num >= 51 ){ while(1){ start.num <- start.num + 1 quote.table <- NULL quote.url <- paste('http://info.finance.yahoo.co.jp/history/?code=',x,start,end,'&p=',start.num,'&tm=',substr(time.interval,1,1),sep="") try(r<-htmlParse(quote.url)) if( is.null(r) ) stop(paste("Can not access :", quote.url)) try( quote.table <- xpathApply(r,"//table")[[2]], TRUE ) quote.size<-xmlSize(quote.table) #cat(paste("size:",quote.size)) if(xmlSize(quote.table)<=1){ return (financial.data) } if( is.null(quote.table) ){ if( is.null(financial.data) ){ stop(paste("Can not quote :", x)) }else{ financial.data <- financial.data[order(financial.data$date),] return(financial.data) } } size <- xmlSize(quote.table) for(i in 2:size){ financial.data <- rbind(financial.data,function.financialproduct(quote.table[[i]])) } #result.num <- xmlSize(quote.table) Sys.sleep(1) } financial.data <- financial.data[order(financial.data$date),] return(financial.data) } #convert string formart date to POSIXct object convertToDate <- function(date.string,time.interval) { #data format is different between monthly and dialy or weekly if(any(time.interval==c('d','w'))){ result <- gsub("^([0-9]{4})([^0-9]+)([0-9]{1,2})([^0-9]+)([0-9]{1,2})([^0-9]+)","\\1-\\3-\\5",date.string) }else if(time.interval=='m'){ result <- gsub("^([0-9]{4})([^0-9]+)([0-9]{1,2})([^0-9]+)","\\1-\\3-01",date.string) } return(as.POSIXct(result)) } #convert string to number. as.number <- function(string) { return(as.double(as.character(gsub("[^0-9.]", "",string)))) } #return end of month date. endOfMonth <- function(date.obj) { startOfMonth <- as.Date(format(date.obj,"%Y%m01"),"%Y%m%d") startOfNextMonth <- as.Date(format(startOfMonth+31,"%Y%m01"),"%Y%m%d") return(startOfNextMonth-1) }</pre> |
以下のコードで2016-01-01から現在に至るまでの、指定した株価コードの株価情報を集めてきます。3600銘柄の株価コードを用意する必要があるので、集めてcodelist.csvで保存しておきます。銘柄コードの情報はここでは掲載しません。（適当なサイトにころがっています。）
|
1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 |
source("rfinance_module.R") codelist <- read.csv("codelist.csv") since <- "2016-01-01" stock_price_database <- 0 for(i in 1:nrow(codelist)){ if(i ==1){ dataset <- quoteStockTsData(paste(codelist[i,1],".t",sep = ""),since=since) code <- matrix(codelist[i,1], nrow(dataset), 1) dataset <- cbind(code, dataset) stock_price_database <- dataset } else{ dataset <- quoteStockTsData(paste(codelist[i,1],".t",sep = ""),since=since) code <- matrix(codelist[i,1], nrow(dataset), 1) dataset <- cbind(code, dataset) stock_price_database <- rbind(stock_price_database,dataset) } } |
以下のコードを実行して、元のデータを調整済み終値について、縦が日付で横が銘柄のテーブルに整形します。
|
1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 |
#make stock price data table library(dplyr) code_list <- data.frame(code=distinct(stock_price_database,code)$code) date_list <- data.frame(date=distinct(stock_price_database,date)$date) stock_dataset <- date_list stock_price_database_adj <- stock_price_database[,which(colnames(stock_price_database) %in% c("code","date","adj_close"))] for( i in 1:nrow(code_list)){ #nrow(code_list) each_stock <- filter(stock_price_database_adj, code == code_list$code[i]) each_stock <- each_stock[,-1] each_stock <- distinct(each_stock, date,.keep_all = TRUE) colnames(each_stock)[2] <- paste("sdj",code_list$code[i],sep = "") stock_dataset <- left_join(stock_dataset, each_stock, by = "date") } |
成長率
指定した期間での成長率を計算します。
|
1 2 3 4 5 6 7 8 9 10 11 12 13 14 |
#calculate growth rate #start_date start_date <- "2016-01-04" #end_date end_date <- "2017-02-17" #return rate datatable return_data <- filter(stock_dataset, date==start_date |date==end_date) return_data <- return_data[,-1] return_data <- rbind(return_data, (as.numeric(return_data[1,]) / as.numeric(return_data[2,])-1)*100) return_data <- data.frame(code=colnames(return_data),t(return_data)) colnames(return_data) <- c("code",end_date,start_date,"return_rate") return_data <- return_data[order(return_data$return_rate,decreasing = TRUE),] |
標準偏差
ここではシンプルに標準偏差の計算を行います。
|
1 2 3 4 5 6 7 8 9 10 11 |
#volatility volatility_data <- stock_dataset[,-1] volatility_data <- log(volatility_data[-nrow(volatility_data),]) - log(volatility_data[-1,]) volatility <- data.frame(matrix(0,1,2)) colnames(volatility) <- c("code", "stdev") for(i in 1:length(volatility_data)){ volatility[i,1] <- colnames(volatility_data)[i] volatility[i,2] <- sd(volatility_data[,i],na.rm = TRUE) } |
自己回帰による決定係数
自己回帰でフィッティングの良い銘柄はどれなのかを単純に知りたいので、3600銘柄分計算させます。
|
1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 |
#AR(1) estimation ar1_estimation <- stock_dataset[,-1] ar1_rsuares <- data.frame(matrix(0,1,2)) colnames(ar1_rsuares) <- c("code", "adj_rsquare") for(i in 1:length(ar1_estimation)){ #objective variable ar1_y <- ar1_estimation[-nrow(ar1_estimation),i] #input variable ar1_x <- ar1_estimation[-1,i] #combine data ar1_data <- data.frame(cbind(ar1_y, ar1_x)) #estimation if(class(try(lm(data = ar1_data, formula = ar1_y ~. )))=="lm"){ model <- lm(data = ar1_data, formula = ar1_y ~. ) model_sum <- summary(model) ar1_rsuares[i,1] <- colnames(ar1_estimation)[i] ar1_rsuares[i,2] <- model_sum$adj.r.squared } else{ ar1_rsuares[i,1] <- colnames(ar1_estimation)[i] ar1_rsuares[i,2] <- "NULL" } } |
以上より、見てみたかった指標の3600銘柄分を一気に計算することができました。
可視化
試しに一番成長率の高かった株価をプロットしてみます。
|
1 2 3 4 5 6 7 8 9 10 11 12 13 |
#visualization library(ggplot2) stock_code <- "sdj7612" df <- data.frame(stock_dataset[,which(colnames(stock_dataset) %in% c("date", stock_code))]) g <- ggplot(df,aes(x= date, y= df[,2]))+ geom_line(colour = "magenta",linetype = 2,size = 0.5)+ geom_smooth (method = "lm")+ xlab("Date")+ ylab("Stock Price")+ ggtitle("Stock Price") plot(g) |
おわりに
株価データにカジュアルにアクセスできる環境を手に入れることができたので、今後は様々な指標や時系列解析・業界ごとの分析を試してみたいと思います。