Stocks SP500 [USA]

chart_us.r

#!/usr/bin/env Rscript
# ═══════════════════════════════════════════════════════════════
# US Equity (S&P 500 TR / VTI) — Log-Trend & Super-Cycle Analysis
# Charts: 1_ Price+Trend, 2_ Deviation+Supertrend, 3_ Normalized
# Data: ERN SWR Toolbox (1871–2025) + Yahoo Finance (recent)
# ═══════════════════════════════════════════════════════════════

library(tidyverse)
library(quantmod)
library(ggrepel)

sf <- 1.5
width <- 600 * sf
height <- 335 * sf
options(vsc.dev.args = list(width = width, height = height, res = 72 * sf))

# ── Configuration ─────────────────────────────────────────────
LABEL  <- "Stocks USA (S&P 500 TR / VTI)"
PREFIX <- "chart_us"
START  <- as.Date("1871-01-01")

# ERN SWR Toolbox — column 6 = S&P 500 Total Return index
ERN_SHEET_ID <- "1QGrMm6XSGWBVLI8I_DOAeJV5whoCnSdmaR8toQB2Jz8"
ERN_COL      <- 6
YAHOO_TICKER <- "SPY"  # adjusted close ≈ total return proxy

RED <- "#C62828"; ORANGE <- "#E65100"; AMBER <- "#F9A825"
GREEN <- "#2E7D32"; BLUE <- "#1565C0"; GRAY <- "#757575"; TEAL <- "#00695C"

# Supertrend
SUPER_CYCLE_YEARS_DEFAULT <- 40
LOESS_SPAN_MULT   <- 1.8
LOESS_MIN_SPAN    <- 0.12
LOESS_MAX_SPAN    <- 0.98
EDGE_DAMPEN_YEARS <- 12
EDGE_MIN_WEIGHT   <- 0.35
Z_BAND            <- 1.96

# ── Helper: parse ERN CSV numbers (comma thousands) ──────────
parse_num <- function(x) suppressWarnings(as.numeric(gsub(",", "", x)))

# ── Download ERN historical data ─────────────────────────────
cat("→ Downloading ERN SWR Toolbox (Asset Returns)...\n")
ern_url <- paste0(
  "https://docs.google.com/spreadsheets/d/", ERN_SHEET_ID,
  "/gviz/tq?tqx=out:csv&sheet=Asset+Returns"
)

ern_raw <- read_csv(ern_url, skip = 2, col_names = FALSE,
                    col_types = cols(.default = "c"), show_col_types = FALSE)

ern_df <- tibble(
  month = parse_num(ern_raw[[1]]),
  year  = parse_num(ern_raw[[2]]),
  close = parse_num(ern_raw[[ERN_COL]])
) |>
  filter(!is.na(month), !is.na(year), !is.na(close), close > 0) |>
  mutate(date = as.Date(sprintf("%04d-%02d-01",
                                as.integer(year), as.integer(month)))) |>
  filter(date >= START) |>
  transmute(date, close) |>
  arrange(date)

cat("  ERN data:", nrow(ern_df), "months,",
    format(min(ern_df$date)), "to", format(max(ern_df$date)), "\n")

# ── Supplement with recent Yahoo Finance data ─────────────────
cat("→ Supplementing with Yahoo Finance (", YAHOO_TICKER, ")...\n")
last_ern <- tail(ern_df, 1)

tryCatch({
  raw <- getSymbols(YAHOO_TICKER, from = last_ern$date - 60,
                    to = Sys.Date(), auto.assign = FALSE)
  yahoo_m <- tibble(
    date  = as.Date(index(raw)),
    close = as.numeric(Ad(raw))
  ) |>
    filter(!is.na(close), close > 0) |>
    mutate(ym = floor_date(date, "month")) |>
    group_by(ym) |> slice_tail(n = 1) |> ungroup() |>
    transmute(date = ym, close) |> arrange(date)

  # Scale Yahoo to ERN at overlap month (try exact, then nearest)
  overlap <- yahoo_m |> filter(date == last_ern$date)
  if (nrow(overlap) == 0) {
    overlap <- yahoo_m |>
      mutate(diff = abs(as.numeric(date - last_ern$date))) |>
      filter(diff <= 45) |> arrange(diff) |> slice(1) |> select(-diff)
  }
  if (nrow(overlap) == 1) {
    scale_f <- last_ern$close / overlap$close
    yahoo_new <- yahoo_m |>
      filter(date > last_ern$date) |>
      mutate(close = close * scale_f)
    if (nrow(yahoo_new) > 0) {
      ern_df <- bind_rows(ern_df, yahoo_new)
      cat("  Added", nrow(yahoo_new), "months from Yahoo.\n")
    }
  }
}, error = \(e) cat("  Yahoo supplement failed:", e$message, "\n"))

