Stochastic volatility models are a class of models for non-stationary time series data. The log of the variance of the series is assumed to evolve according to a latent Gaussian process. The equation below shows a univariate stochastic volatility model with an AR(1) latent state:
\[\begin{align*} Y_t &= \varepsilon_t\exp(\alpha_t/2), &\varepsilon_t &\sim \mathcal{N}(0, 1), \\ \alpha_t &= \mu + \phi(\alpha_{t-1} - \mu) + \eta_t, &\eta_t &\sim \mathcal{N}(0, \sigma_\eta^2).\numberthis \label{eq:usv-spec} \end{align*}\]\(Y_t\), \(t = 1,\dots,T\) are discrete observations of a time series process. The latent state represents the log-volatility and \(\alpha_t\) follows an AR(1) process with mean \(\mu\).
The figure below shows a simulation from this model with fixed values of the parameters, \(\phi = 0.8, \mu = 1, \sigma_\eta = 0.1\), the initial state of the log-volatility is the stationary solution of the AR(1) process: \(\alpha_0 \sim \mathcal{N}\left(\mu, \frac{\sigma_\eta^2}{1 - \phi^2}\right)\).
To define this model in Scala and simulate observations, first import the required libraries and use the method simulate which initialises a Breeze MarkovChain
object with the initial state equivalent to the stationary solution of the AR(1) process.
scala> import com.github.jonnylaw.dlm._
<console>:12: warning: Unused import
import com.github.jonnylaw.dlm._
^
import com.github.jonnylaw.dlm._
scala> import cats.implicits._
<console>:10: warning: Unused import
import com.github.jonnylaw.dlm._
^
<console>:15: warning: Unused import
import cats.implicits._
^
import cats.implicits._
scala> val p = SvParameters(phi = 0.8, mu = 2.0, sigmaEta = 0.3)
<console>:13: warning: Unused import
import cats.implicits._
^
p: com.github.jonnylaw.dlm.SvParameters = SvParameters(0.8,2.0,0.3)
scala> val sims = StochasticVolatility.simulate(p).
| steps.take(1000).toVector.tail.map { case (t, y, a) => (t, y)}
<console>:13: warning: Unused import
import cats.implicits._
^
sims: scala.collection.immutable.Vector[(Double, Option[Double])] = Vector((2.0,Some(3.846583926103539)), (3.0,Some(-0.23095146212092024)), (4.0,Some(-3.4143942394698206)), (5.0,Some(1.7572309122071414)), (6.0,Some(2.1047649586787083)), (7.0,Some(-2.2661717069966496)), (8.0,Some(0.4066468240993743)), (9.0,Some(-0.3843542564748422)), (10.0,Some(-4.271175785394365)), (11.0,Some(2.2641525966108533)), (12.0,Some(0.832053197800046)), (13.0,Some(-1.071886957638012)), (14.0,Some(-4.587727837914248)), (15.0,Some(0.4561115567367286)), (16.0,Some(0.7214564076353294)), (17.0,Some(-1.2356909193260661)), (18.0,Some(3.9836016044778764)), (19.0,Some(1.2202853102451159)), (20.0,Some(0.562178388757327)), (21.0,Some(5.0632597721313735)), (22.0,Some(-1.0518501203753887)), (23.0,Some(0.7311813047530124)), ...
To determine the posterior distribution of the parameters specify the prior distributions of the parameters:
scala> import breeze.stats.distributions.Gaussian
<console>:10: warning: Unused import
import com.github.jonnylaw.dlm._
^
<console>:13: warning: Unused import
import cats.implicits._
^
<console>:18: warning: Unused import
import breeze.stats.distributions.Gaussian
^
import breeze.stats.distributions.Gaussian
scala> val priorPhi = Gaussian(0.8, 0.1)
<console>:10: warning: Unused import
import com.github.jonnylaw.dlm._
^
<console>:13: warning: Unused import
import cats.implicits._
^
priorPhi: breeze.stats.distributions.Gaussian = Gaussian(0.8, 0.1)
scala> val priorMu = Gaussian(2.0, 1.0)
<console>:10: warning: Unused import
import com.github.jonnylaw.dlm._
^
<console>:13: warning: Unused import
import cats.implicits._
^
priorMu: breeze.stats.distributions.Gaussian = Gaussian(2.0, 1.0)
scala> val priorSigma = InverseGamma(2.0, 2.0)
<console>:13: warning: Unused import
import cats.implicits._
^
<console>:15: warning: Unused import
import breeze.stats.distributions.Gaussian
^
priorSigma: com.github.jonnylaw.dlm.InverseGamma = InverseGamma(2.0,2.0)
There are two inference algorithms defined for the stochastic volatility model, either one can be used and should produce identical inferences:
y1. Likelihood Analysis of Non-Gaussian Measurement Time Series: http://www.jstor.org/stable/2337586 2. Stochastic volatility: likelihood inference and comparison with ARCH models: https://faculty.mccombs.utexas.edu/carlos.carvalho/teaching/KSB98-SV.pdf
scala> val iters1 = StochasticVolatilityKnots.sampleParametersAr(priorPhi,
| priorMu,
| priorSigma,
| sims)
<console>:13: warning: Unused import
import cats.implicits._
^
<console>:16: warning: Unused import
import breeze.stats.distributions.Gaussian
^
iters1: breeze.stats.distributions.Process[com.github.jonnylaw.dlm.StochVolState] = breeze.stats.distributions.MarkovChain$$anon$1@dcf30df
scala> val iters2 = StochasticVolatility.sampleUni(priorPhi, priorMu, priorSigma, sims)
<console>:13: warning: Unused import
import cats.implicits._
^
<console>:16: warning: Unused import
import breeze.stats.distributions.Gaussian
^
iters2: breeze.stats.distributions.Process[com.github.jonnylaw.dlm.StochVolState] = breeze.stats.distributions.MarkovChain$$anon$1@87ff76c
The figure below shows the diagnostic plots for the posterior distribution of the parameters of the univariate stochastic volatility model.