Global COE Hi-Stat Discussion Paper Series 072

Research Unit for Statisticaland Empirical Analysis in Social Sciences (Hi-Stat)

Bayesian Analysis of Time-Varying Parameter Vector
Autoregressive Model for the Japanese Economy andMonetary Policy

Hi-Stat Discussion Paper

Jouchi Nakajima
Munehisa Kasuya
Toshiaki Watanabe

May 2009

Hi-Stat
Institute of Economic Research
Hitotsubashi University
2-1 Naka, Kunitatchi Tokyo, 186-8601 Japan
http://gcoe.ier.hit-u.ac.jp

Bayesian Analysis of Time-Varying Parameter Vector
Autoregressive Model for the Japanese Economy andMonetary Policy∗
Jouchi Nakajima
Institute for Monetary and Economic Studies, Bank of Japan

Munehisa Kasuya
Research and Statistics Department, Bank of Japan

Toshiaki Watanabe
Institute of Economic Research, Hitotsubashi University

May 2009
Abstract
This paper analyzes the time-varying parameter vector autoregressive (TVP-VAR)
model for the Japanese economy and monetary policy. The time-varying parameters are
estimated via the Markov chain Monte Carlo method and the posterior estimates of parameters reveal the time-varying structure of the Japanese economy and monetary policy
during the period from 1981 to 2008. The marginal likelihoods of the TVP-VAR modeland other VAR models are also estimated. The estimated marginal likelihoods indicatethat the TVP-VAR model best ﬁts the Japanese economic data.
Key words: Bayesian inference, Markov chain Monte Carlo, Monetary policy, State spacemodel, Structural vector autoregressive model, Stochastic volatility, Time-varying parameter.

∗

The authors would like to thank Tomoyoshi Yabu and seminar participants at the Third Japanese-European
Bayesian Econometrics and Statistics Meeting (JEuBES 2008) and the Macro-money workshop at Hitotsubashi
University, and the staﬀ of the Institute for Monetary and Economic Studies (IMES), the Bank of Japan,
for their useful comments. Financial support from the Ministry of Education, Culture, Sports, Science andTechnology of the Japanese Government through the Global COE program “Research Unit for Statistical andEmpirical Analysis in Social Sciences” at Hitotsubashi University is gratefully acknowledged. Views expressedin this paper are those of the authors and do not necessarily reﬂect the oﬃcial views of the Bank of Japan.

11Introduction

The Japanese economy has experienced several distinct periods of macroeconomic activity inrecent decades, resulting in many of Japan’s macroeconomic variables exhibiting changing behaviors over time. Since Miyao (2000, 2002) analyzed the Japanese economy using a vector autoregressive (VAR) model, the time variation of the relations among Japanese macroeconomic
variables has been investigated in several studies (e.g., Fujiwara (2006), Inoue and Okimoto(2008) using a Markov-switching VAR model, and Kimura et al. (2003) using a VAR modelwith time-varying coeﬃcients). In these studies, the changes in the coeﬃcients in the VAR
system are well studied, although the variance of the structural shocks is assumed constantover the sample period or subsample period. This paper estimates a time-varying parameter
vector autoregressive (TVP-VAR) model for the Japanese economy and monetary policy, whichallows both the coeﬃcients and the variance of structural shock to vary over time.
From a methodological viewpoint, the TVP-VAR model has recently become increasingly
popular in the macroeconomics literature following the use of this estimation technique by
Primiceri (2005) using data for the US economy. The spirit of its speciﬁcation is derived
from Cogly and Sargent (2005). Benati and Mumtaz (2005) provide empirical results for the
TVP-VAR model for the UK economy and Baumeister et al. (2008) for the Euro economy.
D’Agostino et al. (2008) show the superior forecasting performance of the TVP-VAR speciﬁcation over other VAR models using US macroeconomic data. We apply their method to theJapanese economy with a slight modiﬁcation. Their scheme of sampling from the posterior
distribution of the stochastic volatility of the TVP-VAR model uses a mixture sampler, which
was originally developed by Kim et al. (1998) in the context of the stochastic volatility modelin ﬁnancial econometrics. The mixture sampler draws sample from the approximated posterior density and its approximation error is small enough to implement the overall model, as
discussed by Omori et al. (2007). However, the multimove sampler proposed by Shephard andPitt (1997) and modiﬁed by Watanabe and Omori (2004) can draw sample from the exact posterior density of the stochastic volatility, and this method is incorporated into the TVP-VAR
model in this paper. Furthermore, Yano and Yoshino (2008) estimate the TVP-VAR modelusing a Monte Carlo particle ﬁltering approach.
In our empirical analysis using Japanese data, a four-variable VAR system is estimated. The
model includes the inﬂation rate, industrial production, nominal short-term interest rate, andmoney supply. The stochastic volatilities and time-varying impulse responses of the macroeco-

2

nomic variables are shown over time. The marginal likelihoods of the TVP-VAR speciﬁcationand other VAR models are also estimated under diﬀerent estimation conditions. The estimatedmarginal likelihood indicates the good performance of the TVP-VAR model.
The paper is organized as follows. In Section 2, we introduce the TVP-VAR speciﬁcation.
Section 3 illustrates the estimation procedure for the TVP-VAR model. Section 4 presents ourempirical results. Finally, Section 5 concludes.

2

Structural VAR models

2.1

Preliminary

We begin with a basic structural VAR model deﬁned as
Ayt = F1 yt−1 + · · · + Fs yt−s + ut ,t = s + 1, . . . , n,(1)where yt is an k × 1 vector of observed variables, A, F1 , . . . , Fs are k × k matrices of coeﬃcients,
and ut is a k ×1 structural shock. We specify the simultaneous relations of the structural shockby recursive identiﬁcation, assuming that A is lower-triangular,

1


0

.. 
. 
.

0 


···
..
.
..
.

0
..



.
 a21
A = 
 ..
..
 .
.

ak1 · · ·ak,k−1 1

We rewrite model (1) as the following reduced form VAR model,
yt = B1 yt−1 + · · · + Bs yt−s + A−1 Σεt ,εt ∼ N (0, Ik ),where Bi = A−1 Fi , for i = 1, . . . , s, and




Σ = 




σ1
0
..
.
0

.

···
..
.
..
.

···

0

0
..

.

..

3


0

.. 
. 
.

0 

σk

The σi (i = 1, . . . , k) is the standard deviation of the structural shock. Stacking the elements
 , . . . , y
in the rows of the Bi ’s to form β (k 2 s × 1 vector), and deﬁning Xt = Ik ⊗ (yt−1t−k ), themodel can be written as
yt = Xt β + A−1 Σ εt .

(2)

All parameters in equation (2) are time-invariant. In the next section, we construct the modelby allowing these parameters to vary over time.

2.2

Time-Varying Parameter VAR