df <- ern_df |> arrange(date)
cat("  Final:", nrow(df), "months,",
    format(min(df$date)), "to", format(max(df$date)), "\n")

# ── Log-linear trend ─────────────────────────────────────────
df <- df |> mutate(date_num = as.numeric(date), log_close = log(close))
fit <- lm(log_close ~ date_num, data = df)
cagr <- (exp(coef(fit)[2] * 365.25) - 1) * 100

df <- df |> mutate(
  trend   = exp(predict(fit, newdata = pick(everything()))),
  pct_dev = (close - trend) / trend * 100
)
cat("  Trend CAGR:", round(cagr, 2), "%/yr\n")

# ── Supertrend (edge-adaptive LOESS) ─────────────────────────
data_years <- as.numeric(difftime(max(df$date), min(df$date),
                                  units = "days")) / 365.25
cycle_years <- if (data_years >= 50) {
  SUPER_CYCLE_YEARS_DEFAULT
} else {
  max(8, round(data_years * 0.4))
}
k_months <- as.integer(cycle_years * 12)
n <- nrow(df)
span <- min(LOESS_MAX_SPAN, max(LOESS_MIN_SPAN,
                                LOESS_SPAN_MULT * k_months / n))

t_idx <- seq_len(n)
edge_n <- min(as.integer(EDGE_DAMPEN_YEARS * 12), floor(n / 2))
w <- rep(1, n)
if (edge_n >= 2) {
  ramp <- seq(EDGE_MIN_WEIGHT, 1, length.out = edge_n)
  w[seq_len(edge_n)] <- ramp
  w[(n - edge_n + 1):n] <- rev(ramp)
}

lo <- loess(df$pct_dev ~ t_idx, span = span, degree = 1,
            family = "symmetric", weights = w,
            control = loess.control(surface = "direct"))
df$supertrend <- as.numeric(predict(lo, data.frame(t_idx = t_idx)))
df$norm_dev   <- df$pct_dev - df$supertrend

cat("  Supertrend:", cycle_years, "yr cycle, span =", round(span, 3), "\n")

# ── Label helper ──────────────────────────────────────────────
decade_labels <- function(dates, values, threshold = 20) {
  tbl <- tibble(date = dates, val = values) |>
    mutate(decade = floor(year(date) / 10) * 10) |>
    group_by(decade) |>
    summarize(
      max_d = date[which.max(val)], max_v = max(val),
      min_d = date[which.min(val)], min_v = min(val),
      .groups = "drop"
    )
  bind_rows(
    tbl |> filter(abs(max_v) >= threshold) |>
      transmute(date = max_d, val = max_v),
    tbl |> filter(abs(min_v) >= threshold) |>
      transmute(date = min_d, val = min_v)
  ) |> distinct(date, .keep_all = TRUE)
}

last_str <- format(max(df$date), "%b %Y")
cur <- tail(df, 1)

# ═══════════════════════════════════════════════════════════════
# Chart 1_: Price (log scale) with Log-Trend
# ═══════════════════════════════════════════════════════════════
cat("→ Chart 1: Price + Trend\n")

extrema1 <- bind_rows(
  decade_labels(df$date, df$pct_dev, 20),
  tibble(date = cur$date, val = cur$pct_dev)
) |> distinct(date, .keep_all = TRUE) |>
  left_join(df |> select(date, close), by = "date")

