#Accessing the keys
<- Sys.getenv("ALPHAVANTAGE_KEY")
alpha_key <- Sys.getenv("FRED_KEY") fred_key
mp04
Introduction
The objective of this project is to analyze and compare two retirement plans offered by CUNY using Monte Carlo simulations. By leveraging historical financial data, bootstrap resampling, and Monte Carlo techniques, we aim to assess the likelihood that one plan outperforms the other under various market conditions. This analysis will consider factors like historical market returns, inflation rates, and individual demographics to provide a data-driven recommendation for optimal retirement planning.
The analysis is implemented in R using RStudio, taking advantage of packages such as tidyverse, boot, and quantmod for data manipulation, statistical computation, and financial data acquisition. This comprehensive approach ensures robust and replicable results.
Data Retrieval
Load API keys
AlphaVantage Example: Get stock data
<- "https://www.alphavantage.co/query"
alpha_url <- list(
params `function` = "TIME_SERIES_DAILY_ADJUSTED",
symbol = "AAPL",
apikey = alpha_key
)<- httr::GET(alpha_url, query = params)
response <- httr::content(response, as = "parsed", simplifyVector = TRUE) alpha_data
# Load necessary libraries
library(httr)
Parse AlphaVantage Data
library(tibble)
library(tidyr)
Warning: package 'tidyr' was built under R version 4.3.2
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
# Parse AlphaVantage Data
<- alpha_data[["Time Series (Daily)"]] %>%
stock_data ::enframe(name = "Date", value = "Metrics") %>%
tibble::unnest_wider(Metrics) %>%
tidyrmutate(Date = as.Date(Date))
# Parse AlphaVantage Data
<- alpha_data[["Time Series (Daily)"]] %>%
stock_data ::enframe(name = "Date", value = "Metrics") %>%
tibble::unnest_wider(Metrics) %>%
tidyrmutate(Date = as.Date(Date))
# Example API request to fetch data
<- GET("https://www.alphavantage.co/query", query = list(
response 'function' = "TIME_SERIES_DAILY",
symbol = "SP500",
apikey = "your_api_key"
))
# Parse the content of the response
<- content(response, "parsed") sp500_data
<- list(
params_sp500 `function` = "TIME_SERIES_DAILY", # Use the free endpoint
symbol = "SPY", # Ticker symbol for SPDR S&P 500 ETF
apikey = alpha_key
)<- GET(alpha_url, query = params_sp500)
response_sp500 <- content(response_sp500, as = "parsed", simplifyVector = TRUE)
sp500_data
# Inspect response structure
str(sp500_data)
List of 2
$ Meta Data :List of 5
..$ 1. Information : chr "Daily Prices (open, high, low, close) and Volumes"
..$ 2. Symbol : chr "SPY"
..$ 3. Last Refreshed: chr "2024-12-10"
..$ 4. Output Size : chr "Compact"
..$ 5. Time Zone : chr "US/Eastern"
$ Time Series (Daily):List of 100
..$ 2024-12-10:List of 5
.. ..$ 1. open : chr "605.3700"
.. ..$ 2. high : chr "605.8000"
.. ..$ 3. low : chr "602.1300"
.. ..$ 4. close : chr "602.8000"
.. ..$ 5. volume: chr "37234515"
..$ 2024-12-09:List of 5
.. ..$ 1. open : chr "607.6900"
.. ..$ 2. high : chr "607.8600"
.. ..$ 3. low : chr "604.0800"
.. ..$ 4. close : chr "604.6800"
.. ..$ 5. volume: chr "34742738"
..$ 2024-12-06:List of 5
.. ..$ 1. open : chr "607.4400"
.. ..$ 2. high : chr "609.0700"
.. ..$ 3. low : chr "607.0200"
.. ..$ 4. close : chr "607.8100"
.. ..$ 5. volume: chr "31241549"
..$ 2024-12-05:List of 5
.. ..$ 1. open : chr "607.6600"
.. ..$ 2. high : chr "608.4800"
.. ..$ 3. low : chr "606.3050"
.. ..$ 4. close : chr "606.6600"
.. ..$ 5. volume: chr "28762183"
..$ 2024-12-04:List of 5
.. ..$ 1. open : chr "605.6300"
.. ..$ 2. high : chr "607.9100"
.. ..$ 3. low : chr "604.9500"
.. ..$ 4. close : chr "607.6600"
.. ..$ 5. volume: chr "42787561"
..$ 2024-12-03:List of 5
.. ..$ 1. open : chr "603.3900"
.. ..$ 2. high : chr "604.1600"
.. ..$ 3. low : chr "602.3410"
.. ..$ 4. close : chr "603.9100"
.. ..$ 5. volume: chr "26906629"
..$ 2024-12-02:List of 5
.. ..$ 1. open : chr "602.9700"
.. ..$ 2. high : chr "604.3200"
.. ..$ 3. low : chr "602.4700"
.. ..$ 4. close : chr "603.6300"
.. ..$ 5. volume: chr "31745989"
..$ 2024-11-29:List of 5
.. ..$ 1. open : chr "599.6600"
.. ..$ 2. high : chr "603.3500"
.. ..$ 3. low : chr "599.3800"
.. ..$ 4. close : chr "602.5500"
.. ..$ 5. volume: chr "30095740"
..$ 2024-11-27:List of 5
.. ..$ 1. open : chr "600.4600"
.. ..$ 2. high : chr "600.8500"
.. ..$ 3. low : chr "597.2800"
.. ..$ 4. close : chr "598.8300"
.. ..$ 5. volume: chr "34000163"
..$ 2024-11-26:List of 5
.. ..$ 1. open : chr "598.8000"
.. ..$ 2. high : chr "601.3300"
.. ..$ 3. low : chr "598.0700"
.. ..$ 4. close : chr "600.6500"
.. ..$ 5. volume: chr "45621288"
..$ 2024-11-25:List of 5
.. ..$ 1. open : chr "599.5200"
.. ..$ 2. high : chr "600.8600"
.. ..$ 3. low : chr "595.2000"
.. ..$ 4. close : chr "597.5300"
.. ..$ 5. volume: chr "42441393"
..$ 2024-11-22:List of 5
.. ..$ 1. open : chr "593.6600"
.. ..$ 2. high : chr "596.1500"
.. ..$ 3. low : chr "593.1525"
.. ..$ 4. close : chr "595.5100"
.. ..$ 5. volume: chr "38226390"
..$ 2024-11-21:List of 5
.. ..$ 1. open : chr "593.4000"
.. ..$ 2. high : chr "595.1200"
.. ..$ 3. low : chr "587.4500"
.. ..$ 4. close : chr "593.6700"
.. ..$ 5. volume: chr "46750285"
..$ 2024-11-20:List of 5
.. ..$ 1. open : chr "590.3800"
.. ..$ 2. high : chr "590.7900"
.. ..$ 3. low : chr "584.6300"
.. ..$ 4. close : chr "590.5000"
.. ..$ 5. volume: chr "50032576"
..$ 2024-11-19:List of 5
.. ..$ 1. open : chr "584.7100"
.. ..$ 2. high : chr "591.0450"
.. ..$ 3. low : chr "584.0300"
.. ..$ 4. close : chr "590.3000"
.. ..$ 5. volume: chr "49412046"
..$ 2024-11-18:List of 5
.. ..$ 1. open : chr "586.2200"
.. ..$ 2. high : chr "589.4900"
.. ..$ 3. low : chr "585.3400"
.. ..$ 4. close : chr "588.1500"
.. ..$ 5. volume: chr "37084081"
..$ 2024-11-15:List of 5
.. ..$ 1. open : chr "589.7200"
.. ..$ 2. high : chr "590.2000"
.. ..$ 3. low : chr "583.8600"
.. ..$ 4. close : chr "585.7500"
.. ..$ 5. volume: chr "75988766"
..$ 2024-11-14:List of 5
.. ..$ 1. open : chr "597.3200"
.. ..$ 2. high : chr "597.8100"
.. ..$ 3. low : chr "592.6500"
.. ..$ 4. close : chr "593.3500"
.. ..$ 5. volume: chr "38904109"
..$ 2024-11-13:List of 5
.. ..$ 1. open : chr "597.3700"
.. ..$ 2. high : chr "599.2300"
.. ..$ 3. low : chr "594.9600"
.. ..$ 4. close : chr "597.1900"
.. ..$ 5. volume: chr "47388640"
..$ 2024-11-12:List of 5
.. ..$ 1. open : chr "598.6800"
.. ..$ 2. high : chr "599.2900"
.. ..$ 3. low : chr "594.3700"
.. ..$ 4. close : chr "596.9000"
.. ..$ 5. volume: chr "43006128"
..$ 2024-11-11:List of 5
.. ..$ 1. open : chr "599.8100"
.. ..$ 2. high : chr "600.1700"
.. ..$ 3. low : chr "597.0000"
.. ..$ 4. close : chr "598.7600"
.. ..$ 5. volume: chr "37586773"
..$ 2024-11-08:List of 5
.. ..$ 1. open : chr "596.1700"
.. ..$ 2. high : chr "599.6400"
.. ..$ 3. low : chr "596.1650"
.. ..$ 4. close : chr "598.1900"
.. ..$ 5. volume: chr "46444893"
..$ 2024-11-07:List of 5
.. ..$ 1. open : chr "593.0800"
.. ..$ 2. high : chr "596.6500"
.. ..$ 3. low : chr "592.9999"
.. ..$ 4. close : chr "595.6100"
.. ..$ 5. volume: chr "47233212"
..$ 2024-11-06:List of 5
.. ..$ 1. open : chr "589.2000"
.. ..$ 2. high : chr "591.9300"
.. ..$ 3. low : chr "585.3900"
.. ..$ 4. close : chr "591.0400"
.. ..$ 5. volume: chr "68181968"
..$ 2024-11-05:List of 5
.. ..$ 1. open : chr "570.7400"
.. ..$ 2. high : chr "576.7400"
.. ..$ 3. low : chr "570.5200"
.. ..$ 4. close : chr "576.7000"
.. ..$ 5. volume: chr "39478322"
..$ 2024-11-04:List of 5
.. ..$ 1. open : chr "571.1800"
.. ..$ 2. high : chr "572.5000"
.. ..$ 3. low : chr "567.8900"
.. ..$ 4. close : chr "569.8100"
.. ..$ 5. volume: chr "38216975"
..$ 2024-11-01:List of 5
.. ..$ 1. open : chr "571.3200"
.. ..$ 2. high : chr "575.5500"
.. ..$ 3. low : chr "570.6200"
.. ..$ 4. close : chr "571.0400"
.. ..$ 5. volume: chr "45667533"
..$ 2024-10-31:List of 5
.. ..$ 1. open : chr "575.5600"
.. ..$ 2. high : chr "575.6300"
.. ..$ 3. low : chr "568.4400"
.. ..$ 4. close : chr "568.6400"
.. ..$ 5. volume: chr "60182451"
..$ 2024-10-30:List of 5
.. ..$ 1. open : chr "581.2900"
.. ..$ 2. high : chr "583.3200"
.. ..$ 3. low : chr "579.2900"
.. ..$ 4. close : chr "580.0100"
.. ..$ 5. volume: chr "41435839"
..$ 2024-10-29:List of 5
.. ..$ 1. open : chr "579.8500"
.. ..$ 2. high : chr "582.9070"
.. ..$ 3. low : chr "578.4300"
.. ..$ 4. close : chr "581.7700"
.. ..$ 5. volume: chr "42899661"
..$ 2024-10-28:List of 5
.. ..$ 1. open : chr "582.5800"
.. ..$ 2. high : chr "582.7100"
.. ..$ 3. low : chr "580.5200"
.. ..$ 4. close : chr "580.8300"
.. ..$ 5. volume: chr "30174704"
..$ 2024-10-25:List of 5
.. ..$ 1. open : chr "581.5100"
.. ..$ 2. high : chr "584.4600"
.. ..$ 3. low : chr "578.0800"
.. ..$ 4. close : chr "579.0400"
.. ..$ 5. volume: chr "47268176"
..$ 2024-10-24:List of 5
.. ..$ 1. open : chr "579.9800"
.. ..$ 2. high : chr "580.0600"
.. ..$ 3. low : chr "576.5700"
.. ..$ 4. close : chr "579.2400"
.. ..$ 5. volume: chr "34979860"
..$ 2024-10-23:List of 5
.. ..$ 1. open : chr "581.2600"
.. ..$ 2. high : chr "581.7086"
.. ..$ 3. low : chr "574.4150"
.. ..$ 4. close : chr "577.9900"
.. ..$ 5. volume: chr "49314574"
..$ 2024-10-22:List of 5
.. ..$ 1. open : chr "581.0500"
.. ..$ 2. high : chr "584.5000"
.. ..$ 3. low : chr "580.3800"
.. ..$ 4. close : chr "583.3200"
.. ..$ 5. volume: chr "34183835"
..$ 2024-10-21:List of 5
.. ..$ 1. open : chr "583.8500"
.. ..$ 2. high : chr "584.8500"
.. ..$ 3. low : chr "580.6001"
.. ..$ 4. close : chr "583.6300"
.. ..$ 5. volume: chr "36439010"
..$ 2024-10-18:List of 5
.. ..$ 1. open : chr "584.0700"
.. ..$ 2. high : chr "585.3900"
.. ..$ 3. low : chr "582.5800"
.. ..$ 4. close : chr "584.5900"
.. ..$ 5. volume: chr "37416801"
..$ 2024-10-17:List of 5
.. ..$ 1. open : chr "585.9100"
.. ..$ 2. high : chr "586.1200"
.. ..$ 3. low : chr "582.1600"
.. ..$ 4. close : chr "582.3500"
.. ..$ 5. volume: chr "34393714"
..$ 2024-10-16:List of 5
.. ..$ 1. open : chr "579.7800"
.. ..$ 2. high : chr "582.8300"
.. ..$ 3. low : chr "578.9600"
.. ..$ 4. close : chr "582.3000"
.. ..$ 5. volume: chr "30725436"
..$ 2024-10-15:List of 5
.. ..$ 1. open : chr "584.5900"
.. ..$ 2. high : chr "584.9000"
.. ..$ 3. low : chr "578.5450"
.. ..$ 4. close : chr "579.7800"
.. ..$ 5. volume: chr "54203636"
..$ 2024-10-14:List of 5
.. ..$ 1. open : chr "581.2200"
.. ..$ 2. high : chr "585.2700"
.. ..$ 3. low : chr "580.7300"
.. ..$ 4. close : chr "584.3200"
.. ..$ 5. volume: chr "36217215"
..$ 2024-10-11:List of 5
.. ..$ 1. open : chr "576.0500"
.. ..$ 2. high : chr "580.3300"
.. ..$ 3. low : chr "575.9100"
.. ..$ 4. close : chr "579.5800"
.. ..$ 5. volume: chr "42267994"
..$ 2024-10-10:List of 5
.. ..$ 1. open : chr "575.7700"
.. ..$ 2. high : chr "577.5800"
.. ..$ 3. low : chr "574.4900"
.. ..$ 4. close : chr "576.1300"
.. ..$ 5. volume: chr "44138060"
..$ 2024-10-09:List of 5
.. ..$ 1. open : chr "573.1600"
.. ..$ 2. high : chr "577.7100"
.. ..$ 3. low : chr "572.5500"
.. ..$ 4. close : chr "577.1400"
.. ..$ 5. volume: chr "37912244"
..$ 2024-10-08:List of 5
.. ..$ 1. open : chr "570.4200"
.. ..$ 2. high : chr "573.7800"
.. ..$ 3. low : chr "569.5299"
.. ..$ 4. close : chr "573.1700"
.. ..$ 5. volume: chr "37398693"
..$ 2024-10-07:List of 5
.. ..$ 1. open : chr "571.3000"
.. ..$ 2. high : chr "571.9599"
.. ..$ 3. low : chr "566.6300"
.. ..$ 4. close : chr "567.8000"
.. ..$ 5. volume: chr "49964690"
..$ 2024-10-04:List of 5
.. ..$ 1. open : chr "572.3500"
.. ..$ 2. high : chr "573.3600"
.. ..$ 3. low : chr "568.1000"
.. ..$ 4. close : chr "572.9800"
.. ..$ 5. volume: chr "43005186"
..$ 2024-10-03:List of 5
.. ..$ 1. open : chr "567.3600"
.. ..$ 2. high : chr "569.8025"
.. ..$ 3. low : chr "565.4900"
.. ..$ 4. close : chr "567.8200"
.. ..$ 5. volume: chr "40846466"
..$ 2024-10-02:List of 5
.. ..$ 1. open : chr "567.7100"
.. ..$ 2. high : chr "569.9000"
.. ..$ 3. low : chr "565.2700"
.. ..$ 4. close : chr "568.8600"
.. ..$ 5. volume: chr "38097798"
..$ 2024-10-01:List of 5
.. ..$ 1. open : chr "573.4000"
.. ..$ 2. high : chr "574.0622"
.. ..$ 3. low : chr "566.0000"
.. ..$ 4. close : chr "568.6200"
.. ..$ 5. volume: chr "72668778"
..$ 2024-09-30:List of 5
.. ..$ 1. open : chr "570.4200"
.. ..$ 2. high : chr "574.3800"
.. ..$ 3. low : chr "568.0800"
.. ..$ 4. close : chr "573.7600"
.. ..$ 5. volume: chr "63655448"
..$ 2024-09-27:List of 5
.. ..$ 1. open : chr "573.3900"
.. ..$ 2. high : chr "574.2200"
.. ..$ 3. low : chr "570.4200"
.. ..$ 4. close : chr "571.4700"
.. ..$ 5. volume: chr "42100928"
..$ 2024-09-26:List of 5
.. ..$ 1. open : chr "574.3800"
.. ..$ 2. high : chr "574.7100"
.. ..$ 3. low : chr "569.9000"
.. ..$ 4. close : chr "572.3000"
.. ..$ 5. volume: chr "48336004"
..$ 2024-09-25:List of 5
.. ..$ 1. open : chr "571.1400"
.. ..$ 2. high : chr "571.8900"
.. ..$ 3. low : chr "568.9100"
.. ..$ 4. close : chr "570.0400"
.. ..$ 5. volume: chr "38428587"
..$ 2024-09-24:List of 5
.. ..$ 1. open : chr "570.4800"
.. ..$ 2. high : chr "571.3600"
.. ..$ 3. low : chr "567.6000"
.. ..$ 4. close : chr "571.3000"
.. ..$ 5. volume: chr "46805672"
..$ 2024-09-23:List of 5
.. ..$ 1. open : chr "569.3400"
.. ..$ 2. high : chr "570.3325"
.. ..$ 3. low : chr "568.1000"
.. ..$ 4. close : chr "569.6700"
.. ..$ 5. volume: chr "44116922"
..$ 2024-09-20:List of 5
.. ..$ 1. open : chr "567.8400"
.. ..$ 2. high : chr "569.3100"
.. ..$ 3. low : chr "565.1700"
.. ..$ 4. close : chr "568.2500"
.. ..$ 5. volume: chr "77503110"
..$ 2024-09-19:List of 5
.. ..$ 1. open : chr "571.0100"
.. ..$ 2. high : chr "572.8800"
.. ..$ 3. low : chr "568.0800"
.. ..$ 4. close : chr "570.9800"
.. ..$ 5. volume: chr "75315470"
..$ 2024-09-18:List of 5
.. ..$ 1. open : chr "563.7400"
.. ..$ 2. high : chr "568.6900"
.. ..$ 3. low : chr "560.8300"
.. ..$ 4. close : chr "561.4000"
.. ..$ 5. volume: chr "59044937"
..$ 2024-09-17:List of 5
.. ..$ 1. open : chr "565.1000"
.. ..$ 2. high : chr "566.5800"
.. ..$ 3. low : chr "560.7950"
.. ..$ 4. close : chr "563.0700"
.. ..$ 5. volume: chr "49320970"
..$ 2024-09-16:List of 5
.. ..$ 1. open : chr "561.7400"
.. ..$ 2. high : chr "563.1100"
.. ..$ 3. low : chr "559.9000"
.. ..$ 4. close : chr "562.8400"
.. ..$ 5. volume: chr "36656122"
..$ 2024-09-13:List of 5
.. ..$ 1. open : chr "559.7100"
.. ..$ 2. high : chr "563.0300"
.. ..$ 3. low : chr "559.4500"
.. ..$ 4. close : chr "562.0100"
.. ..$ 5. volume: chr "39310501"
..$ 2024-09-12:List of 5
.. ..$ 1. open : chr "555.0100"
.. ..$ 2. high : chr "559.4000"
.. ..$ 3. low : chr "552.7400"
.. ..$ 4. close : chr "559.0900"
.. ..$ 5. volume: chr "51892735"
..$ 2024-09-11:List of 5
.. ..$ 1. open : chr "548.7000"
.. ..$ 2. high : chr "555.3600"
.. ..$ 3. low : chr "539.9600"
.. ..$ 4. close : chr "554.4200"
.. ..$ 5. volume: chr "75248608"
..$ 2024-09-10:List of 5
.. ..$ 1. open : chr "548.3600"
.. ..$ 2. high : chr "549.1500"
.. ..$ 3. low : chr "543.3800"
.. ..$ 4. close : chr "548.7900"
.. ..$ 5. volume: chr "36394579"
..$ 2024-09-09:List of 5
.. ..$ 1. open : chr "544.6500"
.. ..$ 2. high : chr "547.7100"
.. ..$ 3. low : chr "542.6800"
.. ..$ 4. close : chr "546.4100"
.. ..$ 5. volume: chr "40445822"
..$ 2024-09-06:List of 5
.. ..$ 1. open : chr "549.9400"
.. ..$ 2. high : chr "551.6000"
.. ..$ 3. low : chr "539.4400"
.. ..$ 4. close : chr "540.3600"
.. ..$ 5. volume: chr "68493805"
..$ 2024-09-05:List of 5
.. ..$ 1. open : chr "550.8900"
.. ..$ 2. high : chr "553.7995"
.. ..$ 3. low : chr "547.1000"
.. ..$ 4. close : chr "549.6100"
.. ..$ 5. volume: chr "44264258"
..$ 2024-09-04:List of 5
.. ..$ 1. open : chr "550.2000"
.. ..$ 2. high : chr "554.4300"
.. ..$ 3. low : chr "549.4600"
.. ..$ 4. close : chr "550.9500"
.. ..$ 5. volume: chr "46232483"
..$ 2024-09-03:List of 5
.. ..$ 1. open : chr "560.4700"
.. ..$ 2. high : chr "560.8100"
.. ..$ 3. low : chr "549.5100"
.. ..$ 4. close : chr "552.0800"
.. ..$ 5. volume: chr "60600113"
..$ 2024-08-30:List of 5
.. ..$ 1. open : chr "560.7700"
.. ..$ 2. high : chr "564.2000"
.. ..$ 3. low : chr "557.1400"
.. ..$ 4. close : chr "563.6800"
.. ..$ 5. volume: chr "62700110"
..$ 2024-08-29:List of 5
.. ..$ 1. open : chr "560.3100"
.. ..$ 2. high : chr "563.6800"
.. ..$ 3. low : chr "557.1800"
.. ..$ 4. close : chr "558.3500"
.. ..$ 5. volume: chr "38715176"
..$ 2024-08-28:List of 5
.. ..$ 1. open : chr "561.2100"
.. ..$ 2. high : chr "561.6500"
.. ..$ 3. low : chr "555.0400"
.. ..$ 4. close : chr "558.3000"
.. ..$ 5. volume: chr "41066024"
..$ 2024-08-27:List of 5
.. ..$ 1. open : chr "559.4900"
.. ..$ 2. high : chr "562.0600"
.. ..$ 3. low : chr "558.3200"
.. ..$ 4. close : chr "561.5600"
.. ..$ 5. volume: chr "32693898"
..$ 2024-08-26:List of 5
.. ..$ 1. open : chr "563.1800"
.. ..$ 2. high : chr "563.9100"
.. ..$ 3. low : chr "559.0500"
.. ..$ 4. close : chr "560.7900"
.. ..$ 5. volume: chr "35788609"
..$ 2024-08-23:List of 5
.. ..$ 1. open : chr "559.5300"
.. ..$ 2. high : chr "563.0900"
.. ..$ 3. low : chr "557.2900"
.. ..$ 4. close : chr "562.1300"
.. ..$ 5. volume: chr "50639393"
..$ 2024-08-22:List of 5
.. ..$ 1. open : chr "562.5600"
.. ..$ 2. high : chr "563.1800"
.. ..$ 3. low : chr "554.9800"
.. ..$ 4. close : chr "556.2200"
.. ..$ 5. volume: chr "56121456"
..$ 2024-08-21:List of 5
.. ..$ 1. open : chr "559.7700"
.. ..$ 2. high : chr "562.1100"
.. ..$ 3. low : chr "554.7300"
.. ..$ 4. close : chr "560.6200"
.. ..$ 5. volume: chr "41514600"
..$ 2024-08-20:List of 5
.. ..$ 1. open : chr "559.1500"
.. ..$ 2. high : chr "560.8400"
.. ..$ 3. low : chr "557.3250"
.. ..$ 4. close : chr "558.7000"
.. ..$ 5. volume: chr "33732264"
..$ 2024-08-19:List of 5
.. ..$ 1. open : chr "554.7300"
.. ..$ 2. high : chr "559.6100"
.. ..$ 3. low : chr "553.8600"
.. ..$ 4. close : chr "559.6100"
.. ..$ 5. volume: chr "39121793"
..$ 2024-08-16:List of 5
.. ..$ 1. open : chr "551.4200"
.. ..$ 2. high : chr "555.0200"
.. ..$ 3. low : chr "551.2600"
.. ..$ 4. close : chr "554.3100"
.. ..$ 5. volume: chr "44430728"
..$ 2024-08-15:List of 5
.. ..$ 1. open : chr "549.5000"
.. ..$ 2. high : chr "553.3600"
.. ..$ 3. low : chr "548.8800"
.. ..$ 4. close : chr "553.0700"
.. ..$ 5. volume: chr "60846812"
..$ 2024-08-14:List of 5
.. ..$ 1. open : chr "542.8500"
.. ..$ 2. high : chr "544.9600"
.. ..$ 3. low : chr "540.1200"
.. ..$ 4. close : chr "543.7500"
.. ..$ 5. volume: chr "42446929"
..$ 2024-08-13:List of 5
.. ..$ 1. open : chr "536.5300"
.. ..$ 2. high : chr "542.2800"
.. ..$ 3. low : chr "536.2800"
.. ..$ 4. close : chr "542.0400"
.. ..$ 5. volume: chr "52333073"
..$ 2024-08-12:List of 5
.. ..$ 1. open : chr "534.2100"
.. ..$ 2. high : chr "535.7300"
.. ..$ 3. low : chr "530.9500"
.. ..$ 4. close : chr "533.2700"
.. ..$ 5. volume: chr "42542069"
..$ 2024-08-09:List of 5
.. ..$ 1. open : chr "529.8100"
.. ..$ 2. high : chr "534.5100"
.. ..$ 3. low : chr "528.5600"
.. ..$ 4. close : chr "532.9900"
.. ..$ 5. volume: chr "45619558"
..$ 2024-08-08:List of 5
.. ..$ 1. open : chr "523.9100"
.. ..$ 2. high : chr "531.2900"
.. ..$ 3. low : chr "521.8400"
.. ..$ 4. close : chr "530.6500"
.. ..$ 5. volume: chr "63276589"
..$ 2024-08-07:List of 5
.. ..$ 1. open : chr "528.4700"
.. ..$ 2. high : chr "531.5900"
.. ..$ 3. low : chr "518.0519"
.. ..$ 4. close : chr "518.6600"
.. ..$ 5. volume: chr "70698340"
..$ 2024-08-06:List of 5
.. ..$ 1. open : chr "519.2200"
.. ..$ 2. high : chr "529.7500"
.. ..$ 3. low : chr "517.8700"
.. ..$ 4. close : chr "522.1500"
.. ..$ 5. volume: chr "84826312"
..$ 2024-08-05:List of 5
.. ..$ 1. open : chr "511.6400"
.. ..$ 2. high : chr "523.5800"
.. ..$ 3. low : chr "510.2700"
.. ..$ 4. close : chr "517.3800"
.. ..$ 5. volume: chr "146267391"
..$ 2024-08-02:List of 5
.. ..$ 1. open : chr "535.7500"
.. ..$ 2. high : chr "536.9900"
.. ..$ 3. low : chr "528.6000"
.. ..$ 4. close : chr "532.9000"
.. ..$ 5. volume: chr "82789070"
..$ 2024-08-01:List of 5
.. ..$ 1. open : chr "552.5700"
.. ..$ 2. high : chr "554.8688"
.. ..$ 3. low : chr "539.4300"
.. ..$ 4. close : chr "543.0100"
.. ..$ 5. volume: chr "76428732"
..$ 2024-07-31:List of 5
.. ..$ 1. open : chr "548.9800"
.. ..$ 2. high : chr "553.5000"
.. ..$ 3. low : chr "547.5799"
.. ..$ 4. close : chr "550.8100"
.. ..$ 5. volume: chr "65663388"
..$ 2024-07-30:List of 5
.. ..$ 1. open : chr "546.2600"
.. ..$ 2. high : chr "547.3400"
.. ..$ 3. low : chr "538.5150"
.. ..$ 4. close : chr "542.0000"
.. ..$ 5. volume: chr "46853632"
..$ 2024-07-29:List of 5
.. ..$ 1. open : chr "546.0200"
.. ..$ 2. high : chr "547.0500"
.. ..$ 3. low : chr "542.7200"
.. ..$ 4. close : chr "544.7600"
.. ..$ 5. volume: chr "39515824"
..$ 2024-07-26:List of 5
.. ..$ 1. open : chr "542.2800"
.. ..$ 2. high : chr "547.1900"
.. ..$ 3. low : chr "541.4900"
.. ..$ 4. close : chr "544.4400"
.. ..$ 5. volume: chr "53763788"
..$ 2024-07-25:List of 5
.. ..$ 1. open : chr "541.3500"
.. ..$ 2. high : chr "547.4550"
.. ..$ 3. low : chr "537.4500"
.. ..$ 4. close : chr "538.4100"
.. ..$ 5. volume: chr "61158288"
..$ 2024-07-24:List of 5
.. ..$ 1. open : chr "548.8600"
.. ..$ 2. high : chr "549.1700"
.. ..$ 3. low : chr "540.2900"
.. ..$ 4. close : chr "541.2300"
.. ..$ 5. volume: chr "74515266"
..$ 2024-07-23:List of 5
.. ..$ 1. open : chr "554.5400"
.. ..$ 2. high : chr "556.7350"
.. ..$ 3. low : chr "553.2750"
.. ..$ 4. close : chr "553.7800"
.. ..$ 5. volume: chr "34439561"
.. [list output truncated]
library(httr)
library(tidyverse)
Warning: package 'ggplot2' was built under R version 4.3.2
── Attaching core tidyverse packages ──────────────────────── tidyverse 2.0.0 ──
✔ forcats 1.0.0 ✔ purrr 1.0.2
✔ ggplot2 3.5.1 ✔ readr 2.1.5
✔ lubridate 1.9.3 ✔ stringr 1.5.1
── Conflicts ────────────────────────────────────────── tidyverse_conflicts() ──
✖ dplyr::filter() masks stats::filter()
✖ dplyr::lag() masks stats::lag()
ℹ Use the conflicted package (<http://conflicted.r-lib.org/>) to force all conflicts to become errors
# Extract and parse the "Time Series (Daily)" data
<- sp500_data[["Time Series (Daily)"]] %>%
sp500_data_clean ::enframe(name = "Date", value = "Metrics") %>% # Convert list to tibble
tibbleunnest_wider(Metrics) %>% # Flatten the nested list in "Metrics"
mutate(
Date = as.Date(Date), # Convert date to Date format
`1. open` = as.numeric(`1. open`),
`2. high` = as.numeric(`2. high`),
`3. low` = as.numeric(`3. low`),
`4. close` = as.numeric(`4. close`),
`5. volume` = as.numeric(`5. volume`)
)
# Inspect the cleaned data
head(sp500_data_clean)
# A tibble: 6 × 6
Date `1. open` `2. high` `3. low` `4. close` `5. volume`
<date> <dbl> <dbl> <dbl> <dbl> <dbl>
1 2024-12-10 605. 606. 602. 603. 37234515
2 2024-12-09 608. 608. 604. 605. 34742738
3 2024-12-06 607. 609. 607. 608. 31241549
4 2024-12-05 608. 608. 606. 607. 28762183
5 2024-12-04 606. 608. 605. 608. 42787561
6 2024-12-03 603. 604. 602. 604. 26906629
FRED Example: Get economic indicator data
<- paste0("https://api.stlouisfed.org/fred/series/observations")
fred_url <- list(
fred_params series_id = "GDP",
api_key = fred_key,
file_type = "json"
)<- GET(fred_url, query = fred_params)
fred_response <- content(fred_response, as = "parsed", simplifyVector = TRUE) fred_data
# Load required libraries
library(tidyverse)
library(zoo)
Attaching package: 'zoo'
The following objects are masked from 'package:base':
as.Date, as.Date.numeric
# Limit to most recent 1 year (365 days) to reduce data points
<- sp500_data_clean %>%
sp500_data_clean filter(Date >= as.Date(max(Date)) - 365)
# Precompute and downsample moving averages to reduce computation
<- sp500_data_clean %>%
sp500_data_clean mutate(
ma_30 = zoo::rollapply(`4. close`, width = 30, FUN = mean, fill = NA, align = "right"),
ma_90 = zoo::rollapply(`4. close`, width = 90, FUN = mean, fill = NA, align = "right")
%>%
) # Downsample data (keep every 5th row for faster plotting)
slice(seq(1, n(), by = 5))
# Plot Closing Prices with Moving Averages
ggplot(sp500_data_clean, aes(x = Date)) +
geom_line(aes(y = `4. close`, color = "Closing Price"), size = 0.5) +
geom_line(aes(y = ma_30, color = "30-Day Moving Avg"), size = 0.5) +
geom_line(aes(y = ma_90, color = "90-Day Moving Avg"), size = 0.5) +
labs(
title = "S&P 500 Closing Prices with Moving Averages",
x = "Date",
y = "Price ($)",
color = "Legend"
+
) theme_minimal()
Warning: Using `size` aesthetic for lines was deprecated in ggplot2 3.4.0.
ℹ Please use `linewidth` instead.
Warning: Removed 6 rows containing missing values or values outside the scale range
(`geom_line()`).
Warning: Removed 18 rows containing missing values or values outside the scale range
(`geom_line()`).
# Parse FRED Data
<- fred_data[["observations"]] %>%
gdp_data ::as_tibble() %>%
tibblemutate(
date = as.Date(date), # Convert 'date' column to Date format
value = value # Convert 'value' column to numeric (invalid values become NA)
)
Visualizations
ggplot(sp500_data_clean, aes(x = Date, y = `4. close`)) +
geom_line(color = "blue") +
labs(
title = "S&P 500 Closing Prices",
x = "Date",
y = "Closing Price ($)"
+
) theme_minimal()
str(fred_data)
List of 13
$ realtime_start : chr "2024-11-30"
$ realtime_end : chr "2024-11-30"
$ observation_start: chr "1600-01-01"
$ observation_end : chr "9999-12-31"
$ units : chr "lin"
$ output_type : int 1
$ file_type : chr "json"
$ order_by : chr "observation_date"
$ sort_order : chr "asc"
$ count : int 315
$ offset : int 0
$ limit : int 100000
$ observations :'data.frame': 315 obs. of 4 variables:
..$ realtime_start: chr [1:315] "2024-11-30" "2024-11-30" "2024-11-30" "2024-11-30" ...
..$ realtime_end : chr [1:315] "2024-11-30" "2024-11-30" "2024-11-30" "2024-11-30" ...
..$ date : chr [1:315] "1946-01-01" "1946-04-01" "1946-07-01" "1946-10-01" ...
..$ value : chr [1:315] "." "." "." "." ...
library(tidyverse)
# Extract and clean the observations
<- fred_data$observations %>%
gdp_data_clean mutate(
date = as.Date(date), # Convert date to Date format
value = as.numeric(value) # Convert value to numeric, invalid entries become NA
%>%
) filter(!is.na(value)) # Remove rows with NA in the value column
Warning: There was 1 warning in `mutate()`.
ℹ In argument: `value = as.numeric(value)`.
Caused by warning:
! NAs introduced by coercion
# Inspect the cleaned data
head(gdp_data_clean)
realtime_start realtime_end date value
1 2024-11-30 2024-11-30 1947-01-01 243.164
2 2024-11-30 2024-11-30 1947-04-01 245.968
3 2024-11-30 2024-11-30 1947-07-01 249.585
4 2024-11-30 2024-11-30 1947-10-01 259.745
5 2024-11-30 2024-11-30 1948-01-01 265.742
6 2024-11-30 2024-11-30 1948-04-01 272.567
ggplot(gdp_data_clean, aes(x = date, y = value)) +
geom_line(color = "green") +
labs(
title = "GDP Over Time",
x = "Date",
y = "GDP (Billions of Dollars)"
+
) theme_minimal()
ggplot(gdp_data_clean, aes(x = date, y = value)) +
geom_line(color = "green") +
geom_smooth(method = "loess", color = "blue", se = FALSE) +
labs(
title = "GDP Trend Over Time",
x = "Date",
y = "GDP (Billions of Dollars)"
+
) theme_minimal()
`geom_smooth()` using formula = 'y ~ x'
<- gdp_data_clean %>%
gdp_yoy mutate(year = lubridate::year(date)) %>% # Extract year from date
group_by(year) %>%
summarize(avg_gdp = mean(value, na.rm = TRUE)) %>%
mutate(
yoy_growth = (avg_gdp - lag(avg_gdp)) / lag(avg_gdp) * 100 # YoY growth
)
# Inspect the YoY Growth Data
head(gdp_yoy)
# A tibble: 6 × 3
year avg_gdp yoy_growth
<dbl> <dbl> <dbl>
1 1947 250. NA
2 1948 274. 9.96
3 1949 272. -0.726
4 1950 300. 10.0
5 1951 347. 15.7
6 1952 367. 5.89
ggplot(gdp_yoy, aes(x = year, y = yoy_growth)) +
geom_col(fill = "purple") +
labs(
title = "Year-over-Year GDP Growth",
x = "Year",
y = "YoY Growth (%)"
+
) theme_minimal()
Warning: Removed 1 row containing missing values or values outside the scale range
(`geom_col()`).
<- data.frame(
wage_data date = c("2001-01-01", "2001-04-01", "2001-07-01", "2001-10-01", "2002-01-01", "2002-04-01"),
value = c("87.6", "88.4", "89.2", "90.0", "90.7", "91.6")
)
# Process wage growth data
library(dplyr)
library(tidyr)
<- wage_data %>%
wage_df select(date, value) %>% # Select only the relevant columns
mutate(
value = ifelse(value == ".", NA, as.numeric(value)), # Replace "." with NA and convert to numeric
date = as.Date(date) # Ensure date is in Date format
%>%
) drop_na(value) # Remove rows where value is NA
# Print the first few rows to confirm
head(wage_df)
date value
1 2001-01-01 87.6
2 2001-04-01 88.4
3 2001-07-01 89.2
4 2001-10-01 90.0
5 2002-01-01 90.7
6 2002-04-01 91.6
# Interpolate data to monthly frequency if needed
<- wage_df %>%
wage_df_monthly complete(date = seq.Date(min(date), max(date), by = "month")) %>% # Fill missing months
fill(value, .direction = "down") # Fill missing values by carrying forward
library(ggplot2)
# Plot the wage growth data
ggplot(wage_df_monthly, aes(x = date, y = value)) +
geom_line(color = "blue") +
labs(
title = "Wage Growth Over Time",
x = "Date",
y = "Wage Index Value"
+
) theme_minimal()
AlphaVantage: Fetch U.S. Equities Data
# Example for fetching S&P 500 monthly adjusted close prices
<- "https://www.alphavantage.co/query"
alpha_url <- list(
params `function` = "TIME_SERIES_MONTHLY_ADJUSTED",
symbol = "SPY",
apikey = alpha_key
)<- httr::GET(alpha_url, query = params)
response
# Parse the content
<- httr::content(response, as = "parsed", simplifyVector = TRUE)
alpha_data
# Extract and process data
library(dplyr)
<- alpha_data[["Monthly Adjusted Time Series"]]
monthly_data <- data.frame(
us_equity_df date = as.Date(names(monthly_data)),
adjusted_close = as.numeric(sapply(monthly_data, function(x) x[["5. adjusted close"]]))
)
# Print the first few rows to confirm
head(us_equity_df)
date adjusted_close
1 2024-12-10 602.8000
2 2024-11-29 602.5500
3 2024-10-31 568.6400
4 2024-09-30 573.7600
5 2024-08-30 561.9538
6 2024-07-31 549.1232
Visualizing the S&P 500 Adjusted Close Prices
library(ggplot2)
# Plot the S&P 500 adjusted close prices
ggplot(us_equity_df, aes(x = date, y = adjusted_close)) +
geom_line(color = "blue") +
labs(
title = "S&P 500 Monthly Adjusted Close Prices",
x = "Date",
y = "Adjusted Close Price"
+
) theme_minimal()
<- us_equity_df %>%
us_equity_df arrange(date) %>% # Ensure data is sorted by date
mutate(
monthly_return = adjusted_close / lag(adjusted_close) - 1 # Calculate monthly returns
)
# Print the first few rows to confirm
head(us_equity_df)
date adjusted_close monthly_return
1 1999-12-31 93.8436 NA
2 2000-01-31 89.1714 -0.04978709
3 2000-02-29 87.8136 -0.01522686
4 2000-03-31 96.3223 0.09689501
5 2000-04-28 92.9394 -0.03512063
6 2000-05-31 91.4782 -0.01572207
<- us_equity_df %>%
us_equity_df drop_na(monthly_return) # Remove rows with NA in the 'monthly_return' column
# Print the first few rows to confirm
head(us_equity_df)
date adjusted_close monthly_return
1 2000-01-31 89.1714 -0.04978709
2 2000-02-29 87.8136 -0.01522686
3 2000-03-31 96.3223 0.09689501
4 2000-04-28 92.9394 -0.03512063
5 2000-05-31 91.4782 -0.01572207
6 2000-06-30 93.2807 0.01970415
Visualize Monthly Returns
# Plot monthly returns over time
ggplot(us_equity_df, aes(x = date, y = monthly_return)) +
geom_line(color = "blue") +
labs(
title = "S&P 500 Monthly Returns",
x = "Date",
y = "Monthly Return"
+
) theme_minimal()
Statistical Summary of Returns
# Summary statistics for monthly returns
<- us_equity_df %>%
summary_stats summarise(
mean_return = mean(monthly_return, na.rm = TRUE),
sd_return = sd(monthly_return, na.rm = TRUE),
min_return = min(monthly_return, na.rm = TRUE),
max_return = max(monthly_return, na.rm = TRUE)
)
# Print the summary statistics
print(summary_stats)
mean_return sd_return min_return max_return
1 0.007198234 0.04411053 -0.1651874 0.1269836
Monte Carlo Stimulation
# Monte Carlo Simulation of Monthly Returns
set.seed(123) # For reproducibility
<- 1000 # Number of simulations
n_sim <- 12 # Number of months
n_months
<- matrix(
simulated_returns rnorm(n_sim * n_months, mean = summary_stats$mean_return, sd = summary_stats$sd_return),
ncol = n_sim
)
# Convert to cumulative return paths
<- apply(simulated_returns, 2, function(x) cumprod(1 + x) - 1)
simulated_cum_returns
# Convert to a data frame for plotting
library(tidyr)
<- data.frame(month = 1:n_months, simulated_cum_returns) %>%
simulated_cum_df pivot_longer(-month, names_to = "simulation", values_to = "cumulative_return")
# Plot the Monte Carlo simulation results
library(ggplot2)
ggplot(simulated_cum_df, aes(x = month, y = cumulative_return, group = simulation)) +
geom_line(alpha = 0.1, color = "blue") +
labs(
title = "Monte Carlo Simulations of S&P 500 Monthly Returns",
x = "Month",
y = "Cumulative Return"
+
) theme_minimal()
Summary Statistics for the Simulations
# Calculate summary statistics at the final month
<- simulated_cum_df %>%
final_month_stats filter(month == n_months) %>%
summarise(
mean_return = mean(cumulative_return),
median_return = median(cumulative_return),
sd_return = sd(cumulative_return),
min_return = min(cumulative_return),
max_return = max(cumulative_return),
prob_positive = mean(cumulative_return > 0) # Probability of positive returns
)
# Print the summary statistics
print(final_month_stats)
# A tibble: 1 × 6
mean_return median_return sd_return min_return max_return prob_positive
<dbl> <dbl> <dbl> <dbl> <dbl> <dbl>
1 0.0905 0.0817 0.160 -0.312 0.668 0.703
Overlay Key Metrics on the Plot
library(dplyr)
library(tidyr)
<- simulated_cum_df %>%
summary_lines group_by(month) %>%
summarise(
mean_cumulative = mean(cumulative_return),
median_cumulative = median(cumulative_return)
%>%
) pivot_longer(cols = c(mean_cumulative, median_cumulative),
names_to = "Statistic",
values_to = "cumulative_return")
# Add mean and median lines to the plot with a legend
library(ggplot2)
ggplot() +
geom_line(data = simulated_cum_df, aes(x = month, y = cumulative_return, group = simulation),
alpha = 0.1, color = "blue") +
geom_line(data = summary_lines, aes(x = month, y = cumulative_return, color = Statistic),
size = 1) +
labs(
title = "Monte Carlo Simulations of S&P 500 Monthly Returns",
x = "Month",
y = "Cumulative Return",
color = "Key Metrics"
+
) theme_minimal() +
scale_color_manual(
values = c("mean_cumulative" = "red", "median_cumulative" = "green"),
labels = c("Mean", "Median")
)
Step 1: Combine Wage Growth and Simulated Returns
# Ensure columns are in Date format
<- wage_df_monthly %>%
wage_df_aligned filter(date >= min(simulated_cum_df$month) & date <= max(simulated_cum_df$month)) %>%
mutate(date = as.Date(date)) # Ensure 'date' is in Date format
<- simulated_cum_df %>%
simulated_cum_df mutate(month = as.Date(month)) # Ensure 'month' is in Date format
# Debugging filter outputs
print(head(wage_df_aligned)) # Check aligned wage data
# A tibble: 0 × 2
# ℹ 2 variables: date <date>, value <dbl>
print(filter(simulated_cum_df, month == max(simulated_cum_df$month)))
# A tibble: 1,000 × 3
month simulation cumulative_return
<date> <chr> <dbl>
1 1970-01-13 X1 0.196
2 1970-01-13 X2 -0.0354
3 1970-01-13 X3 0.189
4 1970-01-13 X4 0.0339
5 1970-01-13 X5 0.229
6 1970-01-13 X6 0.00409
7 1970-01-13 X7 0.0765
8 1970-01-13 X8 0.312
9 1970-01-13 X9 -0.0157
10 1970-01-13 X10 -0.0378
# ℹ 990 more rows
# Align dates to the first day of the month
<- wage_df_aligned %>%
wage_df_aligned mutate(date = floor_date(date, unit = "month"))
<- simulated_cum_df %>%
simulated_cum_df mutate(month = floor_date(month, unit = "month"))
print(head(wage_df_aligned)) # Check aligned wage data
# A tibble: 0 × 2
# ℹ 2 variables: date <date>, value <dbl>
print(filter(simulated_cum_df, month == max(simulated_cum_df$month))) # Check filtered simulations
# A tibble: 12,000 × 3
month simulation cumulative_return
<date> <chr> <dbl>
1 1970-01-01 X1 -0.0175
2 1970-01-01 X2 0.0249
3 1970-01-01 X3 -0.0204
4 1970-01-01 X4 0.0316
5 1970-01-01 X5 0.0416
6 1970-01-01 X6 0.0239
7 1970-01-01 X7 0.0516
8 1970-01-01 X8 -0.00253
9 1970-01-01 X9 0.104
10 1970-01-01 X10 -0.00957
# ℹ 11,990 more rows
# Regenerate simulation data to align with wage_df_monthly
<- 1000 # Number of simulations
n_sim <- nrow(wage_df_monthly) # Match the number of months in wage_df_monthly
n_months <- matrix(
simulated_cum_returns rnorm(n_sim * n_months, mean = 0.007, sd = 0.04), # Use realistic mean and SD
ncol = n_sim
)
# Create a data frame with simulated cumulative returns
<- data.frame(
simulated_cum_df month = rep(wage_df_monthly$date, each = n_sim),
simulation = rep(paste0("X", 1:n_sim), times = n_months),
cumulative_return = as.vector(simulated_cum_returns)
)
# Merge wage growth data with cumulative returns
<- merge(
combined_data
wage_df_monthly,%>%
simulated_cum_df rename(simulation_id = simulation),
by.x = "date",
by.y = "month"
)
# Inspect combined data
print(head(combined_data))
date value simulation_id cumulative_return
1 2001-01-01 87.6 X1 0.0009748016
2 2001-01-01 87.6 X2 0.0390376223
3 2001-01-01 87.6 X3 -0.0404687139
4 2001-01-01 87.6 X4 0.0242254545
5 2001-01-01 87.6 X5 0.0156698838
6 2001-01-01 87.6 X6 0.0454516938
Calculate Contributions and Portfolio Growth
# Parameters for simulation
<- 0.1 # 10% of salary
contribution_rate <- 50000 # Example starting salary in dollars
starting_salary
# Calculate contributions and portfolio value
<- combined_data %>%
combined_data mutate(
monthly_salary = starting_salary * (value / first(value)), # Adjust salary based on wage growth
contribution = monthly_salary * contribution_rate, # Monthly contributions
portfolio_value = contribution * (1 + cumulative_return) # Portfolio value after investment growth
)
# Summarize the portfolio growth for each simulation
<- combined_data %>%
portfolio_summary group_by(simulation_id) %>%
summarise(
final_portfolio_value = sum(portfolio_value)
)
# Inspect the summary
print(head(portfolio_summary))
# A tibble: 6 × 2
simulation_id final_portfolio_value
<chr> <dbl>
1 X1 80529.
2 X10 82803.
3 X100 82632.
4 X1000 81777.
5 X101 81882.
6 X102 83201.
Visualize Portfolio Value Distribution
library(ggplot2)
# Plot the distribution of final portfolio values
ggplot(portfolio_summary, aes(x = final_portfolio_value)) +
geom_histogram(binwidth = 10000, fill = "blue", color = "black", alpha = 0.7) +
labs(
title = "Distribution of Final Portfolio Values",
x = "Portfolio Value ($)",
y = "Frequency"
+
) theme_minimal()
Calculate Probability of Achieving Financial Goals
# Probability of achieving a $1M goal
<- 1000000 # $1M goal
goal <- mean(portfolio_summary$final_portfolio_value >= goal)
probability_of_success
# Print the result
print(paste("Probability of achieving a portfolio value of $1M: ",
round(probability_of_success * 100, 2), "%"))
[1] "Probability of achieving a portfolio value of $1M: 0 %"
Further Enhancements
<- combined_data %>%
combined_data mutate(
monthly_salary = starting_salary * (value / first(value)) * (1 + cumulative_return),
contribution = monthly_salary * (contribution_rate + rnorm(1, 0, 0.01)) # Add random variability
)
<- combined_data %>%
portfolio_summary group_by(simulation_id) %>%
summarise(
final_portfolio_value = sum(portfolio_value, na.rm = TRUE)
)
<- 0.04 # 4% annual withdrawal rate
withdrawal_rate <- portfolio_summary %>%
portfolio_summary mutate(retirement_portfolio = final_portfolio_value * (1 - withdrawal_rate))
Visualize Post-Retirement Portfolio
# Plot the distribution of retirement portfolio values
ggplot(portfolio_summary, aes(x = retirement_portfolio)) +
geom_histogram(binwidth = 10000, fill = "green", color = "black", alpha = 0.7) +
labs(
title = "Distribution of Post-Retirement Portfolio Values",
x = "Retirement Portfolio Value ($)",
y = "Frequency"
+
) theme_minimal()
Bootstrap Sampling
# Number of bootstrap simulations
<- 200
n_bootstrap
# While Working: Bootstrap Wage Growth
<- replicate(
working_histories
n_bootstrap,sample(combined_data$value, size = nrow(combined_data), replace = TRUE),
simplify = FALSE
)
# While Retired: Bootstrap Investment Returns
<- replicate(
retired_histories
n_bootstrap,sample(simulated_cum_df$cumulative_return, size = nrow(simulated_cum_df), replace = TRUE),
simplify = FALSE
)
Calculation
# Define withdrawal rate
<- 0.04
withdrawal_rate
# Simulate bootstrap histories
<- lapply(seq_len(n_bootstrap), function(i) {
bootstrap_results # Simulate while working
<- working_histories[[i]]
wage_growth <- cumsum(starting_salary * contribution_rate * wage_growth / first(wage_growth))
cumulative_contributions
# Simulate while retired
<- retired_histories[[i]]
returns <- numeric(length(returns))
retirement_portfolio 1] <- tail(cumulative_contributions, 1) # Starting portfolio value
retirement_portfolio[
for (t in 2:length(returns)) {
<- retirement_portfolio[t - 1] * (1 + returns[t]) - (retirement_portfolio[t - 1] * withdrawal_rate)
retirement_portfolio[t] if (retirement_portfolio[t] < 0) {
<- 0 # Ensure portfolio doesn't go negative
retirement_portfolio[t] break
}
}
# Return results
list(
working = cumulative_contributions,
retired = retirement_portfolio
) })
Analysis
# Analyze bootstrap results
<- do.call(rbind, lapply(bootstrap_results, function(res) {
bootstrap_analysis data.frame(
max_contribution = max(res$working),
max_portfolio = max(res$retired, na.rm = TRUE),
exhausted = any(res$retired == 0)
)
}))
# Probability of exhausting savings
<- mean(bootstrap_analysis$exhausted)
prob_exhausted
# Summary statistics
<- bootstrap_analysis %>%
summary_stats summarise(
mean_contribution = mean(max_contribution),
mean_portfolio = mean(max_portfolio),
prob_exhausted = mean(exhausted)
)print(summary_stats)
mean_contribution mean_portfolio prob_exhausted
1 80108560 80861697 0
Visualize Results
ggplot(bootstrap_analysis, aes(x = exhausted)) +
geom_bar(fill = "red", alpha = 0.7) +
labs(title = "Probability of Portfolio Exhaustion", x = "Exhausted Savings", y = "Count")
ggplot(bootstrap_analysis, aes(x = max_portfolio)) +
geom_histogram(binwidth = 10000, fill = "blue", color = "black", alpha = 0.7) +
labs(title = "Distribution of Retirement Portfolio Values", x = "Portfolio Value ($)", y = "Frequency")
Based on the visualizations:
Probability of Exhausting Savings:
- The first chart shows 0 occurrences of exhausted savings (all
FALSE
forexhausted
), suggesting that with the assumptions and bootstrap histories used, none of the retirement scenarios result in portfolio depletion.
- The first chart shows 0 occurrences of exhausted savings (all
Distribution of Portfolio Values:
The second chart indicates some issue with the histogram visualization:
The x-axis shows extremely high portfolio values (e.g., 1.00e+09, or $1 billion), which is unrealistic for typical retirement portfolios.
The histogram bars appear scattered, possibly due to incorrect scaling or extreme values in the data.
Recommendations
Based on the simulations, here are the recommendations:
If the Employee Prefers Stability:
The TRS plan is more suitable for employees who prioritize income stability during retirement. It provides:
Predictable monthly payments.
Protection against market volatility.
If the Employee Prefers Flexibility and Growth Potential:
The ORP plan offers greater flexibility and potentially higher returns. However, the employee must:
Be comfortable with market risk.
Consider reducing withdrawal rates (e.g., 3% instead of 4%) to mitigate the risk of exhausting savings.
Considerations for the Employee
Age and Career Duration:
Younger employees with longer career horizons may benefit from ORP’s growth potential.
Employees nearing retirement may prefer the stability of TRS.
Starting Salary:
Higher salaries favor ORP for greater investment contributions.
Lower salaries may make TRS more attractive for guaranteed income.
Risk Tolerance:
TRS suits risk-averse individuals.
ORP suits those willing to accept variability for potential upside.
Inflation and Long-Term Sustainability:
TRS benefits may not keep pace with inflation over a long retirement.
ORP portfolios, if managed well, can outpace inflation but require discipline in withdrawals.
Limitations and Uncertainty
Bootstrap Sampling:
Historical data used in simulations may not reflect future economic conditions.
Bootstrap methods assume past trends are indicative of future outcomes.
Simplified Assumptions:
Fixed withdrawal rates and contribution percentages may not reflect real-life variability.
Employee preferences and external factors (e.g., tax changes) were not modeled.
Market Risks:
- Simulations cannot predict black swan events (e.g., 2008 financial crisis).
Conclusion
The choice between TRS and ORP depends on the employee’s individual priorities and circumstances. For those seeking stability and simplicity, TRS is a safer option. For those aiming to maximize retirement wealth and are comfortable with risk, ORP can offer higher potential returns, provided they follow a disciplined withdrawal strategy.
Final Recommendation: Align the decision with personal goals, risk tolerance, and financial needs. Regularly review and adjust the chosen plan as circumstances evolve.