We consider a time-varying parameter VAR (TVP-VAR) model speciﬁed by
yt = Xt βt + A−1t Σt εt ,t = s + 1, . . . , n,(3)where the coeﬃcients βt , and the parameters At , and Σt are all time varying. There would
be many ways to model the process for these time-varying parameters. Let at be a stacked
2 , forvector of the lower-triangular elements in At and ht = (h1t , . . . , hkt ) with hjt = log σjt

j = 1, . . . , k, t = s+1, . . . , n. As suggested by Primiceri (2005), we assume that the parametersin (3) follow a random walk process as follows:
βt+1 = βt + uβt ,at+1 = at + uat ,ht+1 = ht + uht ,






I

εt

O

O

O


 




 

  O Σβ O O 
 uβt  ∼ N 0, 
 ,


 


  O O Σa O  uat 


 

uht
O O O Σhfor t = s + 1, . . . , n, where βs+1 ∼ N (µβ0 , Σβ0 ), as+1 ∼ N (µa0 , Σa0 ) and hs+1 ∼ N (µh0 , Σh0 ).
The shocks to the innovations of the time-varying parameters are assumed uncorrelated among
the parameters βt , at and ht . We further assume that Σβ , Σa and Σh is all diagonal matrices.
The drifting coeﬃcients and parameters are modeled to fully capture possible changes of the
VAR structure over time. Our dynamic speciﬁcation is adequate to permit the parameters tovary even if the shocks in the processes driving the time-varying parameters are uncorrelated.
Note that the log of variance (σt2 ) for the structural shocks is modeled to follow a random
walk process, which belongs to the class of stochastic volatility (e.g., Shephard (2005)). In
ﬁnancial literature, the log volatility (ht ) is often formulated to follow a stationary process such4

as the ﬁrst-order autoregressive process. The random walk process is non-stationary and it
would be undesirable to analyze long-series behavior of market products such as ﬁnancial daily
data, while the purpose of our study is the empirical analysis for the quarterly macroeconomics
data whose sample size is around one hundred. Moreover, as discussed by Primiceri (2005),
the random-walk assumption can capture possible gradual (or sudden) structural changes.

3

Bayesian inference

The TVP-VAR model includes a number of parameters, while the estimation procedure can beconstructed using the MCMC methods. In a Bayesian inference, under certain prior probability
distributions, the MCMC algorithm produces the sample drawn from a high-dimensional posterior distribution of parameters including unobserved latent variables (see e.g., Chib (2001)).
Using the time-varying parameters in our model as latent variables, the model forms a state
space speciﬁcation. The key to constructing an eﬃcient sampling scheme for the TVP-VAR
model is the joint sampling of β = {βt }nt=s+1 (and in turn, a = {at }nt=s+1 , h = {ht }nt=s+1 )conditioned on the rest of the parameters, which is better than the approaches that rely on
one-at-a-time sampling. To accomplish this strategy, the simulation smoother (de Jong andShephard (1995), Durbin and Koopman (2002)) is suitable for sampling the time-varying coeﬃcient β and parameter a because the model can be written in a linear Gaussian state spaceform.
Regarding the stochastic volatility h, the model forms a non-linear non-Gaussian state
space form; thus, we need more technical methods of sampling. One way to sample stochastic
volatility is the approach of Kim et al. (1998), called the mixture sampler. This method has
been widely used in the ﬁnancial and macroeconomics literature (Omori et al. (2007), Primiceri(2005)). The other way is the multimove sampler of Shephard and Pitt (1997) and Watanabe
and Omori (2004). The former method deals with the approximated linear state space model.
As studied by Kim et al. (1998) and Omori et al. (2007), its approximation error is smallenough to implement the original model and can be corrected by the compensation step (seethe details in Kim et al. (1998)). On the other hand, the latter algorithm approaches themodel by drawing from the exact conditional posterior density. Both methods are appropriate
to implement the volatility part in the TVP-VAR model, while we use the latter method inthis paper for its direct sampling from the original form of the model.

5

3.1

MCMC algorithm

Let y = {yt }nt=1 and ω = (Σβ , Σa , Σh ). We set the prior probability density as π(ω) forω. Given the data y, we draw sample from the posterior distribution, π(β, a, h, ω|y), by the
MCMC method. There are several ways of sampling the posterior distribution; we apply thefollowing MCMC algorithm:
1. Initialize β, a, h, ω.
2. Sample β|a, h, Σβ , y.
3. Sample Σβ |β.
4. Sample a|β, h, Σa , y.
5. Sample Σa |a.
6. Sample h|β, a, Σh , y.
7. Sample Σh |h.
8. Go to 2.
As mentioned above, Steps 2 and 4 are conducted with the help of the simulation smoother,and Step 6 requires the multi-move sampler for the stochastic volatility. The assumption of the
diagonal matrix of Σh makes the conditional posterior distribution of {hjt }nt=s+1 independent
with respect to the series j (= 1, . . . , k) and the sampling algorithm for h becomes simple.
Steps 3, 5, and 7 are straightforward, drawing sample from a Wishart or Gamma distribution
under conjugate priors. The details of the procedure are illustrated in Appendix.

3.2

PriorsThe priors should be carefully chosen because the TVP-VAR model has many state variables inthe VAR speciﬁcation and their process is modeled as a non-stationary random walk process.
The TVP-VAR model is so ﬂexible that the state variables can capture both gradual and suddenchanges of the underlying economic structure. On the other hand, allowing time variation in
every parameter in the VAR may cause an over-identiﬁcation problem in the ﬁnite period
data sample. As mentioned by Primiceri (2005), the tight prior for the covariance matrix ofthe disturbance in the random walk process avoids implausible behaviors of the time-varying

6

parameters. A tighter prior should sometimes be avoided in empirical econometrics, although
the TVP-VAR model needs slightly tighter priors to provide reasonable identiﬁcation.
In this paper, the time-varying coeﬃcient (β) needs a tighter prior than the simultaneous
relations (a) and the volatility (h) of the structural shock for the variance of the disturbancein their time-varying process. The structural shock we consider in the model unexpectedly
hits the economic system and its size would more widely ﬂuctuate over time than the possible
change of the autoregressive system of the economic variables speciﬁed by VAR coeﬃcients.
Through the estimation procedure in the following sections, we set a slightly tighter priorfor Σβ and a rather diﬀuse prior for Σa and Σh . Of course, a prior sensitivity analysis isnecessary in this situation to check the robustness of the empirical result with respect to theprior tightness. Section 4.3 reports the robustness check.
An additional remark is required for the prior of the initial state of the time-varying parameters. When the time series to estimate is modeled as a stationary process, we often assume
the initial state follows a stationary distribution of the process. However, our time-varying
parameters are random walks; thus, we specify the prior of the normal distribution for the
initial state of each time-varying parameter. Following Primiceri (2005), we determine themean of these normal priors with the estimates of a time-invariance VAR model using the
pre-sample period. It is possible to specify a ﬂat prior for the initial state to draw the sampleof the posterior fully using information from data. However, we consider it is reasonable touse the economic structure estimated from the pre-sample period up to the start of the main
sample data.