p1 <- ggplot(df, aes(x = date)) +
  # Green: below trend (undervalued)
  geom_ribbon(aes(ymin = pmin(close, trend), ymax = trend),
              fill = GREEN, alpha = 0.15) +
  # Amber: above trend up to 60% over
  geom_ribbon(aes(ymin = trend,
                  ymax = pmin(pmax(close, trend), trend * 1.6)),
              fill = AMBER, alpha = 0.12) +
  # Red: above 60% over trend
  geom_ribbon(aes(ymin = trend * 1.6,
                  ymax = pmax(close, trend * 1.6)),
              fill = RED, alpha = 0.18) +
  geom_line(aes(y = close), color = BLUE, linewidth = 0.5) +
  geom_line(aes(y = trend), color = "black", linetype = "dashed",
            linewidth = 0.5) +
  geom_point(data = extrema1, aes(x = date, y = close),
             color = ORANGE, size = 2) +
  geom_text_repel(data = extrema1,
    aes(x = date, y = close, label = sprintf("%+.0f%%", val)),
    color = BLUE, size = 2.8, nudge_y = 0.15,
    max.overlaps = 15, seed = 42) +
  scale_x_date(date_breaks = "10 year", date_labels = "%Y") +
  scale_y_continuous(trans = "log2",
    labels = scales::label_number(scale_cut = scales::cut_short_scale())) +
  labs(
    title = paste0(LABEL, " (Log Scale)"),
    subtitle = paste0("Since ", format(min(df$date), "%Y"),
      " | Trend CAGR ~ ", round(cagr, 1), "%/yr"),
    x = NULL, y = "Index Level",
    caption = "Source: ERN SWR Toolbox + Yahoo Finance"
  ) +
  theme_bw() +
  theme(axis.text.x = element_text(angle = 30, hjust = 1),
        legend.position = "none")

ggsave(paste0(PREFIX, "_1_price.png"), p1,
       width = width / 72 / sf, height = height / 72 / sf, dpi = 72 * sf)

# ═══════════════════════════════════════════════════════════════
# Chart 2_: Deviation from Log-Trend + Supertrend
# ═══════════════════════════════════════════════════════════════
cat("→ Chart 2: Deviation + Supertrend\n")

lbl2 <- bind_rows(
  decade_labels(df$date, df$pct_dev, 20),
  tibble(date = cur$date, val = cur$pct_dev)
) |> distinct(date, .keep_all = TRUE)

p2 <- ggplot(df, aes(date, pct_dev)) +
  geom_ribbon(aes(ymin = pmin(pct_dev, 0), ymax = 0),
              fill = GREEN, alpha = 0.15) +
  geom_ribbon(aes(ymin = 0, ymax = pmin(pmax(pct_dev, 0), 60)),
              fill = AMBER, alpha = 0.12) +
  geom_ribbon(aes(ymin = 60, ymax = pmax(pct_dev, 60)),
              fill = RED, alpha = 0.18) +
  geom_hline(yintercept = 0, linetype = "dashed", linewidth = 0.4) +
  geom_line(color = BLUE, linewidth = 0.4, alpha = 0.6) +
  geom_line(aes(y = supertrend), color = TEAL, linewidth = 0.8) +
  geom_point(data = lbl2, aes(y = val), color = ORANGE, size = 1.8) +
  geom_text_repel(
    data = lbl2 |> mutate(
      lbl = if_else(date == max(date),
        paste0(last_str, ": ", sprintf("%+.0f%%", val)),
        sprintf("%+.0f%%", val))),
    aes(y = val, label = lbl),
    color = BLUE, size = 2.8, fontface = "bold",
    nudge_y = 5, box.padding = 0.5, max.overlaps = 15, seed = 42
  ) +
  scale_x_date(date_breaks = "10 year", date_labels = "%Y") +
  scale_y_continuous(labels = \(x) paste0(x, "%"),
                     breaks = seq(-60, 120, 20)) +
  labs(
    title = paste0(LABEL, ": Deviation from Log-Trend + Supertrend"),
    subtitle = paste0(
      "Trend CAGR ~ ", round(cagr, 1), "%/yr  |  Teal = ",
      cycle_years, "yr supertrend (edge-adaptive LOESS)"),
    x = NULL, y = "Deviation from trend (%)",
    caption = "Source: ERN SWR Toolbox + Yahoo Finance. Supertrend: robust LOESS with edge down-weighting."
  ) +
  theme_bw() +
  theme(axis.text.x = element_text(angle = 30, hjust = 1),
        legend.position = "none")

ggsave(paste0(PREFIX, "_2_deviation.png"), p2,
       width = width / 72 / sf, height = height / 72 / sf, dpi = 72 * sf)

# ═══════════════════════════════════════════════════════════════
# Chart 3_: Double-Normalized (deviation minus supertrend)
# ═══════════════════════════════════════════════════════════════
cat("→ Chart 3: Double-Normalized\n")

