分割投資 vs 一括投資

借り入れを行うという条件を許して分割投資と一括投資を比較する。

例えば、現在手持ちの資産は300万円、年間に60万円投資に回せる程度の資金の余裕があるとする。 20年後にアーリーリタイアを目指したい。

A: 毎年75万円ずつ投資して20年かけて投資額(元本)を1500万円にする
B: レバレッジを用いて450万円分借り入れて750万円分投資して、徐々にレバレッジを返済していく。
C: はじめに300万円投資して、毎年60万円ずつ積み立てる。
D: レバレッジを用いて600万円分借り入れて900万円分投資して、徐々にレバレッジを返済していく。
どれがアーリーリタイアできる見込みが高いか?

リスクの測定として、20年間運用した後の総資産額が、S(万円)(3000, 5000, 10000)を上回る確率を計算する。
ただし、金利コスト(=0.03)が20年間一定という前提を設ける。
現金は現在の日本の状況を考え、リターンを産まない(0%)資産として、インフレ率は0%として考慮していません。 投資ポートフォリオの特性として、1年あたりの超過リターンが、年率期待リターンmu, 年率期待標準偏差 sigmaの正規分布に従うi.i.d系列とする。
ポートフォリオリターンをBootstrap法により繰り返し発生させて比較を行う。

運用後の総資産の分布の比較も行う。

#ライブラリ読み込み
library(ggplot2)
library(ggsci)

投資シミュレーションの関数定義

