可以使用
function  na.approx(object, ...) or
function  na.spline(object, ...)
  http://rss.acs.unt.edu/Rdoc/library/zoo/html/na.approx.html
不过顾名思义(Interpolation),你的数据格式似m[,1],是会产生错误的,
要类似m[,4],才能插值.
example:
  set.seed(12)
  m <- matrix(runif(16, 0, 100), nrow = 4)
  missing_values <- sample(16, 7)
  m[missing_values] <- NA
  m
#           [,1]      [,2]     [,3]     [,4]
#[1,]        NA 16.934812       NA 37.62485
#[2,]        NA  3.389562       NA 38.08122
#[3,]        NA 17.878500 39.26972       NA
#[4,]  26.93819        NA 81.38806 43.93343

  na.approx(m[,1])
#Error in approx(x[!na], y[!na], xout, ...) :
#  need at least two non-NA values to interpolate

  na.approx(m[,2])
#[1] 16.934812  3.389562 17.878500
  na.approx(m[,3])
#[1] 39.26972 81.38806
  na.approx(m[,4])
#[1] 37.62485 38.08122 41.00733 43.93343

epoh老师，附件是一篇用贝叶斯估计多变量随机波动的，我想请教下，为何Fig 3和Fig 5的很多密度图都不成正态分布，
作者也没有对这个进行解释，不是说一般要成正态分布吗？

谢谢epoh老师。

先了解下,你为何要用fat tail SV (t-MSV)
MULTIVARIATE VOLATILITY MODELS.pdf
page 17/20已经告诉你
Form the financial point of view t-MSV models are less likely to provide a good fit

麻烦你看一下数据对不对
ywh_msv.txt   

epoh老师，恩，您给出的这个数据是正确的。
我想问下62楼给出的文献，很多密度图不是正态
分布的，这样可以吗？还是有什么其他标准判别？
因为作者对图形解释都没有提到。

非常感谢epoh老师，恩，现在主要是用这个模型且考虑two-way granger
causality test。

通常tau的图形比较不会 smooth
我没用fat-tailed(t-msv).
你的数据跑出来的图形如下,
你先参考看看有啥问题没?
代码我再检查一下,晚点给你.
  

好的，非常感谢epoh老师
您这个用的是文章哪个模型？
还有，可以怎么解释tau一般比较不平滑。

1.c:/Bugs/msv/作为working.directory
  放进两文件 msvmodel6.bug
             ywh_msv.txt(your data)
2.注意你的winbugs安装路径,若不是"D:/WinBUGS14/"
  请依你所安装修改.
#####
3.in R 运行底下程序:
library(R2WinBUGS)
#Data
n=618
yy=read.table(file="c:/Bugs/msv/ywh_msv.txt")
...
...


密度图形解释的问题,
我找找看有哪篇文献解释到.
请注意短信息

好的，epoh老师。
谢谢您的帮忙。

epoh老师，利用R运行Winbugs可以得出核密度图吗？
现在默认的好像只有轨迹图。
ps：利用R运行，方便很多，谢谢。

核密度图,需要另外编程.
我再看一下文献的图后,再给你代码.

好的，非常感谢epoh老师。

好像是我弄混淆了
你是要这篇文献的Model 6 (Heavy-tailed MSV or t-MSV)??
   http://harrisd.net/papers/ARCHSV ... s/YuMeyer2006ER.pdf[/url]
如果是的话,我再加入w好了.

作者的computer太好了,羡慕!!
The CPU time required for 100 iterations ranges from 2.3 seconds to 14.4 seconds.

谢谢epoh老师
是的，要的就是Model 6，但是同时像Model 3一样考虑Granger causality。
即
即假设上述公式的非对角线元素也非0，这样就可以做双向因果检验了。
我看了您修改的程序，好像也是这样的。
哈哈，我也羡慕。
还有请epoh老师看能不能找到对于核密度图不成正态的文章解释。

t-MSV with two-way granger causality:
请注意这次我改同MULTIVARIATE STOCHASTIC VOLATILITY MODELS
page 15/25的试验条件.
  burn-in period of 10,000 iterations and a follow-up period of 100,000,
  we stored every 20th iteration.
msv.sim <- bugs(data, inits, parameters.to.save=parameters,"msvmodel6_t.bug",n.chains=1,
               n.thin=20,n.iter=110000,n.burnin=10000,debug=TRUE,DIC=TRUE,
               bugs.directory="D:/WinBUGS14/",working.directory = "c:/Bugs/msv/")