norm_sd <- sd(df$norm_dev, na.rm = TRUE)
band <- Z_BAND * norm_sd

lbl3 <- bind_rows(
  decade_labels(df$date, df$norm_dev, 15),
  tibble(date = cur$date, val = cur$norm_dev)
) |> distinct(date, .keep_all = TRUE)

p3 <- ggplot(df, aes(date, norm_dev)) +
  geom_ribbon(aes(ymin = pmin(norm_dev, 0), ymax = 0),
              fill = GREEN, alpha = 0.15) +
  geom_ribbon(aes(ymin = 0, ymax = pmin(pmax(norm_dev, 0), band)),
              fill = AMBER, alpha = 0.12) +
  geom_ribbon(aes(ymin = band, ymax = pmax(norm_dev, band)),
              fill = RED, alpha = 0.18) +
  geom_hline(yintercept = 0, linetype = "dashed", linewidth = 0.4) +
  geom_hline(yintercept = c(-band, band), color = GRAY,
             linetype = "dotted", linewidth = 0.35) +
  geom_line(color = BLUE, linewidth = 0.5) +
  geom_point(data = lbl3, aes(y = val), color = ORANGE, size = 1.8) +
  geom_text_repel(
    data = lbl3 |> mutate(
      lbl = if_else(date == max(date),
        paste0(last_str, ": ", sprintf("%+.0f%%", val)),
        sprintf("%+.0f%%", val))),
    aes(y = val, label = lbl),
    color = BLUE, size = 2.8, fontface = "bold",
    nudge_y = 3, box.padding = 0.4, max.overlaps = 20, seed = 42
  ) +
  scale_x_date(date_breaks = "5 year", date_labels = "%Y") +
  scale_y_continuous(labels = \(x) paste0(x, "%")) +
  labs(
    title = paste0(LABEL, ": Double-Normalized (Log + Supertrend)"),
    subtitle = paste0(
      "Deviation minus ", cycle_years,
      "yr supertrend  |  Dotted = +/-", Z_BAND, " SD"),
    x = NULL, y = "Normalized deviation (%)",
    caption = "Source: ERN SWR Toolbox + Yahoo Finance. Cyclical over/undervaluation after removing secular trend."
  ) +
  theme_bw() +
  theme(axis.text.x = element_text(angle = 30, hjust = 1),
        legend.position = "none")

ggsave(paste0(PREFIX, "_3_normalized.png"), p3,
       width = width / 72 / sf, height = height / 72 / sf, dpi = 72 * sf)

cat("\n✓ Saved:", paste0(PREFIX, c("_1_price.png", "_2_deviation.png",
                                   "_3_normalized.png")), "\n")

chart_usa_sp500.r

#!/usr/bin/env Rscript
# ═══════════════════════════════════════════════════════════════
# S&P 500 Log-Trend Analysis
# Original: log-scale with trend line
# New: normalized deviation from trend (flat = on trend)
# ═══════════════════════════════════════════════════════════════

library(tidyverse)
library(tsibble)
library(ggrepel)
library(quantmod)

sf <- 1.5
width <- 600 * sf
height <- 335 * sf
options(vsc.dev.args = list(width = width, height = height, res = 72 * sf))

sdate <- as.Date("1900-01-01")
edate <- Sys.Date()

# ── Helper functions ──────────────────────────────────────────
get_symbol_data <- function(symbol, close_col) {
  data <- getSymbols(symbol, from = sdate, to = edate, auto.assign = FALSE) |>
    as.data.frame() |>
    tibble::rownames_to_column("date") |>
    mutate(date = as.Date(date)) |>
    select(date, close = !!sym(close_col)) |>
    as_tsibble(index = date)
  return(data)
}

find_extrema <- function(df) {
  df <- df |>
    mutate(year = year(date), decade = floor(year / 10) * 10)
  maxima <- df |> group_by(decade) |> slice(which.max(percent_diff)) |> ungroup()
  minima <- df |> group_by(decade) |> slice(which.min(percent_diff)) |> ungroup()
  bind_rows(maxima, minima) |> filter(abs(percent_diff) >= 10)
}

calculate_cagr <- function(data) {
  num_years <- as.numeric(difftime(max(data$date), min(data$date), units = "days")) / 365.25
  (last(data$close) / first(data$close))^(1 / num_years) - 1
}