3.3

Marginal likelihood

In a Bayesian framework, we can compare model ﬁt given data using the posterior probabilitiesof the models. The posterior probability of each model is proportional to the prior probabilityof the model, times the marginal likelihood. The ratio of two posterior probabilities is also wellknown as a Bayes factor. If the prior probabilities are assumed equal, we choose the modelthat yields the largest marginal likelihood. The marginal likelihood is deﬁned as the integralof the likelihood with respect to the prior density of the parameter. There are several ways toestimate the marginal likelihood for the model including state variables; we use the harmonicmean method (e.g., Geweke (1999)). The simulation-based harmonic mean estimator of the7marginal likelihood, denoted by m̂(y), is obtained by1 g(ω (j) )
1
,
=
m̂(y)
M
f (y|ω (j) , ϑ(j) )π(ω (j) )M

j=1

where ϑ = (β, a, h), f (y|ω (j) , ϑ(j) ) and π(ω (j) ) denote the likelihood function and prior density,
respectively, and M is the iteration size of the MCMC. If the fraction g(ω)/f (y|ω, ϑ)π(ω) isbounded above, then the approximation is simulation consistent and the rate of convergenceis likely to be practical. While g(ω) can be any p.d.f. with support contained in the parameter
space of the model, Geweke (1999) recommends the choice of g for the modiﬁed harmonic mean
estimator to guarantee the boundness of this fraction. Consider the normal density with thetails truncated, namely
g(ω

(j)

) =


1
exp − (ω (j) − ω̂) V̂ −1 (ω (j) − ω̂)
p/2
1/2
2
τ (2π) |V̂ |1

× I (ω (j) − ω̂) V̂ −1 (ω (j) − ω̂) ≤ χ2τ (p) ,where I[Ω] is an indicator function that takes the value of one if Ω is satisﬁed and zerootherwise, p is the number of elements in ω, and χ2τ (p) denote the τ percentile of the Chisquare distribution with p degrees of freedom. By cutting oﬀ the tails, sample that drop in
the potentially problematic regions are avoided for the computation of the marginal likelihood.
We set ω̂ and V̂ equal to the sample mean and covariance matrix computed from the posteriordraws {ω (j) }M
j=1 , and τ = 0.99 in this paper. As mentioned by Schorfheide (2000), interestingly,
the estimated marginal likelihood does not so serially depend on τ . We also computed themarginal likelihood of the TVP-VAR model with τ = 0.95, 0.90, and found the estimatedmarginal likelihood is seldom sensitive to the value of τ .

4

Empirical evidences for Japanese economy and monetary
policy4.1

Data and estimation procedure

In this section, we estimate the TVP-VAR model for the Japanese economy. Our main dataset
is quarterly and the sample period is from 1981/1Q to 2008/3Q. The pre-sample period from1970/2Q to 1980/4Q is used for the prior distribution of the initial state in the process of8the time-varying parameters, which is explained in detail below. The model includes four
macroeconomic variables: inﬂation rate, industrial production, nominal short-term interest
rate, and monetary base.1 Recursive identiﬁcation is assumed in this order. These variables
are used for a standard VAR model of the Japanese economy, analyzed by several papers(Miyao (2000, 2002), Fujiwara (2006), Inoue and Okimoto (2008), Yano and Yoshino (2008)).
Most of these studies use monthly data for the estimation, while we change the data to aquarterly base by monthly average. The VAR estimation with monthly data often needs manylags (for example, Miyao (2000) sets 10 lags), because the changes and shocks in the economic
variables are considered to aﬀect the other variables of the system with a delay. As mentioned
above, the TVP-VAR model has so many parameters and our choice of quarterly data reduces
the number of parameters to estimate. The lags are determined by the estimated marginal
likelihood. We estimate the model with one to four lags and choose the lags whose marginal
likelihood is the highest among them. As shown below, the main empirical results are obtainedwith two lags.
The following priors are assumed for the i-th diagonals of the covariance matrices:∼ Gamma(10, 0.01), (Σa )−2∼ Gamma(2, 0.01), (Σh )−2∼ Gamma(2, 0.01).
(Σβ )−2iii(4)For the initial state of the time-varying parameters, µβ0 = β̂0 , µa0 = â0 , µh0 = log σ̂02 , andΣβ0 = Σa0 = Σh0 = 4 × I, where β̂0 , â0 and σ̂0 are the OLS estimators obtained using the
pre-sample period.
We draw M = 10, 000 sample after the initial 1,000 sample are discarded. The computational results are generated using Ox version 5.0 (Doornik (2006)). Figure 1 shows the sample
autocorrelation function, the sample paths and the posterior densities for selected parameters.
After discarding the sample in the burn-in period, the sample paths look stable and the sample
autocorrelations drop stably, indicating our sampling method eﬃciently produces uncorrelated
samples.
Table 1 gives the estimates for posterior means, standard deviations, the 95% credibleintervals, the convergence diagnostics (CD) of Geweke (1992) and ineﬃciency factors. Geweke
1

The inﬂation rate is taken from the CPI (consumer price index, general excluding fresh food, the eﬀects
of the increase in the consumption tax removed, and seasonally adjusted). Industrial production is seasonallyadjusted. The nominal short-term interest rate is the overnight call rate. The monetary base is the average
outstanding, adjusted for the reserve requirement ratio changes, and seasonally adjusted. For the sudden andtemporal increases of the monetary base around December 1999 and February 2002, a linear interpolation isused. Except for the call rate, all the variables are transformed in logarithm, and multiplied by 100. In theestimation, we take the ﬁrst diﬀerence of all variables including the call rate.

9

(1992) suggests a comparison between the ﬁrst n0 draws and the last n1 draws, droppingout the middle draws. The CD statistics is computed by CD = (x̄0 − x̄1 )/ σ̂02 /n0 + σ̂12 /n1 ,
mj +nj −1 (i)(i) is the i-th draw, andx,xσ̂j2 /nj is the standard error of x̄j ,where x̄j = n1j i=mjrespectively, for j = 0, 1. If the sequence of the MCMC sampling is stationary, it converges in
distribution to a standard normal. We set m0 = 1, n0 = 1, 000, m1 = 5, 001, and n1 = 5, 000.
The σ̂j2 is computed using a Parzen window with bandwidth, Bm = 500. The ineﬃciency
factor is deﬁned as 1 + 2