############
density plot可以利用package coda function densplot(),或plot.mcmc()
#coda.out1 <- read.coda("c:/Bugs/msv/coda1.txt", "c:/Bugs/msv/codaIndex.txt")
#densplot(coda.out1)
#plot(coda.out1)
但图形比较乱,我改写代码如下:
M = 5000
chain1 = matrix(scan("c:/Bugs/msv/coda1.txt"),11*M,2,byrow=T)
chain1 = matrix(chain1[,2],M,11)
chain1=apply(chain1,2,sort)

####density.plot
names = c("d","deviance", "mu1","mu2","phi1","phi2","phi12","phi21","rhoep","taua","taub")
par(mfrow=c(3,4))
for (i in 1:11)
{
plot(density(chain1[,i]),main="")
title(main=names)
}

###########
核密度图不成正态,应该不是重点
请看模型检验:
Oil and stock market volatility.pdf  page 7/10
4.5. Model checking
     We observe that the mean of standardized residual in each market
     is close to zero while its standard deviation is close to one.

谢谢epoh老师
恩，看了那句话，意思是指标准化残差的均值为0，标准差
为1，那说明是服从正态分布。
我没有理解，那和参数的核密度图不成正态分布之间关系。
不过这一段的最后，作者提到了
However,it is clear that the non-normality is
less pronounced in the standardized residuals.Thus,both CC-MSV and
DC-MSV models do a fine job in modeling heteroscedasticity in the
data.

这个应该就是类似
Stochastic Volatility and Jumps...pdf
page 29/42,Normal Probability Plot (QQ-plot)
https://bbs.pinggu.org/thread-2126597-1-1.html

是的，Q-Q 图也是一种判别残差是不是服从
正态分布的方法。epoh老师，您修改的程序不是很理解

您修改的程序多定义了一个w，即服从卡方分布
d<-dstar*2;
dstar~dchisqr(4);  
我不是很理解。在您发给我的文章15 / 25里面提到：
For the last three models,there are two sets of parameters:parameters in
the mean equation,and in the factor equation We assume that the parameters are mutually independent.The prior distributions are specified as

这里提到model 7,8,9才有w这个参数，为何现在利用model 6就
需要这个参数？还有，model 7,8,9中的w出现在哪？我都没有看到
这个变量，只是在设置prior时提到这个参数。

哈哈!你忘了你提供的代码,不也是有w

这是在Heavy-tailed才用上
Model 6 (Heavy-tailed MSV or t-MSV)
后三者也只有Model 8才用的着
Model 8 (Heavy-tailed factor-MSV or AFactor-t-MSV).

抱歉啊，epoh老师，现在明白了。
谢谢您的指导。

epoh老师，我想请教下，在您写的代码中
#coda.out1 <- read.coda("c:/Bugs/msv/coda1.txt", "c:/Bugs/msv/codaIndex.txt")
#densplot(coda.out1)
#plot(coda.out1)

plot(coda.out1)这句命令直接同时画出了trace 和density 图，
这个plot只有在library(R2WinBUGS)之后有这个功能，我想问下
如果我现在想分开来画trace和density，所以我想看这个环境下plot
的代码，应该如何查看？

coda.pdf page 26/44
plot(x, trace = TRUE, density = TRUE, smooth = FALSE, bwf,
     auto.layout = TRUE, ask = dev.interactive(), ...)

#plot trace
   plot(x, trace = TRUE, density = F)
#plot density
   plot(x, trace = F, density = TRUE)

太感谢epoh老师，明白了。

epoh老师，您好！
    1、附件中的数据做dccgarch数据太少，您看在R中如何能实现自体抽样扩大
样本数据量，以比较准确地做出dccgarch估计？
    2、如何修改以下程序，用我附件中的dccdata.txt中的数据进行估计？
非常感谢！

# Simulating data from the original DCC-GARCH(1,1) process
  nobs <- 1000; cut <- 1000
  a <- c(0.003, 0.005, 0.001)
  A <- diag(c(0.2,0.3,0.15))
  B <- diag(c(0.75, 0.6, 0.8))
  uncR <- matrix(c(1.0, 0.4, 0.3, 0.4, 1.0, 0.12, 0.3, 0.12, 1.0),3,3)
  dcc.para <- c(0.01,0.98)
  dcc.data <- dcc.sim(nobs, a, A, B, uncR, dcc.para, model="diagonal")

# Estimating a DCC-GARCH(1,1) model
  dcc.results <- dcc.estimation(inia=a, iniA=A, iniB=B, ini.dcc=dcc.para,
        dvar=dcc.data$eps, model="diagonal")

# Parameter estimates and their robust standard errors
  dcc.results$out