Simulation <- function(mu, sigma, rf, year, initialPrice, income, contributeA, contributeB, contributeC, contributeD){
  # mu: 年率の期待リターン(算術平均リターン)
  # sigma: 年率の期待標準偏差
  # year: 投資する年数
  # initialPrice: 開始時の資産額
  
  # Parameters
  N <- year
  bankruptA <- 0
  bankruptB <- 0
  bankruptC <- 0
  bankruptD <- 0
  
  # Return Generation
  yearlyReturn <- rnorm(N, mu, sigma) # 超過年率リターンを正規分布(N(mu, sigma))に従うi.i.d系列として発生
  
  # 分割して投資する(A)
  totalAssetA <- rep(0, N) # 総資産額の格納
  rfAsset <- rep(0, N)    # リスクフリー資産
  investPrice <- rep(0, N) # 投資評価額の格納用
  
  rfAsset[1] <- initialPrice # はじめのcashを設定
  
  # N年分シミュレーション
  for(i in 1:N){
    # 年初に収入を計上して、投資
    rfAsset[i] <- rfAsset[i] + income
    
    investPrice[i] <- investPrice[i] + contributeA
    rfAsset[i] <- rfAsset[i] - contributeA
    
    # リターン計算
    investPrice[i] <- investPrice[i]*(yearlyReturn[i]+1)
    if(investPrice[i]+rfAsset[i]<=0){
      investPrice[i] <- 0
      rfAsset[i] <- 0
      bankruptA <- 1
      } # 0以下(破産)になったら0を代入
    totalAssetA[i] <- investPrice[i] + rfAsset[i]
    if(i == N){break}
    
    # 繰越
    investPrice[i+1] <- investPrice[i]
    rfAsset[i+1] <- rfAsset[i]
  }
  
  # 一括投資して運用し続ける(B)
  totalAssetB <- rep(0, N) # 総資産額の格納
  rfAsset <- rep(0, N)    # リスクフリー資産
  debt <- rep(0, N)    # 借金
  investPrice <- rep(0, N) # 投資評価額の格納用
  
  investPrice[1] <- contributeB
  debt[1] <- initialPrice - contributeB # はじめの借金
  rfAsset[1] <- initialPrice - debt[1] - investPrice[1] 
  # はじめのリスクフリー資産
  for(i in 1:N){
    # 収入計上と借り入れコスト算入
    rfAsset[i] <- rfAsset[i] + income
    if(debt[i]<0){ #借金あれば返済する
      d <- debt[i]
      r <- rfAsset[i]
      debt[i] <- min(0, d + r)
      rfAsset[i] <- max(0, d + r)
    } 
    debt[i] <- debt[i] + min(debt[i],0) * rf # コストを計算
    
    # リターン計算
    investPrice[i] <- investPrice[i] * (1+yearlyReturn[i])
    
    if(investPrice[i]+rfAsset[i]+debt[i]<=0){ #損失が余剰資金よりも大きくなった場合
      rfAsset[i] <- investPrice[i] + rfAsset[i]
      investPrice[i] <- 0
      bankruptB <- 1
        } # 0以下(破産)になったら0を代入
    totalAssetB[i] <- investPrice[i] + rfAsset[i] + debt[i]
    if(i==N){break}
    
    # 繰越
    investPrice[i+1] <- investPrice[i]
    rfAsset[i+1] <- rfAsset[i]
    debt[i+1] <- debt[i]
  }
  
  # 最初に一括投資して残りは積み立てる(C)
  totalAssetC <- rep(0, N) # 総資産額の格納
  rfAsset <- rep(0, N)    # リスクフリー資産(借金)
  investPrice <- rep(0, N) # 投資評価額の格納用
  
  investPrice[1] <- initialPrice
  rfAsset[1] <- 0 # はじめのリスクフリー資産
  for(i in 1:N){
    # 収入計上と積立投資
    rfAsset[i] <- rfAsset[i] + income
    
    investPrice[i] <- investPrice[i] + contributeC
    rfAsset[i] <- rfAsset[i] - contributeC
    
    # リターン計算
    investPrice[i] <- investPrice[i] * (1+yearlyReturn[i])
    if(investPrice[i]+rfAsset[i]<=0){ #損失が余剰資金よりも大きくなった場合
      investPrice[i] <- 0
      rfAsset[i] <- 0
      bankruptC <- 1
        } # 0以下(破産)になったら0を代入
    totalAssetC[i] <- investPrice[i] + rfAsset[i]
    if(i==N){break}
    
    # 繰越
    investPrice[i+1] <- investPrice[i]
    rfAsset[i+1] <- rfAsset[i]
  }
  
  # 一括投資して運用し続ける2つめ(D)
  totalAssetD <- rep(0, N) # 総資産額の格納
  rfAsset <- rep(0, N)    # リスクフリー資産
  debt <- rep(0, N)    # 借金
  investPrice <- rep(0, N) # 投資評価額の格納用
  
  investPrice[1] <- contributeD
  debt[1] <- initialPrice - contributeD # はじめの借金
  rfAsset[1] <- initialPrice - debt[1] - investPrice[1] 
  # はじめのリスクフリー資産
  for(i in 1:N){
    # 収入計上と借り入れコスト算入
    rfAsset[i] <- rfAsset[i] + income
    if(debt[i]<0){ #借金あれば返済する
      d <- debt[i]
      r <- rfAsset[i]
      debt[i] <- min(0, d + r)
      rfAsset[i] <- max(0, d + r)
    } 
      debt[i] <- debt[i] + min(debt[i],0) * rf # コストを計算
  
    # リターン計算
    investPrice[i] <- investPrice[i] * (1+yearlyReturn[i])
    if(investPrice[i]+rfAsset[i]+debt[i]<=0){ #損失が余剰資金よりも大きくなった場合
      rfAsset[i] <- investPrice[i] + rfAsset[i]
      investPrice[i] <- 0
      bankruptD <- 1
        } # 0以下(破産)になったら0を代入
    totalAssetD[i] <- investPrice[i] + rfAsset[i] + debt[i]
    if(i==N){break}
    
    # 繰越
    investPrice[i+1] <- investPrice[i]
    rfAsset[i+1] <- rfAsset[i]
    debt[i+1] <- debt[i]
  }
  
  return(c(totalAssetA[N], totalAssetB[N], totalAssetC[N], totalAssetD[N], bankruptA, bankruptB, bankruptC, bankruptD)) # A, Bの最終資産額を出力
}
Simulation(mu=0.09, sigma=0.2, rf=0.03, year=20, initialPrice=300, income=60, contributeA=75, contributeB=750, contributeC=60, contributeD = 900) # mu=0.09, sigma=0.2としてシミュレーションした結果
## [1] 3112.102 3098.660 3449.346 3379.973    0.000    0.000    0.000    0.000