Bm
s=1 ρs , where ρs is the sample autocorrelation at lag s, which iscomputed to measure how well the MCMC chain mixes (see e.g., Chib (2001)). It is the ratio of
the numerical variance of the posterior sample mean to the variance of the sample mean from
uncorrelated draws. The inverse of the ineﬃciency factor is also known as relative numericaleﬃciency (Geweke (1992)). When the ineﬃciency factor is equal to m, we need to draw m times
as many MCMC sample as uncorrelated samples. In our empirical result, the null hypothesisof the convergence to the posterior distribution is not rejected for the parameters at the 5%
signiﬁcance level based on the CD statistics, and the ineﬃciency factors are very low, whichindicates an eﬃcient sampling for the parameters and state variables.
1 (Σ β )51 (Σ β )15

1 (Σ β )25

1

0

0

0

0

0

250 500(Σ β )50

250 500(Σ β )15

0.003

1000

1 (Σ h )4

0

0

0

10000
(Σ β )1500.025

0 250 500
3 (Σ h )4
2

0.25

1

0.001
10000
(Σ β )25

1000

1000

500

500

0.000 0.003

0.000

0

10000
(Σ a )3

750

0

100

10000

0

10000

0

10000
(Σ h )4(Σ h )2(Σ a )6

1.0

4

500

500

0.000

0.002
0.001

0.001
0

1 (Σ h )2

0

0.002

0
10000
1500 (Σ β )51 (Σ a )6

250 500
0 250 500
0 250 500
0 250 500
(Σ a )6(Σ h )2(Σ β )25 0.004 (Σ a )3
0.050
0.50
0.003
0.002

0.002

0.001

(Σ a )3

0.003

50

0.003

0.000

0.005

0.00

0.5

2

250
0.05

0 0.5 10

2

4

Figure 1: Estimation results for selected parameters in the TVP-VAR model. Sample autocorrelations (top), sample paths (middle) and posterior densities (bottom).

10Parameter(Σβ )5(Σβ )15(Σβ )25(Σa )3(Σa )6(Σh )2(Σh )4Mean
0.10370.11040.11490.16520.93000.17480.8198

Stdev.
0.03550.03870.04130.06100.87590.10470.4306

95% interval[0.0570, 0.1919][0.0579, 0.2061][0.0594, 0.2175][0.0827, 0.3178][0.0174, 3.4204][0.0411, 0.4419][0.2516, 1.9022]

CD
0.974
0.075
0.700
0.922
0.149
0.215
0.622

Ineﬃciency
8.01
10.11
10.90
24.62
68.76
23.41
44.64

Table 1: Estimation results for selected parameters in the TVP-VAR model (the estimates ofΣβ and Σa are multiplied by 100).

4.2

Empirical results

This section provides quantitative empirical results for the Japanese economy and monetary
policy of the TVP-VAR model. Through the sample period from 1981 to 2008, the Japanese
economy experienced several diﬀerent periods. The bubble economy burst in the early 1990s,
and the Bank of Japan introduced the zero interest rate policy from 1999 to 2000 and thequantitative easing policy from 2001 to 2006. Using the TVP-VAR model, we investigate thetime-varying structure of the Japanese economy and monetary policy, as follows.
4.2.1

Time-varing volatility

Figure 2 plots the series for estimated stochastic volatility of the structural shock on four
variables, based on the posterior mean, and the 16-th and 84-th quantile2 intervals of thestandard deviation of the shock, σit = exp(hit /2). It presents the dynamics of the volatility
over time, which diﬀers across variables.
Figure 2(i) plots the time-varying volatility of the inﬂation rate. The period 1981–1985
exhibits a relatively higher volatility because of the second oil shock. The additional reduction
is observed around 1998 to 2006, when the Japanese economy experienced deﬂation. The timevarying volatility of industrial production, displayed by Figure 2(ii), shows that its standard
deviation increased substantially from the mid-1990s to the ﬁrst half of the 2000s comparedwith the preceding period. Sakura et al. (2005) point out that the variance of real GDP became
more volatile during the 1990s compared with the 1980s. They estimate the time-invariant
VAR model by dividing the sample into subsamples and ﬁnd that both the parameters in the
VAR and the shock contribute to the ﬂuctuation of real GDP. One source of the ﬂuctuationis considered to be the higher volatility of ﬁrm investment after the bubble economy burst, as
2
As Primiceri (2005) uses, the 16-th and 84-th quantiles correspond to the one-standard-deviation band
under normality.

11

(i) Inflation rate(ii) Industrial production

0.4

2.5

0.3

2.0

0.2

1.5
1.0

0.1
1985

1990

1995

2000

2005

(iii) Call rate

1985
7.5

1.0

1990

1995

2000

2005

2000

2005

(iv) Monetary base

5.0

0.5

2.5

1985

1990

1995

2000

2005

1985

1990

1995

Figure 2: Posterior estimates for stochastic volatility of structural shock, σit = exp(hit /2).
Posterior mean (solid line), 16-th and 84-th percentiles (dotted line).
discussed by Sakura et al. (2005). The same ﬂuctuation in the shock of industrial production
is observed in our analysis. The TVP-VAR model has an advantage over the constant VAR
models in the sense that it need not divide data into subsamples to conﬁrm the change of thestructure of the model.
Figure 2(iii) and (iv) plots the series of stochastic volatility for the call rate and the monetary base, respectively. The BOJ’s policy instrument of the overnight uncollateralized callrate decreased to 0.25% in 1995. Then, the BOJ implemented the zero interest rate policy3 in1999 and the quantitative easing policy in 2001. During the quantitative easing policy period,
the outstanding balance of the current account is the target of the BOJ’s monetary policy,
and expansions of the monetary base were observed by March 2006, when the BOJ altered thepolicy instrument to the call rate again. The estimated variances of the call rate and monetarybase show these changes in monetary policy; that is, the stochastic volatility of the call rate
3
In our speciﬁcation of the TVP-VAR model, the zero lower bound of the nominal interest rate is notconsidered. While it may be possible to incorporate the lower bound for a certain variable on the TVP-VAR
model, it requires a computational burden. Thus, in the current paper, this issue is left for future work.
Regarding the model ﬁt of the time-varying parameters, the estimation of the subsample period that excludes
the zero interest rate regime is also examined below. One way to solve this econometric issue is well discussed
by Iwata and Wu (2006) with a censored variable.

12drops close to zero around the mid-1990s and the monetary base increases rapidly from the
beginning of the 2000s, and a hike is marked in 2006, which corresponds to the termination of
the quantitative easing policy. The spike around 1986 for the volatility of the call rate would
indicate that the easing policy in this period is measured as the structural shock. After that,the volatility of the call rate is relatively low and looks moderate.
Overall, the time-varying volatility contributes to the VAR estimation, identifying thestructural shock with the appropriate variance of the shock size. For the data we analyze
here, the estimates of the time-invariant VAR model with constant parameters would result inbiases in the covariance matrix for the disturbances and at the same time in the autoregressive
coeﬃcients because of the misspeciﬁcation of the dynamics of the parameters. Again, thenominal interest rate lower bound is not explicitly incorporated into the model, although theestimated volatility for the call rate that stays low enough would help the VAR system toidentify exogenous shocks better than the time-invariant VAR model.
4.2.2

