• 设为首页
  • 点击收藏
  • 手机版
    手机扫一扫访问
    迪恩网络手机版
  • 关注官方公众号
    微信扫一扫关注
    公众号

拓端tecdat|R语言代写估计时变VAR模型时间序列的实证研究分析案例 ...

原作者: [db:作者] 来自: [db:来源] 收藏 邀请

 

原文 http://tecdat.cn/?p=3364

 



上述症状数据集包含在R-package  中,并在加载时自动可用。 加载包后,我们将此数据集中包含的12个心情变量进行子集化:

mood_data <- as.matrix(symptom_data$data[, 1:12]) # Subset variables
mood_labels <- symptom_data$colnames[1:12] # Subset variable labels
colnames(mood_data) <- mood_labels
time_data <- symptom_data$data_time

 

对象mood_data是一个1476×12矩阵,测量了12个心情变量:

> dim(mood_data)
[1] 1476 12
> head(mood_data[,1:7])
Relaxed Down Irritated Satisfied Lonely Anxious Enthusiastic
[1,] 5 -1 1 5 -1 -1 4
[2,] 4 0 3 3 0 0 3
[3,] 4 0 2 3 0 0 4
[4,] 4 0 1 4 0 0 4
[5,] 4 0 2 4 0 0 4
[6,] 5 0 1 4 0 0 3

 

 

time_data包含有关每次测量的时间戳的信息。数据预处理需要此信息。

> head(time_data)
date dayno beepno beeptime resptime_s resptime_e time_norm
1 13/08/12 226 1 08:58 08:58:56 09:00:15 0.000000000
2 14/08/12 227 5 14:32 14:32:09 14:33:25 0.005164874
3 14/08/12 227 6 16:17 16:17:13 16:23:16 0.005470574
4 14/08/12 227 8 18:04 18:04:10 18:06:29 0.005782097
5 14/08/12 227 9 20:57 20:58:23 21:00:18 0.006285774
6 14/08/12 227 10 21:54 21:54:15 21:56:05 0.006451726

 

该数据集中的一些变量是高度偏斜的,这可能导致不可靠的参数估计。 在这里,我们通过计算自举置信区间(KS方法)和可信区间(GAM方法)来处理这个问题,以判断估计的可靠性。 由于本教程的重点是估计时变VAR模型,因此我们不会详细研究变量的偏度。 然而,在实践中,应该在拟合(时变)VAR模型之前始终检查边际分布。




估计时变VAR模型



通过参数lags = 1,我们指定拟合滞后1 VAR模型,并通过lambdaSel =“CV”选择具有交叉验证的参数λ。 最后,使用参数scale = TRUE,我们指定在模型拟合之前,所有变量都应缩放为零和标准差1。 当使用“1正则化”时,建议这样做,因为否则参数惩罚的强度取决于预测变量的方差。 由于交叉验证方案使用随机抽取来定义折叠,因此我们设置种子以确保重现性。

在查看结果之前,我们检查了1476个时间点中有多少用于估算,这在调用控制台中的输出对象时打印的摘要中显示

> tvvar_obj
mgm fit-object
Model class: Time-varying mixed Vector Autoregressive (tv-mVAR) model
Lags: 1
Rows included in VAR design matrix: 876 / 1475 ( 59.39 %)
Nodes: 12
Estimation points: 20

 

估计的VAR系数的绝对值存储在对象tvvar_obj $ wadj中,该对象是维度p×p×滞后×estpoints的数组。


参数估计的可靠性

res_obj <- resample(object = tvvar_obj,
data = mood_data,
nB = 50,
blocks = 10,seeds = 1:50,
quantiles = c(.05, .95))

 

res_obj $ bootParameters包含每个参数的经验采样分布。



计算时变预测误差



函数predict()计算给定mgm模型对象的预测和预测误差。 

 

预测存储在pred_obj $预测中,并且所有时变模型的预测误差组合在pred_obj中:

> pred_obj$errors
Variable Error.RMSE Error.R2
1 Relaxed 0.939 0.155
2 Down 0.825 0.297
3 Irritated 0.942 0.119
4 Satisfied 0.879 0.201
5 Lonely 0.921 0.182
6 Anxious 0.950 0.086
7 Enthusiastic 0.922 0.169
8 Suspicious 0.818 0.247
9 Cheerful 0.889 0.200
10 Guilty 0.928 0.175
11 Doubt 0.871 0.268
12 Strong 0.896 0.195

 

可视化时变VAR模型

可视化上面估计的一部分随时间变化的VAR参数:

# Two Network Plots

# Get layout of mean graph
Q <- qgraph(t(mean_wadj), DoNotPlot=TRUE)
saveRDS(Q$layout, "Tutorials/files/layout_mgm.RDS")

# Plot graph at selected fixed time points
tpSelect <- c(2, 10, 18)

# Switch to colorblind scheme
tvvar_obj$edgecolor[, , , ][tvvar_obj$edgecolor[, , , ] == "darkgreen"] <- c("darkblue")
lty_array <- array(1, dim=c(12, 12, 1, 20))
lty_array[tvvar_obj$edgecolor[, , , ] != "darkblue"] <- 2

for(tp in tpSelect) {
  qgraph(t(tvvar_obj$wadj[, , 1, tp]), 
         layout = Q$layout,
         edge.color = t(tvvar_obj$edgecolor[, , 1, tp]), 
         labels = mood_labels, 
         vsize = 13, 
         esize = 10,
         asize = 10, 
         mar = rep(5, 4), 
         minimum = 0, 
         maximum = .5, 
         lty = t(lty_array[, , 1, tp]),
         pie = pred_obj$tverrors[[tp]][, 3])
}

CIs <- apply(res_obj$bootParameters[par_row[1], par_row[2], 1, , ], 1, function(x) {
    quantile(x, probs = c(.05, .95))
  } )
  
  # Plot shading
  polygon(x = c(1:20, 20:1), y = c(CIs[1,], rev(CIs[2,])), col=alpha(colour = cols[i], alpha = .3), border=FALSE)
  

  
} # end for: i
 

 

图  显示了上面估计的时变VAR参数的一部分。 顶行显示估计点8,15和18的VAR参数的可视化。蓝色实线箭头表示正关系,红色虚线箭头表示负关系。 箭头的宽度与相应参数的绝对值成比例。

 

 

如果您有任何疑问,请在下面发表评论。 

如果您有任何疑问,请在下面发表评论。 


鲜花

握手

雷人

路过

鸡蛋
该文章已有0人参与评论

请发表评论

全部评论

专题导读
上一篇:
R语言画曲线图发布时间:2022-07-18
下一篇:
R语言的ARIMA模型预测发布时间:2022-07-18
热门推荐
热门话题
阅读排行榜

扫描微信二维码

查看手机版网站

随时了解更新最新资讯

139-2527-9053

在线客服(服务时间 9:00~18:00)

在线QQ客服
地址:深圳市南山区西丽大学城创智工业园
电邮:jeky_zhao#qq.com
移动电话:139-2527-9053

Powered by 互联科技 X3.4© 2001-2213 极客世界.|Sitemap