# ── Chart 1: Original log-scale with trend line ──────────────
generate_chart <- function(df, title, subtitle, file_name, x_break = "5 year") {
  df <- df |>
    mutate(date_num = as.numeric(date), log_close = log2(close))

  lm_model <- lm(log_close ~ date_num, data = df)

  df <- df |>
    mutate(
      log_predicted = predict(lm_model, newdata = df),
      predicted_close = 2^(log_predicted),
      percent_diff = (close - predicted_close) / predicted_close * 100
    )

  max_abs <- max(abs(df$percent_diff), na.rm = TRUE)
  percent_band_edges <- seq(-max_abs, max_abs, length.out = 7)
  labels <- c(
    "Negative High", "Negative Medium", "Negative Low",
    "Positive Low", "Positive Medium", "Positive High"
  )

  df <- df |>
    mutate(band = cut(percent_diff, breaks = percent_band_edges,
                      labels = labels, right = FALSE, include.lowest = TRUE))

  df_ribbon <- df |>
    arrange(date) |>
    mutate(
      ymin = pmin(close, predicted_close),
      ymax = pmax(close, predicted_close),
      band_change = (band != lag(band)) | (row_number() == 1),
      group_id = cumsum(ifelse(is.na(band_change), 0, band_change))
    ) |> filter(!is.na(band))

  band_colors <- c(
    "Negative High" = "red", "Negative Medium" = "yellow",
    "Negative Low" = "green", "Positive Low" = "green",
    "Positive Medium" = "yellow", "Positive High" = "red"
  )

  extrema <- bind_rows(
    find_extrema(df) |> as_tibble(),
    df |> filter(date == max(date)) |> as_tibble()
  ) |> distinct(date, .keep_all = TRUE) |>
    as_tsibble(index = date) |>
    mutate(label = sprintf("%.0f%%", percent_diff))

  chart <- ggplot() +
    geom_ribbon(data = df_ribbon,
      aes(x = date, ymin = ymin, ymax = ymax, fill = band, group = group_id),
      alpha = 0.3) +
    scale_fill_manual(values = band_colors) +
    geom_line(data = df, aes(x = date, y = close),
              color = "#5383EC", linewidth = 1) +
    geom_line(data = df, aes(x = date, y = predicted_close),
              color = "black", linetype = "dashed") +
    geom_point(data = extrema, aes(x = date, y = close),
               color = "orange", size = 2) +
    geom_text_repel(data = extrema,
      aes(x = date, y = close, label = label),
      nudge_y = 0.15, color = "blue", size = 3) +
    scale_x_date(date_breaks = x_break, date_labels = "%Y") +
    scale_y_continuous(trans = "log2") +
    labs(title = title, subtitle = subtitle, x = "Date", y = "Points") +
    theme_bw() +
    theme(axis.text.x = element_text(angle = 30, hjust = 0.5, vjust = 0.5),
          legend.position = "none")

  return(chart)
}