老兄,R package BayesPanel
应该适合你参考.
   http://cran.r-project.org/src/contrib/Archive/BayesPanel/
#####
library(BayesPanel)
library(plm)
data(Gasoline)
est <- BayesPanel(lgaspcar~lincomep+lrpmg+lcarpcap|lincomep,
                     data=Gasoline, index = c("country"),
                     control = list(verbose=10, thin=100) )

Summary Statistics for Posterior Distribution:


The Model Specification is:
BayesPanel(formula = lgaspcar ~ lincomep + lrpmg + lcarpcap |
    lincomep, data = Gasoline, index = c("country"), control = list(verbose = 10,
    thin = 100))


Posterior Estimates of Fixed Effect Coefficients:

Iterations = 1001:10901
Thinning interval = 100
Number of chains = 1
Sample size per chain = 100

1. Empirical mean and standard deviation for each variable,
   plus standard error of the mean:

               Mean      SD Naive SE Time-series SE
(Intercept)  2.1868 0.44051 0.044051       0.044051
lincomep     0.3698 0.11068 0.011068       0.011068
lrpmg       -0.3246 0.04666 0.004666       0.004666
lcarpcap    -0.4641 0.04210 0.004210       0.004210

2. Quantiles for each variable:

               2.5%     25%     50%     75%   97.5%
(Intercept)  1.2418  1.9139  2.2340  2.5018  2.9245
lincomep     0.1246  0.3160  0.3750  0.4464  0.5572
lrpmg       -0.4198 -0.3548 -0.3187 -0.2957 -0.2389
lcarpcap    -0.5451 -0.4930 -0.4646 -0.4355 -0.3859



Posterior Mean Estimates of Random Coefficients:
          (Intercept)     lincomep
AUSTRIA   1.910348096  0.333092578
BELGIUM  -0.215327621 -0.001088742
CANADA    1.152187949  0.090922766
DENMARK  -1.838407500 -0.324946283
FRANCE    0.612314329  0.140129525
GERMANY   1.802215674  0.350771888
GREECE    1.312208989  0.180310358
IRELAND   1.587910877  0.229149729
ITALY     0.004301386  0.040926582
JAPAN    -4.191993903 -0.689216994
NETHERLA -2.400792420 -0.383283336
NORWAY   -0.131194057 -0.004403345
SPAIN    -1.114354594 -0.126608625
SWEDEN    0.258225455  0.032832038
SWITZERL  1.005269013  0.171583757
TURKEY   -2.224410353 -0.353748547
U.K.      1.538322575  0.270546584
U.S.A.    1.592090810  0.174681097


Posterior Estimates of Variance Structure for Random Coefficients:

Iterations = 1001:10901
Thinning interval = 100
Number of chains = 1
Sample size per chain = 100

1. Empirical mean and standard deviation for each variable,
   plus standard error of the mean:

                           Mean      SD Naive SE Time-series SE
(Intercept):(Intercept) 2.86597 1.06173 0.106173       0.153698
lincomep:(Intercept)    0.45888 0.16832 0.016832       0.023458
lincomep:lincomep       0.07573 0.02713 0.002713       0.002824

2. Quantiles for each variable:

                           2.5%     25%     50%     75%  97.5%
(Intercept):(Intercept) 1.51928 2.15398 2.58854 3.32357 5.3332
lincomep:(Intercept)    0.23322 0.34598 0.42955 0.53947 0.8216
lincomep:lincomep       0.03904 0.05631 0.07156 0.09101 0.1302



Posterior Estimates of Variance:

Iterations = 1001:10901
Thinning interval = 100
Number of chains = 1
Sample size per chain = 100

1. Empirical mean and standard deviation for each variable,
   plus standard error of the mean:

          Mean             SD       Naive SE Time-series SE
     0.0087293      0.0007850      0.0000785      0.0000785

2. Quantiles for each variable:

    2.5%      25%      50%      75%    97.5%
0.007159 0.008182 0.008735 0.009255 0.010253

epoh老师，谢谢您
这个package我当时看到了，但是好像因为一些
原因被cran移除了。您怎么还可以使用？哈哈
而且我看没有提供下载到本地安装的zip文件，只有
提供原代码文件。

的确只有原代码,
须要自行下载"Gnu Scientific Library" build package
  

好的，非常感谢epoh老师
我先按照您给的codes学习下

epoh老师，您自行build package可以正常安装，但是
在R2.15.2运行library(BayesPanel)时出现一个错误提示框：
没有找到libgsl-0.dll，因此这个应用程序未能启动。重新安装应用
程序可能会修复此问题。

您是用哪个版本的R软件？