Time-varying impulse responses

Impulse response analyses are provided for the time series in the TVP-VAR model. Because thecoeﬃcients are time varying, the impulse responses are calculated at each date over the sample
period. In our study, the shock size for the responses is set equal to the time-series average ofthe stochastic volatility for each series over the sample period. Note that the shock size for the
response is not based on the estimated variance each time. The impulse summarizes the eﬀects
of average-sized experimental structural shocks hitting the VAR system. The impulse responses
are calculated for each iteration of the MCMC with the current draw of the parameters, andthe sample mean and standard deviation of the responses are computed. Figure 3–6 plot theresults of the main responses.
Figure 3 reports the impulse responses of industrial production to a positive interest
rate shock. The time series of the response and the responses in three diﬀerent dates with
one-standard-deviation bands are shown.

The dates for the comparison are (ii)1986/1Q,(iii)1994/1Q, and (iv)2002/1Q, which are chosen arbitrarily, while the overall time variationof the response is summarized in the time-series plot(i). The period around (iii)1994/1Q represents the typical economic conditions after the bubble burst and (iv)2002/1Q under thequantitative easing policy. The interest rate shock has a negative eﬀect on industrial production, although its impact seems to exhibit time variation. It is natural that the response ofindustrial production becomes smaller from the second half of the 1990s to the 2000s because

13(ii) 1986/1Q(i) Time−varying impulse response
0.00
0.00
−0.05

−0.05
1 year
2 years
3 years

−0.10

−0.10

−0.15
1985
1990
(iii) 1994/1Q

1995

2000

2005

0

0.00

0.00

−0.05

−0.05

−0.10

−0.10

−0.15
0

5

10

15

5

10

15

(iv) 2002/1Q

−0.15
5

10

15

0

Figure 3: Impulse responses of industrial production to a positive interest rate shock. (i)Time series of response after one, two and three years, (ii)–(iv) response in each period with
one-standard-deviation bands.
of the environment where the call rate stays close to zero. This behavior of the response isconsistent with the result of Kimura et al. (2003). In addition, the estimation result indicatesthat the bottom of the impact of the interest rate shock is around 2006, when the quantitativeeasing policy is terminated. The time-varying impulse response implies that the increase in
the policy interest rate in 2006 does not decrease industrial production.
Second, Figure 4 shows the responses of inﬂation to a positive interest rate shock. A risein inﬂation after a monetary tightening using VAR estimates is well known as the price puzzle(Sims (1992), Eichenbaum (1992)). In our empirical estimates of the TVP-VAR model, theimpulse response exhibits a diﬀerent shape in each period. For 1981–1987, the impulse response
falls suﬃciently in response to the positive interest rate shock and Figure 4(ii) shows a slightand temporary price puzzle in 1986/1Q. On the other hand, from 1988 to the mid-1990s, theimpulse response keeps rising and stays high even at a three-year horizon. As plotted in 4(iii),
the impulse response is a positive domain over all of 1990/1Q. There are several explanations
for the price puzzle and the VAR estimates would depend on the lags or some omitted variablesin the VAR system to some extent. Figure 4(iv) shows that the impulse response stays near

14

0.04

(i) Time−varying impulse response(ii) 1986/1Q
0.05

0.02
0.00

0.00
−0.02
1 year
2 years
3 years

−0.04
1985

1990

1995

2000

−0.05

2005

0

0.05

0.05

0.00

0.00

−0.05

−0.05

0

5

10

15

5

10

15

(iv) 2002/1Q(iii) 1990/1Q

5

10

15

0

Figure 4: Impulse responses of inﬂation to positive interest rate shock. (i) Time series of
response after one, two and three years, (ii)–(iv) response of each period with one-standarddeviation bands.
zero and the one-standard-deviation bands includes zero for 2002/1Q because the call rate
is close to zero and the structural shock to the interest rate diminishes in this period. Ourestimation result indicates that the TVP-VAR model has an advantage in the sense that it canassess the time-varying dynamics of the relation among economic variables.
Third, the responses of the interest rate to a positive inﬂation rate shock are shown in
Figure 5. This response is one measure of monetary policy activism with respect to inﬂation.
Similar to the previous ﬁgure, the response of the interest rate has time variation and theimpact of the response seems to diminish for the period of the zero interest rate policy and thequantitative easing policy. In Figure 5(i), the time-series line of the response after one year
does not seem to ﬂuctuate for 1981–1998, while the response after two and three years showsa high degree of time variation commencing in the early 1990s. It means that the response
increases smoothly in reaching the three-year horizon of the response. For the other periods,
the response seems to reach the three-year horizon of the response in the ﬁrst year. From
another point of view, the time-varying impulse response shows variation of the weight ofinﬂation in the monetary policy reaction function throughout the sample period.

15(i) Time−varying impulse response(ii) 1986/1Q
0.10

0.15
1 year
2 years
3 years

0.10

0.05

0.05

0.00

0.00

−0.05
1985
1990
(iii) 1994/1Q

1995

2000

2005

0

0.10

0.10

0.05

0.05

0.00

0.00

−0.05

−0.05

0

5

5

10

15

5

10

15

(iv) 2002/1Q

10

15

0

Figure 5: Impulse responses of interest rate to positive inﬂation rate shock. (i) Time series of
response after one, two and three years, (ii)–(iv) response in each period with one-standarddeviation bands.
Fourth, Figure 6 reports the impulse responses of industrial production to a positive monetary base shock. Overall, the monetary base shock has a positive eﬀect on industrial production
during the sample period including the quantitative easing policy period, and in particular, thetime series of the response after one year seems to be stable for 1990 to 2000, while the size ofthe response decreases around the mid-2000s. Throughout the sample period, industrial production seems to reach the three-year horizon response level in two years. Because we observethat the rise in the monetary base is followed by a positive interest rate shock in our estimation
results (not shown here), the monetary base shock can be interpreted as a shock to the money
demand function. This indicates that the money demand shock would have a positive eﬀect
through bank lending or other ﬁnancial market variables. However, when the nominal interest
rate is close to zero, the response becomes smaller than the one in the 1980s and the ﬁrst half ofthe 1990s. For the impulse response of (iv)2002/1Q, the one-standard-deviation bands include
zero. This result indicates that the eﬀect of an increase in the monetary base is uncertain, as
discussed by Kimura et al. (2003) and Fujiwara (2006). Similar to the interest rate shock, theimpulse response of the monetary base shock reaches its bottom around 2006, when the quan-16

(i) Time−varying impulse response
1.5(ii) 1986/1Q
1 year
2 years
3 years

2

1.0

1

0.5

0

1985

1990

1995

2000

2005

0

2

2

1

1

0

0