# ── Chart 2: Normalized deviation from trend (flat) ──────────
generate_deviation_chart <- function(df, title, subtitle, x_break = "10 year") {
  df <- as_tibble(df) |>
    mutate(date_num = as.numeric(date), log_close = log(close))

  lm_model <- lm(log_close ~ date_num, data = df)
  trend_slope <- coef(lm_model)[2]
  trend_cagr <- (exp(trend_slope * 365.25) - 1) * 100

  df <- df |>
    mutate(
      predicted = exp(predict(lm_model, newdata = df)),
      pct_deviation = (close - predicted) / predicted * 100
    )

  # Decade extrema labels (threshold >=20%)
  decade_extrema <- df |>
    mutate(decade = floor(year(date) / 10) * 10) |>
    group_by(decade) |>
    summarize(
      max_date = date[which.max(pct_deviation)],
      max_pct  = max(pct_deviation),
      min_date = date[which.min(pct_deviation)],
      min_pct  = min(pct_deviation),
      .groups  = "drop"
    )

  label_points <- bind_rows(
    decade_extrema |> filter(abs(max_pct) >= 20) |>
      transmute(date = max_date, pct = max_pct),
    decade_extrema |> filter(abs(min_pct) >= 20) |>
      transmute(date = min_date, pct = min_pct)
  ) |> distinct(date, .keep_all = TRUE)

  # Add current point
  current <- tail(df, 1)
  label_points <- bind_rows(
    label_points,
    tibble(date = current$date, pct = current$pct_deviation)
  ) |> distinct(date, .keep_all = TRUE)

  RED <- "#C62828"; ORANGE <- "#E65100"; AMBER <- "#F9A825"
  GREEN <- "#2E7D32"; BLUE <- "#1565C0"

  chart <- df |>
    ggplot(aes(x = date, y = pct_deviation)) +
    geom_ribbon(aes(ymin = pmin(pct_deviation, 0), ymax = 0),
                fill = GREEN, alpha = 0.15) +
    geom_ribbon(aes(ymin = 0, ymax = pmin(pmax(pct_deviation, 0), 60)),
                fill = AMBER, alpha = 0.12) +
    geom_ribbon(aes(ymin = 60, ymax = pmax(pct_deviation, 60)),
                fill = RED, alpha = 0.18) +
    geom_hline(yintercept = 0, color = "black", linetype = "dashed",
               linewidth = 0.4) +
    geom_line(color = BLUE, linewidth = 0.5) +
    geom_point(data = label_points, aes(x = date, y = pct),
               color = ORANGE, size = 1.8) +
    geom_text_repel(data = label_points,
      aes(x = date, y = pct, label = sprintf("%+.0f%%", pct)),
      color = BLUE, size = 2.8, fontface = "bold",
      nudge_y = 5, min.segment.length = 0.3,
      box.padding = 0.4, point.padding = 0.3,
      max.overlaps = 20, seed = 42) +
    scale_x_date(date_breaks = x_break, date_labels = "%Y") +
    scale_y_continuous(labels = \(x) paste0(x, "%"),
                       breaks = seq(-60, 120, 20)) +
    labs(
      title = title,
      subtitle = paste0(
        "Trend CAGR ≈ ", round(trend_cagr, 1),
        "%/yr  |  0% = on trend  |  Deviations are mean-reverting"
      ),
      x = NULL, y = "Deviation from trend (%)",
      caption = "Source: Yahoo Finance (^GSPC)"
    ) +
    theme_bw() +
    theme(axis.text.x = element_text(angle = 30, hjust = 1))

  return(chart)
}

# ═══════════════════════════════════════════════════════════════
# Run
# ═══════════════════════════════════════════════════════════════
sp500_df <- get_symbol_data(symbol = "^GSPC", close_col = "GSPC.Close")
cat("S&P 500 CAGR:", round(calculate_cagr(sp500_df) * 100, 2), "%/yr\n")

sp500_post45 <- sp500_df |> filter(year(date) > 1945)

# Chart 1: Log-scale with trend line (original)
chart1 <- generate_chart(
  df = sp500_post45,
  title = "S&P 500 (Log Scale)",
  subtitle = "Since 1945 · Source: ^GSPC",
  x_break = "5 year"
)

ggsave("chart1.png", chart1,
       width = width / 72, height = height / 72,
       units = "in", dpi = 72 * sf)

# Chart 2: Normalized deviation from trend (flat)
chart2 <- generate_deviation_chart(
  df = sp500_post45,
  title = "S&P 500: Deviation from Log-Linear Trend",
  subtitle = "Since 1945",
  x_break = "10 year"
)

ggsave("chart2.png", chart2,
       width = width / 72, height = height / 72,
       units = "in", dpi = 72 * sf)

# Chart 3: Total Return (if available)
tryCatch({
  sp500tr_df <- get_symbol_data("^SP500TR", "SP500TR.Close")
  cat("S&P 500 TR CAGR:", round(calculate_cagr(sp500tr_df) * 100, 2), "%/yr\n")

  chart3 <- generate_chart(
    df = sp500tr_df,
    title = "S&P 500 Total Return (Log Scale)",
    subtitle = "Source: ^SP500TR",
    x_break = "2 year"
  )

  ggsave("chart3.png", chart3,
         width = width / 72, height = height / 72,
         units = "in", dpi = 72 * sf)

  chart4 <- generate_deviation_chart(
    df = sp500tr_df,
    title = "S&P 500 Total Return: Deviation from Log-Linear Trend",
    subtitle = "Source: ^SP500TR"
  )

  ggsave("chart4.png", chart4,
         width = width / 72, height = height / 72,
         units = "in", dpi = 72 * sf)
}, error = \(e) cat("SP500TR not available:", e$message, "\n"))

https://www.mortality.watch/charts/list.html#stocks-sp500-usa