2013年1月25日金曜日

Rで225とSP500の週次での季節性を見る

日経225とSP500の週次での季節性を見る

最近、季節性の話を友達としてたので、week numbers別のreturnを色々とgraph化してみた。 対象は、日経225とSP500とする。データは、quantmodパッケージでUSのyahooから取得する。 以下のような関数を書いておいて週次ごとのreturnにしておく

library(quantmod)
## 週次return作成関数
get.nweek.return <- function(sym, syear = 2001, eyear = 2012) {
    ## 前年の終わりをstart日付にする
    fday <- as.Date(paste(syear, "-01-01", sep = "")) - 3
    org.data <- getSymbols(sym, from = fday, auto.assign = F)
    data <- adjustOHLC(org.data, use.Adjusted = T)[, 4]
    years <- format(range(index(data)), "%Y")
    d0 <- to.weekly(data)[, 1]
    d0 <- ROC(d0)[-1]
    dt <- do.call(rbind, lapply(syear:eyear, function(year) {
        tmp <- d0[as.character(year)]
        ret <- coredata(tmp)
        data.frame(year = year, nweek = 1:nrow(tmp), ret = ret)
    }))
    names(dt) <- c("year", "nweek", "ret")
    dt
}

でもって、SP500,日経225を取得、週次に加工、一つのデータフレームにまとめる。(gggplotでの表示用に) データの期間は、2001~2012年

## sp500とnikkeiを取得してrbindする
syear = 2001
eyear = 2012
spx <- get.nweek.return("^GSPC", syear, eyear)
spx$idx <- "sp500"
nk <- get.nweek.return("^N225", syear, eyear)
nk$idx <- "n225"
df <- rbind(spx, nk)

グラフ化する。 まず平均をみて、次にboxplotでばらつきをみて、最後に年数を絞って(2010-2012)個別の週次returnをそのまま見る まず、平均。

library(ggplot2)
library(plyr)
df0 <- ddply(df, .(nweek, idx), summarise, ret = mean(ret))
ggplot(df0, aes(factor(nweek), ret)) + geom_bar(stat = "identity") + facet_grid(idx ~ 
    .) + ggtitle(paste("average return by n-week", paste(syear, eyear, sep = "-")))
## Warning: Stacking not well defined when ymin != 0
## Warning: Stacking not well defined when ymin != 0

plot of chunk unnamed-chunk-3

次にboxplotでばらつき

p <- ggplot(df, aes(factor(nweek), ret))
p + geom_boxplot() + facet_grid(idx ~ .) + ggtitle(paste("boxplot of return by n-week", 
    paste(syear, eyear, sep = "-")))

plot of chunk unnamed-chunk-4

最後に個別の年のreturnを2010,2011,2012年で

df1 <- subset(df, year > 2009)
p1 <- ggplot(df1, aes(factor(nweek), ret, color = factor(year)))  ## lineの時はgroup化, group=factor(year)))
p1 + geom_point() + facet_grid(idx ~ .) + ggtitle("return by n-week year(2010-2012)")

plot of chunk unnamed-chunk-5

ドンチャンルールのバックテスト

ドンチャンルールのバックテスト

Rで簡単なバックテストです。

(* 以下のURLを参考にしました http://d.hatena.ne.jp/teramonagi/20110329/1301352620) (ただ、ここでの間違いがあれば、それは私です。)

1.定数を設定

ここでは、120日高値抜けでLONG、30日安値抜けでCLOSE POSを-1にすれば、SHORT

D1 <- 120
D2 <- 30
POS <- 1

2.データをload。ここでは日経225の現物を使う.indicatorとpostionをcbind(column bind)する

suppressPackageStartupMessages(library(quantmod))
org.data <- getSymbols("^N225", from = "2000-01-01", auto.assign = F)
data <- adjustOHLC(org.data, use.Adjusted = T)[, 1:4]
n1 <- ifelse(POS > 0, 1, 3)
n2 <- ifelse(POS > 0, 3, 1)
data <- cbind(data, DonchianChannel(data[, 4], D1)[, n1])
data <- cbind(data, DonchianChannel(data[, 4], D2)[, n2])
data <- cbind(data, data[, 4] == data[, 5])
data <- cbind(data, data[, 4] == data[, 6])
names(data)[5:8] <- c("ind1", "ind2", "sigEntry", "sigClose")
data$pos <- 0
data <- na.omit(data)

3.for-loopで、前日までのpositionを見ながら、当日のpositionを設定

これで、トレードした(positionをとった)ことになる。

## library(foreach) for-loopが順番依存なので、parallelは使えない??
for (i in 1:(nrow(data) - 1)) {
    data$pos[i + 1] <- if (data$sigEntry[i] == 1) {
        POS
    } else if (data$sigClose[i] == 1) {
        0
    } else {
        data$pos[i]
    }
}

4.損益をグラフかする。ROCでreturnを計算。,returnを一個lagする。

suppressPackageStartupMessages(library(PerformanceAnalytics))
chart.CumReturns(cbind(ROC(data[, 4]) * (POS), ROC(data[, 4]) * lag(data$pos)), 
    wealth.index = TRUE, main = paste(ifelse(POS > 0, "BuyHold", "SellHold"), 
        "=RED, TAAproach=BLUE", sep = ""))

plot of chunk unnamed-chunk-4

5.月別returnを見てみる

cat("monthly return\n")
## monthly return
ret.yearmon <- aggregate(coredata(ROC(data[, 4])) * lag(data$pos) * POS, as.Date(as.yearmon(index(data))), 
    function(x) exp(cumsum(x[1])))
hist(coredata(ret.yearmon), breaks = 50)

plot of chunk unnamed-chunk-5

これは、未実現損益含む。ほとんどの場合、ちょい負け。
じゃなくて、ポジションもたないとき(return=0)がでてるだけか、、、
なんか、プラスがおおい、、どこかおかしいかも、、

6.関連数値(サンプル期間、エントリー数、個別トレード期間、損益、対象物のreturn)

PerformanceAnalyticsの関数を使ったほうがいいが、勉強不足

cat("total return\n")
## total return
rl <- rle(as.vector(coredata(data$pos)))
cat(paste(paste("date range", paste(range(index(data)), collapse = "-"), sep = ":"), 
    paste("number entries:", sum(rl$values), "  position-day/all-day:", sum(rl$length[rl$values == 
        1]), "/", (nrow(data) - 1), sep = ""), paste("trade periods", paste(sort(rl$length[rl$values == 
        1], decreasing = T), collapse = " "), sep = ":"), paste("return", round(exp(sum(na.omit(coredata(ROC(data[, 
        4])) * lag(data$pos)) * POS)), 3), sep = ":"), paste("B&H return", round(exp(sum(na.omit(ROC(data[, 
        4]) * POS))), 3), sep = ":"), sep = "\n"), "\n")
## date range:2000-06-26-2013-01-22
## number entries:14  position-day/all-day:816/3084
## trade periods:189 106 69 65 62 51 46 43 39 37 37 25 24 23
## return:1.133
## B&H return:0.607