0

5

10

15

5

10

15

(iv) 2002/1Q(iii) 1994/1Q

5

10

15

0

Figure 6: Impulse responses of industrial production to positive monetary base shock. (i)Time series of response after one, two and three years, (ii)–(iv) response of each period with
one-standard-deviation bands.
titative easing policy is terminated and the monetary base decreases. The estimation results
show that the rapid decrease in the monetary base does not decrease industrial production
after 2006. It implies that this period would have been an appropriate time to terminate thequantitative easing policy.

4.3

Model selection and robustness

In the previous section, we showed the empirical results of the TVP-VAR model for the main
dataset with the speciﬁed priors. Below, we examine whether the TVP-VAR model also ﬁts thedata better than other VAR models. We estimate the marginal likelihood using the modiﬁed
harmonic mean method in Section 3.3 for diﬀerent model speciﬁcations, lags, and priors. The
competing model speciﬁcations are a constant parameter VAR model (CP-VAR) and a semitime-varying parameter VAR model (STVP-VAR). The CP-VAR model refers to the equation(2). The STVP-VAR model allows only the coeﬃcients (β) and the simultaneous relations ofstructural shock (a) to vary over time, and the volatility of structural shock is time invariant
with Σt = Σ, for all t = s + 1, . . . , n. By comparing the model ﬁt between the TVP-VAR

17

Lag

CP-VARSTVP-VAR

TVP-VARLag(i) Full sample, Prior1
1
2
3
4-519.01-548.76-581.30-621.13-267.08-292.63-315.38-350.11

STVP-VAR

TVP-VAR(ii) Full sample, Prior2-457.41-430.67-404.67-442.36-339.11-311.24-327.12-346.301234(iii) Subsample, Prior1
1
2
3
4

CP-VAR-535.90-577.73-622.73-675.30-515.79-527.24-535.14-567.05-369.87-382.85-397.44-407.87(iv) Subsample, Prior2-266.23-282.30-292.20-315.50-265.31-278.17-279.30-285.771
2
3
4-284.13-321.82-357.18-404.59-298.07-307.57-322.95-347.79-263.43-287.26-305.14-329.82

Table 2: Estimated marginal likelihood for competing VAR models (logarithm scale; CP:constant parameter, STVP: semi time-varying parameter, TVP: time-varying parameter).
and STVP-VAR models, the advantage of the stochastic volatility of the structural shock, asemphasized by Cogly and Sargent (2005), is examined. The lag is set one up to four. We
specify additional priors for the CP-VAR and STVP-VAR models and alternative priors for
the TVP-VAR model, as follows:
Prior1 :(4) and β ∼ N (0, 2 × I),ᾱ ∼ N (0, 2 × I),

σi−1 ∼ Gamma(2, 0.02),

Prior2 :∼ Gamma(20, 0.01), (Σa )−2∼ Gamma(5, 0.01), (Σh )−2∼ Gamma(2, 0.01),(Σβ )−2iiiβ ∼ N (0, 10 × I),ᾱ ∼ N (0, 10 × I),

σi−1 ∼ Gamma(2, 0.02),

where ᾱ denotes the stacked vector of the elements in A.
In addition to the estimation with the full sample period, we estimate the marginal likelihood using the subsample period from 1981/1Q to 1995/4Q. The period of the zero interest
rate policy and the monetary easing policy is omitted in this subsample period. We estimate
three competing models for the subsample data with Prior1 and Prior2.
Table 2 reports the estimated marginal likelihoods for competing models under four conditions: (i) full sample, Prior1, (ii) full sample, Prior2, (iii) subsample, Prior1, and (iv) subsample, Prior2. The iteration size of MCMC is the same as in the main empirical results.
For all the conditions, the TVP-VAR model has the highest marginal likelihood. The STVPVAR model is favored over the CP-VAR model model, while the marginal likelihood of the
TVP-VAR model is much higher than for the STVP-VAR model for the full-sample period.

18Interestingly, the diﬀerence in the marginal likelihood between the competing models is
smaller for the subsample period than for the full-sample period under both priors. In particular, the marginal likelihood of the STVP-VAR model becomes close to the one for the
TVP-VAR model in the subsample estimation. One possible reason for this result is that the
full-sample period includes the zero interest rate policy and the monetary easing policy; thus,
the posterior distribution of constant volatility estimated through the full-sample period would
be biased and unable to follow the volatility dynamics, especially that seen for the call rate
and money supply, as shown in Figure 2. We understand that the TVP-VAR speciﬁcation is avery ﬂexible model for analyzing the economic variables including the period when the nominalinterest rate is close to zero.
Overall, the time-varying parameter contributes to the VAR model based on the estimated marginal likelihoods. The time-varying coeﬃcients increase the marginal likelihood,
although not by enough, and the time-varying simultaneous relations and stochastic volatilityof structural shock contributes the model ﬁt over our data set, especially for the sample period
including the zero interest rate policy. These results are robust as shown by the estimations
under diﬀerent priors and lags.

5

Conclusion

This paper analyzes the TVP-VAR models of the Japanese economy and monetary policy. The
time-varying parameters are estimated via the Markov chain Monte Carlo method and theposterior estimates of parameters reveal the time-varying structure of the Japanese economy
and monetary policy during the period from 1981 to 2008. The marginal likelihoods of the
TVP-VAR model and other VAR models are estimated under diﬀerent priors, lags, and sample
periods. The estimated marginal likelihoods indicate that the TVP-VAR model best ﬁts theJapanese economic data over; in particular, it marks much higher marginal likelihoods for the
sample period that includes the zero interest rate policy.

19

Appendix. MCMC algorithm for the TVP-VAR model
A.1

Sampling β

To sample β from the conditional posterior distribution π(β|a, h, Σβ , y), we write the model inthe state space form as
yt = Xt βt + A−1t Σt εt ,βt+1 = βt + uβt ,t = s + 1, . . . , n,t = s, . . . , n − 1,

where βs = µβ0 , and uβs ∼ N (0, Σβ0 ). It is possible to draw sample from the joint posterior
distribution π(βs+1 , . . . , βn |a, h, Σβ , y). Following de Jong and Shephard (1995), we show thealgorithm of the simulation smoother on the state space modelyt = Zt αt + Gt ut ,t = s + 1, . . . , n,αt+1 = Tt αt + Ht ut ,t = s, . . . , n − 1,

where ut ∼ N (0, I) and Gt Ht = O. The simulation smoother draws η = (ηs , . . . , ηn−1 ) ∼π(η|y, θ) where ηt = Ht ut for t = s, . . . , n − 1, and θ denotes the rest of the parameters in themodel. We run the Kalman ﬁlter:
et = yt − Zt at ,Dt = Zt Pt Zt + Gt Gt ,Kt = Tt Pt Zt Dt−1 ,Lt = Tt − Kt Zt ,at+1 = Tt at + Kt et ,Pt+1 = Tt Pt Lt + Ht Ht ,for t = s + 1, . . . , n, where as+1 = Ts αs and Ps+1 = Hs Hs . Then, letting Λt = Ht Ht , we run
the simulation smoother:
Ct = Λt − Λt Ut Λt ,

