如何用R语言hmmTMB包对观测状态做一步超前预测？

问题背景 #

基于训练好的hmmTMB隐马尔可夫模型，对测试集的油价进行一步超前预测，核心是利用模型的状态转移规律和观测分布参数，结合滤波后的状态概率生成预测值。

完整代码示例（含预测部分） #

# Load necessary libraries
library(hmmTMB)
library(tidyverse)

# Set seed for reproducibility
set.seed(123)

# Example dataset preparation (extended for 1 year)
df <- data.frame(
    day = seq.Date(from = as.Date("2020-01-01"), to = as.Date("2020-12-31"), by = "day"),
    OILPRICE = rnorm(366, 100, 10)
    ) %>%
    as_tibble()

# Splitting into training and test sets
split_index <- round(nrow(df) * 0.8)
training_set <- df[1:split_index, ] %>% as.data.frame()
testing_set <- df[(split_index + 1):nrow(df), ] %>% as.data.frame()

# Define a hidden Markov model with 2 states for the training set
hid1 <- MarkovChain$new(data = training_set, n_states = 2)
dists <- list(OILPRICE = "norm")
par0 <- list(OILPRICE = list(mean = c(110, 90), sd = c(1, 1)))
obs1 <- Observation$new(data = training_set, n_states = 2, dists = dists, par = par0)
par0 <- obs1$suggest_initial()
obs1 <- Observation$new(data = training_set, n_states = 2, dists = dists, par = par0)
hmm1 <- HMM$new(obs = obs1, hid = hid1)
hmm1$fit(silent = TRUE)

# Rename and prepare visualization data
data_plot <- training_set
pal <- c("#1f77b4", "#ff7f0e") # 定义状态颜色
data_plot$state <- factor(paste0("State ", hmm1$viterbi()))

# Plot training data with Viterbi states
ggplot(data_plot, aes(day, OILPRICE, col = state)) +
    geom_point() +
    scale_color_manual(values = pal, name = NULL) +
    theme_minimal()

# --------------------------
# 方法1：手动计算一步超前预测
# --------------------------
# 获取训练好的模型参数
trans_mat <- hmm1$hid$par$transition # 状态转移矩阵
obs_params <- hmm1$obs$par$OILPRICE # 观测分布的均值和标准差
mean_vec <- obs_params$mean
sd_vec <- obs_params$sd

# 获取训练集最后时刻的滤波状态概率（作为预测的初始状态）
train_filter <- hmm1$filter()
initial_state_probs <- train_filter$alpha[nrow(train_filter$alpha), ]

# 初始化存储变量
testing_set$predicted_price_manual <- NA
current_state_probs <- initial_state_probs

# 逐步计算预测值
for (i in 1:nrow(testing_set)) {
  # 预测下一时刻的状态概率
  next_state_probs <- current_state_probs %*% trans_mat
  # 计算油价预测值（基于状态概率的加权期望）
  predicted_val <- sum(next_state_probs * mean_vec)
  testing_set$predicted_price_manual[i] <- predicted_val
  
  # 更新状态概率：用当前观测值滤波
  obs_probs <- dnorm(testing_set$OILPRICE[i], mean = mean_vec, sd = sd_vec)
  updated_probs <- next_state_probs * obs_probs
  current_state_probs <- updated_probs / sum(updated_probs)
}

# --------------------------
# 方法2：使用hmmTMB内置predict函数
# --------------------------
# 合并训练与测试数据，让模型从训练数据结尾开始预测
full_data <- rbind(training_set, testing_set)
# 生成一步超前预测
predictions <- predict(hmm1, newdata = full_data, type = "response", one_step_ahead = TRUE)
# 提取测试集对应的预测值
testing_set$predicted_price_auto <- predictions$OILPRICE[(split_index + 1):nrow(full_data)]

# --------------------------
# 可视化预测结果
# --------------------------
ggplot(testing_set, aes(x = day)) +
  geom_point(aes(y = OILPRICE, color = "实际油价"), size = 2) +
  geom_point(aes(y = predicted_price_manual, color = "手动预测"), size = 2, shape = 17) +
  geom_point(aes(y = predicted_price_auto, color = "内置函数预测"), size = 2, shape = 15) +
  scale_color_manual(values = c("实际油价" = "black", "手动预测" = "#ff7f0e", "内置函数预测" = "#1f77b4")) +
  labs(title = "油价一步超前预测对比", x = "日期", y = "油价", color = NULL) +
  theme_minimal()

关键说明 #

• 手动计算逻辑：先通过训练集最后一步的滤波状态概率，结合转移矩阵得到下一状态概率，再用状态概率加权观测分布的均值得到预测值，之后用当前观测值更新状态概率，循环处理整个测试集。
• 内置函数预测：predict方法的one_step_ahead=TRUE参数会生成每个时间步的一步超前预测，需要将训练集和测试集合并，确保模型从训练数据的状态延续过来。

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