シミュレーションの反復試行

Bootstrap法による反復試行は、1組の(mu, sigma)あたり100,000回とする。

repeatN <- 100000 #1組のパラメータあたりの試行数
resultA <- rep(0, repeatN) # Aの最終資産額(対数変換後)格納用
resultB <- rep(0, repeatN) # Bの最終資産額(対数変換後)格納用
resultC <- rep(0, repeatN) # Bの最終資産額(対数変換後)格納用
resultD <- rep(0, repeatN) # Bの最終資産額(対数変換後)格納用

# 以下のmu, sigmaの総当りの組み合わせで検定を行う。
m <- c(0.03, 0.06, 0.09, 0.12, 0.15) # mu(vector)
s <- c(0.05, 0.1, 0.15, 0.2, 0.25, 0.3) # sigma(vector)
rf <-  0.03

# 結果の格納用
results <- data.frame(matrix(rep(NA, 4), nrow=1))[numeric(0), ]
colnames(results) <- c("method", "mu", "sigma", "value")
bankrupt_df <- data.frame(matrix(rep(NA, 4), nrow=1))[numeric(0), ]
colnames(bankrupt_df) <- c("method", "mu", "sigma", "value")
success_df <- data.frame(matrix(rep(NA, 4), nrow=1))[numeric(0), ]
colnames(success_df) <- c("method", "mu", "sigma", "value")

# ループ
for(k in 1:length(m)){ # muについてのループ
for(j in 1:length(s)){ # sigmaについてのループ
  bankrupt <- c(0,0,0,0)
  for(i in 1:repeatN){ # bootstrap
    # 特定の(mu, sigma)についてのシミュレーション
    res <- Simulation(m[k], s[j], rf, 20, 3000000, 600000, 750000, 7500000, 600000, 9000000) # 1試行での結果(最終資産額)
    resultA[i] <- log10(max(res[1],1)) # A: 積立投資
    resultB[i] <- log10(max(res[2],1)) # B: 一括投資
    resultC[i] <- log10(max(res[3],1)) # C: 一括投資 + 積立投資
    resultD[i] <- log10(max(res[4],1)) # C: 一括投資 + 積立投資
    bankrupt <- bankrupt + res[5:8]
  }
  
  # 結果の格納
  results <- rbind(results, 
                 data.frame(method=rep("Dollar Cost",repeatN), 
                            mu=rep(m[k], repeatN), 
                            sigma=rep(s[j], repeatN), 
                            value=resultA),
                 data.frame(method=rep("Bulk",repeatN), 
                            mu=rep(m[k], repeatN), 
                            sigma=rep(s[j], repeatN), 
                            value=resultB),
                 data.frame(method=rep("Bulk & Dollar Cost",repeatN), 
                            mu=rep(m[k], repeatN), 
                            sigma=rep(s[j], repeatN), 
                            value=resultC),
                 data.frame(method=rep("Bulk-High Risk",repeatN), 
                            mu=rep(m[k], repeatN), 
                            sigma=rep(s[j], repeatN), 
                            value=resultD))
  bankrupt_df <- rbind(bankrupt_df, 
                 data.frame(method="Dollar Cost", 
                            mu=m[k], 
                            sigma=s[j], 
                            value=bankrupt[1]/repeatN),
                 data.frame(method="Bulk", 
                            mu=m[k], 
                            sigma=s[j], 
                            value=bankrupt[2]/repeatN),
                 data.frame(method="Bulk & Dollar Cost", 
                            mu=m[k], 
                            sigma=s[j], 
                            value=bankrupt[3]/repeatN),
                 data.frame(method="Bulk-High Risk", 
                            mu=m[k], 
                            sigma=s[j], 
                            value=bankrupt[4]/repeatN))
}
}
library(ggplot2)
library(ggsci)
# 箱ひげ図でそれぞれのsigmaについての中央値、四分位を示す。
plot_est_ci_bp <- function(results, sigma){
  title = paste("Final wealth distribution after 20 years in sigma = ", sigma, sep="")
  df_plot <- results[results$sigma==sigma,]
  g <- ggplot(df_plot, aes(x=as.factor(mu), y=value, fill=method))+
    geom_boxplot() +
    scale_color_nejm() +
    xlab("mu") + ylab("Wealth") +
    labs(title=title) +
    ylim(c(6,NA))
  return(g)
}