ηt = Λt rt + εt ,rt−1 = Zt Dt−1 et + Lt rt − Vt Ct−1 εt ,εt ∼ N (0, Ct ),Vt = Λt Ut Lt ,Ut−1 = Zt Dt−1 Zt + Lt Ut Lt + Vt Ct−1 Vt ,for t = n, n−1, . . . , s+1, with rn = Un = 0. Finally, we can draw ηs = Λs rs +εs , εs ∼ N (0, Cs )with Cs = Λs − Λs Us Λs . We construct the sample of {αt }nt=s+1 via the state equation usingn−1{ηt }t=s
drawn through the simulation smoother. In the case of sampling β, we coordinate
1/2

the parameters as Zt = Xt , Tt = I, Gt = (A−1t Σt , O), Ht = (O, Σβ ), for t = s + 1, . . . , n,1/2

Ts αs = µβ0 , and Hs = (O, Σβ0 ).

20

A.2

Sampling a

For sampling a from its conditional posterior distribution π(a|β, h, Σa , y), the expression of thestate space form with respect to a is beneﬁcial to see the implementation, namely,ŷt = X̂t at + Σt εt ,t = s + 1, . . . , n,t = s, . . . , n − 1,

at+1 = at + uat ,where as = µa0 , uas ∼ N (0, Σa0 ), ŷt = yt − Xt βt , and


0

···


 −ŷ1t
0
0
···


 0
0
···
−ŷ1t −ŷ2tX̂t =  0
00−ŷ1t · · ·

 .
..
 ..
.

0
···
0

0
..
.

0

···

−ŷ1t · · ·

0








,






−ŷk−1,t

for t = s + 1, . . . , n. Similar to sampling β, we use the simulation smoother with Zt = X̂t ,
1/2

1/2

Tt = I, Gt = (Σt , O), Ht = (O, Σa ) for t = s + 1, . . . , n, Ts αs = µa0 , and Hs = (O, Σa0 ).

A.3

Sampling h

The state space equations for state variable h belong to a non-linear state space class. Because
we assume Σh and Σh0 are diagonal matrices, we make the inference for {hjt }nt=s+1 separately∗ denote the i-th element of A ŷ . Then, we can write
for j (= 1, . . . , k). Let yit
t t
∗
yit
= exp(hit /2)εit ,t = s + 1, . . . , n,hi,t+1 = hit + ηit , t = s, . . . , n − 1,


 
εit
1 0

 ∼ N 0, 
 ,
ηit
0 vi2


where ηis ∼ N (0, vi20 ), and vi2 and vi20 are the i-th diagonal elements of Σh and Σh0 , respectively,and ηit is the i-th element of uht . We sample (hi,s+1 , . . . , hin ) from its conditional posterior
density by running the multi-move sampler for the non-linear Gaussian state space model,
developed by Shephard and Pitt (1997) and Watanabe and Omori (2004).

21

We divide (hi,s+1 , . . . , hin ) into K + 1 blocks and sample each block from its joint posterior
density conditioned on the other blocks and parameters. Now we focus on sampling a typical
block (hir , . . . , hi,r+d ) from its joint posterior density, (note that r ≥ s + 1, d ≥ 1, r + d ≤ n).
The conditional posterior distribution is given byπ(ηi,r−1 , . . . , ηi,r+d−1 |θi ) ∝ r+d−1

∗2

yit
×
exp
−
f (ηit ) × f (hi,r+d ),hhit /2
it
2e
e
t=r
t=r−1

r+d


1

(5)

where



2

η

is

(if t = s), exp − 22vi0
f (ηit ) =


2

ηit


 exp − 2
(if t ≥ s + 1),
2vi


2

 exp − (hi,r+d+1 − hi,r+d )(if r + d < n),2vi2
f (hi,r+d ) =
1(if r + d = n),∗ , . . . , y∗and θi = (hi,r−1 , hi,r+d+1 , yiri,r+d , vi , vi0 ). The posterior draw of (hir , . . . , hi,r+d ) can beobtained by running the state equation with this draw of (ηi,r−1 , . . . , ηi,r+d−1 ) given hi,r−1 .
We sample (ηi,r−1 , . . . , ηi,r+d−1 ) from (5) using the AR-MH algorithm (Tierney (1994),
Chib and Greenberg (1995)) with the following proposal distribution. Our construction of the
proposal density begins with the second-order Taylor expansion ofg(hit ) ≡ −

hit
y ∗2
− ith ,
22e itaround the certain point ĥt discussed later, namely,
1g(hit ) ≈ g(ĥit ) + g  (ĥit )(hit − ĥit ) + g  (ĥit )(hit − ĥit )2
2
2



1 
g (ĥit )g (ĥit ) hit − ĥit −∝,2g  (ĥit )
where
y ∗2
1
g  (ĥit ) = − + ith ,2 2e itg  (ĥit ) = −

22

∗2
yit
.
2ehit

We use the proposal density formed as(h∗ − hit )2
exp − it ∗2q(ηi,r−1 , . . . , ηi,r+d−1 |θi ) ∝
2σit
t=r
r+d


×

r+d−1


f (ηit ),

t=r−1

where
∗2
=−
σit

1
g  (ĥit ),∗2 h∗it = ĥit + σitg (ĥit ),(6)

for t = r, . . . , r + d − 1 and t = r + d (when r + d = n). For t = r + d (when r + d < n),∗2=σi,r+dh∗i,r+d1,−g  (ĥi,r+d ) + 1/vi2∗2g  (ĥi,r+d ) − g  (ĥi,r+d )ĥi,r+d + hi,r+d+1 /vi2 .
= σi,r+d(7)(8)

The choice of this proposal density is derived from its correspondence to the state space modelh∗it = hit + ζit ,t = r, . . . , r + d,

hi,t+1 = hit + ηit , t = r − 1, . . . , r + d − 1,

 

∗2
0
ζ
σ
 it  ∼ N 0,  it
 , t = r, . . . , r + d,
ηit
0 vi2

(9)



