R语言使用ggplot2绘制时间序列结构断点分段回归趋势图

我尝试使用ggplot2复现下图的图表结构，但尚未找到实现方法：

我持有时间序列数据，需要为单位根检验分析中计算得到的每个结构断点添加对应的回归趋势线。

现有代码如下：

autoplot(merval_usd) +
  labs(x = "Año", y = "Log_Merval en USD")+
  theme_bw()

我的数据为MERVAL价格指数，共4888条观测值，前20条数据结构如下：

> dput(head(merval_usd,20))
structure(c(6.31324002782979, 6.25952410104302, 6.27792086863394, 
6.26998603927151, 6.25789744652059, 6.25377111760352, 6.25553888702416, 
6.25132568135736, 6.32806177591735, 6.33873541097142, 6.34166356744438, 
6.35940084332373, 6.36564536517638, 6.36090516203878, 6.3456363608286, 
6.34242066971575, 6.33416733347983, 6.36629861166262, 6.36406205159053, 
6.34344097026516), class = c("xts", "zoo"), index = structure(c(946857600, 
946944000, 947030400, 947116800, 947203200, 947462400, 947548800, 
947635200, 947721600, 947808000, 948067200, 948153600, 948240000, 
948326400, 948412800, 948672000, 948758400, 948844800, 948931200, 
949017600), tzone = "UTC", tclass = "Date"), .Dim = c(20L, 1L
), .Dimnames = list(NULL, "MERV"))

前10条数据预览：

> head(merval_usd, 10)
           MERV
2000-01-03 6.313240
2000-01-04 6.259524
2000-01-05 6.277921
2000-01-06 6.269986
2000-01-07 6.257897
2000-01-10 6.253771
2000-01-11 6.255539
2000-01-12 6.251326
2000-01-13 6.328062
2000-01-14 6.338735

我使用修改自urca::ur.za的ur.ka函数计算断点。
断点计算与绘图的现有代码如下：

for (i in 1){
  inicio<- Sys.time()
  urka <- ur.ka(merval_usd$MERV, model = 'trend', bp = 5) #bp5 y 12 lag, trend, se rechaza, con 24 rezagos tmbn es "bueno" el ajuste
  #bp5, 24, both
  fin<-Sys.time()
  print(urka)
  print(stringr::str_c("tiempo de ejecucion ", fin-inicio))
}

plot(urka$testreg$residuals)

plot<- ggplot(merval_usd, aes(x = Index, y = MERV))+params
bp<-vector("list", length = length(urka$bpoints))
count<-0
for(i in urka$bpoints){
  print(merval_usd[i])
  count = count+1
  bp[[count]]<- geom_vline(xintercept = as.Date(index(merval_usd$MERV[i])), color = "black", lwd = 0.5, lty=2)
  print(plot+bp)
}

params<- list(geom_line(linetype=1, lwd=1, colour="steelblue"),
              labs(title = "Log Indice Merval - período 2000-2020"),
              theme_minimal())

我目前已生成的图表如下：

我此前使用autoplot做了数据预览，但希望用ggplot2完成自定义绘图，目前已跑完模型，得到了每个断点对应的截距和趋势系数，请问如何实现分段回归趋势线的绘制？

首先将xts格式的时间序列转换为data.frame格式，方便后续分组和绘图：

library(tidyverse)
library(zoo)
# 转换数据格式
merval_df <- as.data.frame(merval_usd) %>%
  rownames_to_column("date") %>%
  mutate(date = as.Date(date),
         time_num = as.numeric(date)) # 数值型时间索引，用于回归计算

方案1：使用已有模型的截距和斜率绘制 #

如果已经从urka结果中提取到每个分段的截距、斜率，以及对应分段的时间范围，直接绘制即可：

# 先提取所有断点的日期
break_dates <- index(merval_usd)[urka$bpoints] %>% as.Date()
# 拼接首尾日期生成分段区间
seg_dates <- c(min(merval_df$date), break_dates, max(merval_df$date))

# 替换下面的数值为你从模型中提取的对应截距、斜率
coef_df <- data.frame(
  seg_id = 1:(length(break_dates)+1),
  intercept = c(你的截距数值列表),
  slope = c(你的斜率数值列表)
)

# 绘图
ggplot(merval_df, aes(x = date, y = MERV)) +
  geom_line(color = "steelblue", linewidth = 1) +
  # 加断点竖线
  geom_vline(xintercept = break_dates, linetype = 2, linewidth = 0.5) +
  # 加每段趋势线
  pmap(coef_df, function(seg_id, intercept, slope) {
    x_start <- seg_dates[seg_id]
    x_end <- seg_dates[seg_id + 1]
    geom_abline(intercept = intercept, slope = slope, 
                xlim = c(x_start, x_end), color = "red", linewidth = 1)
  }) +
  labs(title = "Log Indice Merval - período 2000-2020",
       x = "Año", y = "Log_Merval en USD") +
  theme_minimal()

方案2：直接按分段自动拟合趋势线 #

如果不需要严格复用urka输出的系数，只想按断点分段拟合线性趋势，可以先给数据打分段标签再分组拟合：

# 给数据打分段标签
merval_df <- merval_df %>%
  mutate(seg_group = cut(date, breaks = seg_dates, include.lowest = T, labels = F))

# 绘图
ggplot(merval_df, aes(x = date, y = MERV)) +
  geom_line(color = "steelblue", linewidth = 1) +
  geom_vline(xintercept = break_dates, linetype = 2, linewidth = 0.5) +
  # 按分组拟合线性趋势，关闭置信区间
  geom_smooth(aes(group = seg_group), method = "lm", se = F, color = "red", linewidth = 1) +
  labs(title = "Log Indice Merval - período 2000-2020",
       x = "Año", y = "Log_Merval en USD") +
  theme_minimal()

内容的提问来源于stack exchange，提问作者Lucas Guzman