plot_est_ci_bp(results, 0.05)

plot_est_ci_bp(results, 0.1)

plot_est_ci_bp(results, 0.2)
## Warning: Removed 2 rows containing non-finite values (stat_boxplot).

plot_est_ci_bp(results, 0.3)
## Warning: Removed 1065 rows containing non-finite values (stat_boxplot).

# 破産した頻度を棒グラフで示す
plot_bankrupt <- function(results, sigma){
  title = paste("Bankrupt frequency in sigma = ", sigma, sep="")
  df_plot <- results[results$sigma==sigma,]
  g <- ggplot(df_plot, aes(x=as.factor(mu), y=value*100-0.1, fill=method))+
    geom_bar(stat="identity", position = "dodge") +
    scale_color_nejm() +
    xlab("mu") + ylab("Bankrupt (%)") +
    labs(title=title) +
    ylim(c(-0.1, 35))
  return(g)
}

plot_bankrupt(bankrupt_df, sigma=0.05)

plot_bankrupt(bankrupt_df, sigma=0.1)

plot_bankrupt(bankrupt_df, sigma=0.2)

plot_bankrupt(bankrupt_df, sigma=0.3)

目標金額を達成する可能性をプロット

library(dplyr)
## 
## Attaching package: 'dplyr'
## The following objects are masked from 'package:stats':
## 
##     filter, lag
## The following objects are masked from 'package:base':
## 
##     intersect, setdiff, setequal, union
# 箱ひげ図でそれぞれのsigmaについての中央値、四分位を示す。
plot_achievement <- function(results, sigma, thres){
  title = paste("Success rate for earning ", paste(as.character(thres/10e6),"0,000,000", sep=""), " in sigma = ", sigma, sep="")
  df_plot <- results[results$sigma==sigma,]
  g <- ggplot(df_plot, aes(x=as.factor(mu), y=value*100, fill=method))+
    geom_bar(stat="identity", position = "dodge") +
    scale_color_nejm() +
    xlab("mu") + ylab("Success (%)") +
    ylim(c(0, 100)) +
    labs(title=title)
  return(g)
}

# 目標額 1億円
thres <- 100000000
success_df <- results %>% group_by(method, mu, sigma) %>%
  mutate(success = value > log10(thres)) %>%
  summarize(value = sum(success)/repeatN)

plot_achievement(success_df, sigma=0.05, thres)

plot_achievement(success_df, sigma=0.1, thres)

plot_achievement(success_df, sigma=0.2, thres)

plot_achievement(success_df, sigma=0.3, thres)

# 目標額 5000万円
thres <- 50000000
success_df <- results %>% group_by(method, mu, sigma) %>%
  mutate(success = value > log10(thres)) %>%
  summarize(value = sum(success)/repeatN)

plot_achievement(success_df, sigma=0.05, thres)

plot_achievement(success_df, sigma=0.1, thres)

plot_achievement(success_df, sigma=0.2, thres)

plot_achievement(success_df, sigma=0.3, thres)

# 目標額 3000万円
thres <- 30000000
success_df <- results %>% group_by(method, mu, sigma) %>%
  mutate(success = value > log10(thres)) %>%
  summarize(value = sum(success)/repeatN)

plot_achievement(success_df, sigma=0.05, thres)

plot_achievement(success_df, sigma=0.1, thres)

plot_achievement(success_df, sigma=0.2, thres)

plot_achievement(success_df, sigma=0.3, thres)