with ηi,r−1 ∼ N (0, vi2 ) when r ≥ 2, and ηis ∼ N (0, vi20 ). Given θi we draw the candidate pointof (ηi,r−1 , . . . , ηi,r+d−1 ) for the AR-MH algorithm by running the simulation smoother over thestate space representation (9).
Now we come to ﬁnd (ĥir , . . . , ĥi,r+d ) which is near the mode of the posterior density for
an eﬃcient sampling. We loop the following steps several times to reach near the mode:
1. Initialize (ĥir , . . . , ĥi,r+d ).
2∗ , . . . , σ 2∗ ) by (6) and (8).
2. Compute (h∗ir , . . . , h∗i,r+d ), (σiri,r+d2∗ , . . . , σ 2∗ ) on (9) and3. Run the simulation smoother using current (h∗ir , . . . , h∗i,r+d ), (σiri,r+dobtain ĥ∗it ≡ E(hit |θi ) for t = r, . . . , r + d.
4. Replace (ĥir , . . . , ĥi,r+d ) by (ĥ∗ir , . . . , ĥ∗i,r+d ).
5. Go to 2.

23

Note that the E(hit |θi ) is the product in the simulation smoother as Λt rt with εt = 0.
Finally, we remark the selection of the blocks. A typical block can be denoted as (hi,kj−1 +1 ,. . . , hi,kj ) for j = 1, . . . , K + 1 with k0 = s and kK+1 = n and Shephard and Pitt (1997)
suggests the determination of (k1 , . . . , kK ), called stochastic knots. It leads ton(j + Uj ),kj = int
K +2


for j = 1, . . . , K, where Uj is the random sample from the uniform distribution U [0, 1]. We
randomly choose (k1 , . . . , kK ) for every iteration of MCMC sampling to obtain a ﬂexible drawof (hi,s+1 , . . . , hin ).

A.4

Sampling ωGiven β, we can derive the conditional posterior density of Σβ . Let σβi denote the i-th
diagonal element of Σβ . Because we assume that Σβ is a diagonal matrix, we sample σβi
independently for i = 1 . . . , k. If we specify the prior as σβ−2
∼ Gamma(sβ0 /2, Sβ0 /2), wei|β ∼ Gamma(ŝβi /2, Ŝβi /2) whereobtain the conditional posterior distribution as σβ−2iŝβi = sβ0 + n − s − 1,Ŝβi = Sβ0 +n−1

(βi,t+1 − βit )2 ,t=s+1

and βit is the i-th element of βt . The gamma prior is conjugate in this case and the posteriordraw is straightforward. We sample the diagonal elements of Σa |a and Σh |h in the same way.

24

References
Baumeister, C., E. Durinck, and G. Peersman (2008). Liquidity, inﬂation and asset prices in a timevarying framework for the Euro area. Manuscript.
Benati, L. and H. Mumtaz (2005). The ‘Great Stability’ in the U.K.: good policy or good luck?
Manuscript, Bank of England.
Chib, S. (2001). Markov chain Monte Carlo methods: Computation and inference. In J. J. Heckman andE. Leamer (Eds.), Handbook of Econometrics, Volume 5, pp. 3569–3649. Amsterdam: North-Holland.
Chib, S. and E. Greenberg (1995). Understanding the Metropolis-Hastings algorithm. The American
Statistician 49, 327–335.
Cogly, T. and T. J. Sargent (2005). Drifts and volatilities: monetary policies and outcomes in the post
WWII U.S. Review of Economic Dynamics 8, 262–302.
D’Agostino, A., L. Gambetti, and D. Giannone (2008). Macroeconomic forecasting and structuralchange. Manuscript.
de Jong, P. and N. Shephard (1995). The simulation smoother for time series models. Biometrika 82,
339–350.
Doornik, J. (2006). Ox: Object Oriented Matrix Programming. London: Timberlake Consultants Press.
Durbin, J. and S. J. Koopman (2002). Simple and eﬃcient simulation smoother for state space timeseries analysis. Biometrika 89, 603–616.
Eichenbaum, M. (1992). Interpreting the macroeconomic time series facts: the eﬀects of monetary
policy – comments. European Economic Review 36, 1001–1011.
Fujiwara, I. (2006). Evaluating monetary policy when nominal interest rates are almost zero. Journalof Japanese and International Economies 20, 434–453.
Geweke, J. (1992). Evaluating the accuracy of sampling-based approaches to the calculation of posterior
moments. In J. M. Bernardo, J. O. Berger, A. P. Dawid, and A. F. M. Smith (Eds.), Bayesian
Statistics, Volume 4, pp. 169–188. New York: Oxford University Press.
Geweke, J. (1999). Using simulation methods for Bayesian econometric models: inference, development
and communication. Econometric Reviews 18, 1–126.
Inoue, T. and T. Okimoto (2008). Were there structural breaks in the eﬀects of Japanese monetary policy? Re-evaluating policy eﬀects of the lost decade. Journal of Japanese and International
Economies 22, 320–342.

25

Iwata, S. and S. Wu (2006). Estimating monetary policy eﬀects when interest rates are close to zero.
Journal of Monetary Economics 53, 1395–1408.
Kim, S., N. Shephard, and S. Chib (1998). Stochastic volatility: likelihood inference and comparison
with ARCH models. Review of Economic Studies 65, 361–393.
Kimura, T., H. Kobayashi, J. Muranaga, and H. Ugai (2003). The eﬀect of the increase in monetarybase on Japan’s economy at zero interest rates: an empirical analysis. In Bank for International
Settlement (Ed.), Monetary policy in a changing environment, 276–312.
Miyao, R. (2000). The role of monetary policy in Japan: a break in the 1990s? Journal of Japaneseand Interenational Economies 14, 366–384.
Miyao, R. (2002). The eﬀects of monetary policy in Japan. Journal of Money, Credit, and Banking 34,
376–392.
Omori, Y., S. Chib, N. Shephard, and J. Nakajima (2007). Stochastic volatility with leverage: fast
likelihood inference. Journal of Econometrics 140, 425–449.
Primiceri, G. E. (2005). Time varying structural vector autoregressions and monetary policy. Review
of Economic Studies 72, 821–852.
Sakura, K., H. Sasaki, and M. Higo (2005). Japan’s economic ﬂuctuation since the 1990s: fact ﬁnding(in Japanese). Bank of Japan Working Paper Series, 05-J-10.
Schorfheide, F. (2000). Loss function-based evaluation of DSGE models. Journal of Applied Econometrics 15, 645–670.
Shephard, N. (Ed.) (2005). Stochastic Volatility: Selected Readings. Oxford: Oxford University Press.
Shephard, N. and M. Pitt (1997). Likelihood analysis of non-Gaussian measurement time series.
Biometrika 84, 653–667.
Sims, C. A. (1992). Interpreting the macroeconomic time series facts: the eﬀects of monetary policy.
European Economic Review 36, 975–1001.
Tierney, L. (1994). Markov chains for exploring posterior distributions. Annals of Statistics 21, 1701–
1762.
Watanabe, T. and Y. Omori (2004). A multi-move sampler for estimating non-Gaussian time seriesmodels: comments on Shephard and Pitt (1997). Biometrika 91, 246–248.
Yano, K. and N. Yoshino (2008). Japanese monetary policy reaction function and time-varying structuralvector autoregressions: a Monte Carlo particle ﬁlering approach. Manuscript.

26