Title: | Inference, Goodness-of-Fit and Forecast for Univariate Gaussian Hidden Markov Models |
---|---|
Description: | Inference, goodness-of-fit test, and prediction densities and intervals for univariate Gaussian Hidden Markov Models (HMM). The goodness-of-fit is based on a Cramer-von Mises statistic and uses parametric bootstrap to estimate the p-value. The description of the methodology is taken from Chapter 10.2 of Remillard (2013) <doi:10.1201/b14285>. |
Authors: | Bouchra R. Nasri [aut, cre, cph], Bruno N Remillard [aut, ctb, cph] |
Maintainer: | Bouchra R. Nasri <[email protected]> |
License: | GPL (>= 2) |
Version: | 1.1.2 |
Built: | 2025-02-05 16:42:46 UTC |
Source: | https://github.com/cran/GaussianHMM1d |
This function estimates parameters (mu, sigma, Q) of a univariate Hidden Markov Model. It computes also the probability of being in each regime, given the past observations (eta) and the whole series (lambda). The conditional distribution given past observations is applied to obtains pseudo-observations W that should be uniformly distributed under the null hypothesis. A Cramér-von Mises test statistic is then computed.
EstHMM1d(y, reg, max_iter = 10000, eps = 1e-04)
EstHMM1d(y, reg, max_iter = 10000, eps = 1e-04)
y |
(nx1) vector of data |
reg |
number of regimes |
max_iter |
maximum number of iterations of the EM algorithm; suggestion 10 000 |
eps |
precision (stopping criteria); suggestion 0.0001. |
mu |
estimated mean for each regime |
sigma |
stimated standard deviation for each regime |
Q |
(reg x reg) estimated transition matrix |
eta |
(n x reg) probabilities of being in regime k at time t given observations up to time t |
lambda |
(n x reg) probabilities of being in regime k at time t given all observations |
cvm |
Cramér-von Mises statistic for the goodness-of-fit test |
U |
Pseudo-observations that should be uniformly distributed under the null hypothesis of a Gaussian HMM |
LL |
Log-likelihood |
Bouchra R Nasri and Bruno N Rémillard, January 31, 2019
Chapter 10.2 of B. Rémillard (2013). Statistical Methods for Financial Engineering, Chapman and Hall/CRC Financial Mathematics Series, Taylor & Francis.
Q <- matrix(c(0.8, 0.3, 0.2, 0.7),2,2); mu <- c(-0.3 ,0.7) ; sigma <- c(0.15,0.05) data <- Sim.HMM.Gaussian.1d(mu,sigma,Q,eta0=1,100)$x est <- EstHMM1d(data, 2, max_iter=10000, eps=0.0001)
Q <- matrix(c(0.8, 0.3, 0.2, 0.7),2,2); mu <- c(-0.3 ,0.7) ; sigma <- c(0.15,0.05) data <- Sim.HMM.Gaussian.1d(mu,sigma,Q,eta0=1,100)$x est <- EstHMM1d(data, 2, max_iter=10000, eps=0.0001)
This function computes and plots the most likely regime for univariate Gaussian HMM using probabilities of being in regime k at time t given all observations (lambda) and probabilities of being in regime k at time t given observations up to time t (eta).
EstRegime(t, y, lambda, eta)
EstRegime(t, y, lambda, eta)
t |
(nx1) vector of dates (years, ...); if no dates then t=[1:length(y)] |
y |
(nx1) vector of data; |
lambda |
(nxreg) probabilities of being in regime k at time t given all observations; |
eta |
(nxreg) probabilities of being in regime k at time t given observations up to time t; |
A |
Estimated Regime using lambda |
B |
Estimated Regime using eta |
runsA |
Estimated number of runs using lambda |
runsB |
Estimated number of runs using eta |
pA |
Graph for the estimated regime for each observation using lambda |
pB |
Graph for the estimated regime for each observation using eta |
Bouchra R Nasri and Bruno N Rémillard, January 31, 2019
Chapter 10.2 of B. Rémillard (2013). Statistical Methods for Financial Engineering, Chapman and Hall/CRC Financial Mathematics Series, Taylor & Francis.
Q <- matrix(c(0.8, 0.3, 0.2, 0.7),2,2); mu <- c(-0.3 ,0.7) ; sigma <- c(0.15,0.05); data <- Sim.HMM.Gaussian.1d(mu,sigma,Q,eta0=1,100)$x t=c(1:100); est <- EstHMM1d(data, 2) EstRegime(t,data,est$lambda, est$eta)
Q <- matrix(c(0.8, 0.3, 0.2, 0.7),2,2); mu <- c(-0.3 ,0.7) ; sigma <- c(0.15,0.05); data <- Sim.HMM.Gaussian.1d(mu,sigma,Q,eta0=1,100)$x t=c(1:100); est <- EstHMM1d(data, 2) EstRegime(t,data,est$lambda, est$eta)
This function computes the estimated probabilities of the regimes for a Gaussian HMM given new observation after time n. it also computes the associated weight of the Gaussian mixtures that can be used for forecasted density, cdf, or quantile function.
ForecastHMMeta(ynew, mu, sigma, Q, eta)
ForecastHMMeta(ynew, mu, sigma, Q, eta)
ynew |
new observations (mx1); |
mu |
vector of means for each regime (r x 1); |
sigma |
vector of standard deviations for each regime (r x 1); |
Q |
transition probality matrix (r x r); |
eta |
vector of the estimated probability of each regime (r x 1) at time n; |
etanew |
values of the estimated probabilities at times n+1 to n+m, using the new observations |
w |
weights of the mixtures for periods n+1 to n+m |
Bouchra R Nasri and Bruno N Rémillard, January 31, 2019
Chapter 10.2 of B. Rémillard (2013). Statistical Methods for Financial Engineering, Chapman and Hall/CRC Financial Mathematics Series, Taylor & Francis.
mu <- c(-0.3 ,0.7) ; sigma <- c(0.15,0.05); Q <- matrix(c(0.8, 0.3, 0.2, 0.7),2,2); eta <- c(.1,.9); x <- c(0.2,-0.1,0.73) out <- ForecastHMMeta(x,mu,sigma,Q,eta)
mu <- c(-0.3 ,0.7) ; sigma <- c(0.15,0.05); Q <- matrix(c(0.8, 0.3, 0.2, 0.7),2,2); eta <- c(.1,.9); x <- c(0.2,-0.1,0.73) out <- ForecastHMMeta(x,mu,sigma,Q,eta)
This function computes the density function of a Gaussian HMM at time n+k, given observation up to time n.
ForecastHMMPdf(x, mu, sigma, Q, eta, k)
ForecastHMMPdf(x, mu, sigma, Q, eta, k)
x |
points at which the density function is comptuted (mx1); |
mu |
vector of means for each regime (r x 1); |
sigma |
vector of standard deviations for each regime (r x 1); |
Q |
transition probality matrix (r x r); |
eta |
vector of the estimated probability of each regime (r x 1) at time n; |
k |
time of prediction. |
f |
values of the density function at time n+k |
w |
weights of the mixture |
Bouchra R Nasri and Bruno N Rémillard, January 31, 2019
Chapter 10.2 of B. Rémillard (2013). Statistical Methods for Financial Engineering, Chapman and Hall/CRC Financial Mathematics Series, Taylor & Francis.
mu <- c(-0.3 ,0.7) ; sigma <- c(0.15,0.05); Q <- matrix(c(0.8, 0.3, 0.2, 0.7),2,2) ; eta <- c(.9,.1); x <- seq(-1, 1, by = 0.01) out <- ForecastHMMPdf(x,mu,sigma,Q,eta,3) plot(x,out$f,type="l")
mu <- c(-0.3 ,0.7) ; sigma <- c(0.15,0.05); Q <- matrix(c(0.8, 0.3, 0.2, 0.7),2,2) ; eta <- c(.9,.1); x <- seq(-1, 1, by = 0.01) out <- ForecastHMMPdf(x,mu,sigma,Q,eta,3) plot(x,out$f,type="l")
This function computes the distribution function of a mixture of Gaussian univariate distributions
GaussianMixtureCdf(x, mu, sigma, w)
GaussianMixtureCdf(x, mu, sigma, w)
x |
Points at which the distribution function is comptuted (nx1); |
mu |
vector of means for each regime (r x 1); |
sigma |
vector of standard deviations for each regime (r x 1); |
w |
vector of the probability of each regime (r x r). |
F |
values of the distribution function |
Bouchra R Nasri and Bruno N Rémillard, January 31, 2019
mu <- c(-0.3 ,0.7) ; sigma <- c(0.15,0.05); w <-c(0.8, 0.2); x <- seq(-1, 1, by = 0.01) F <- GaussianMixtureCdf(x,mu,sigma,w) plot(x,F,type="l")
mu <- c(-0.3 ,0.7) ; sigma <- c(0.15,0.05); w <-c(0.8, 0.2); x <- seq(-1, 1, by = 0.01) F <- GaussianMixtureCdf(x,mu,sigma,w) plot(x,F,type="l")
This function computes the inverse distribution function of a mixture of Gaussian univariate distributions
GaussianMixtureInv(p, mu, sigma, w)
GaussianMixtureInv(p, mu, sigma, w)
p |
Points in (0,1) at which the distribution function is computed (nx1); |
mu |
vector of means for each regime (r x 1); |
sigma |
vector of standard deviations for each regime (r x 1); |
w |
vector of the probability of each regime (r x 1). |
q |
values of the quantile function |
Bouchra R Nasri and Bruno N Rémillard, January 31, 2019
mu <- c(-0.3 ,0.7) ; sigma <- c(0.15,0.05); w <-c(0.8, 0.2); p <- seq(0.01, 0.99, by = 0.01) q <- GaussianMixtureInv(p,mu,sigma,w) plot(p,q,type="l")
mu <- c(-0.3 ,0.7) ; sigma <- c(0.15,0.05); w <-c(0.8, 0.2); p <- seq(0.01, 0.99, by = 0.01) q <- GaussianMixtureInv(p,mu,sigma,w) plot(p,q,type="l")
This function computes the density function of a mixture of Gaussian univariate distributions
GaussianMixturePdf(x, mu, sigma, w)
GaussianMixturePdf(x, mu, sigma, w)
x |
Points at which the density is comptuted (n x 1); |
mu |
vector of means for each regime (r x 1); |
sigma |
vector of standard deviations for each regime (r x 1); |
w |
vector of the probability of each regime (r x 1). |
f |
Values of the distribution function |
Bouchra R Nasri and Bruno N Rémillard, January 31, 2019
mu <- c(-0.3 ,0.7) ; sigma <- c(0.15,0.05); w <-c(0.8, 0.2); x <- seq(-1, 1, by = 0.01) f <- GaussianMixturePdf(x,mu,sigma,w) plot(x,f,type="l")
mu <- c(-0.3 ,0.7) ; sigma <- c(0.15,0.05); w <-c(0.8, 0.2); x <- seq(-1, 1, by = 0.01) f <- GaussianMixturePdf(x,mu,sigma,w) plot(x,f,type="l")
This function performs a goodness-of-fit test of a Gaussian HMM based on a Cramér-von Mises statistic using parametric bootstrap.
GofHMM1d(y, reg, max_iter = 10000, eps = 1e-04, n_sample = 1000, n_cores)
GofHMM1d(y, reg, max_iter = 10000, eps = 1e-04, n_sample = 1000, n_cores)
y |
(n x 1) data vector |
reg |
number of regimes |
max_iter |
maxmimum number of iterations of the EM algorithm; suggestion 10 000 |
eps |
eps (stopping criteria); suggestion 0.0001 |
n_sample |
number of bootstrap samples; suggestion 1000 |
n_cores |
number of cores to use in the parallel computing |
pvalue |
pvalue of the Cram\'er-von Mises statistic in percent |
mu |
estimated mean for each regime |
sigma |
estimated standard deviation for each regime |
Q |
(reg x reg) estimated transition matrix |
eta |
(n x reg) conditional probabilities of being in regime k at time t given observations up to time t |
lambda |
(n x reg) probabilities of being in regime k at time t given all observations |
cvm |
Cramér-von Mises statistic for the goodness-of-fit test |
W |
Pseudo-observations that should be uniformly distributed under the null hypothesis of a Gaussian HMM |
LL |
Log-likelihood |
Bouchra R Nasri and Bruno N Rémillard, January 31, 2019
Chapter 10.2 of B. Rémillard (2013). Statistical Methods for Financial Engineering, Chapman and Hall/CRC Financial Mathematics Series, Taylor & Francis.
Q <- matrix(c(0.8, 0.3, 0.2, 0.7),2,2); mu <- c(-0.3 ,0.7) ; sigma <- c(0.15,0.05) data <- Sim.HMM.Gaussian.1d(mu,sigma,Q,eta0=1,100)$x gof <- GofHMM1d(data, 2, max_iter=10000, eps=0.0001, n_sample=100,n_cores=2)
Q <- matrix(c(0.8, 0.3, 0.2, 0.7),2,2); mu <- c(-0.3 ,0.7) ; sigma <- c(0.15,0.05) data <- Sim.HMM.Gaussian.1d(mu,sigma,Q,eta0=1,100)$x gof <- GofHMM1d(data, 2, max_iter=10000, eps=0.0001, n_sample=100,n_cores=2)
This function simulates observations from a univariate Gaussian HMM
Sim.HMM.Gaussian.1d(mu, sigma, Q, eta0, n)
Sim.HMM.Gaussian.1d(mu, sigma, Q, eta0, n)
mu |
vector of means for each regime (r x 1); |
sigma |
vector of standard deviations for each regime (r x 1); |
Q |
Transition probality matrix (r x r); |
eta0 |
Initial value for the regime; |
n |
number of simulated observations. |
x |
Simulated Data |
reg |
Markov chain regimes |
Bouchra R Nasri and Bruno N Rémillard, January 31, 2019
Q <- matrix(c(0.8, 0.3, 0.2, 0.7),2,2) ; mu <- c(-0.3 ,0.7) ; sigma <- c(0.15,0.05); sim <- Sim.HMM.Gaussian.1d(mu,sigma,Q,eta0=1,n=100)
Q <- matrix(c(0.8, 0.3, 0.2, 0.7),2,2) ; mu <- c(-0.3 ,0.7) ; sigma <- c(0.15,0.05); sim <- Sim.HMM.Gaussian.1d(mu,sigma,Q,eta0=1,n=100)
This function generates a Markov chain X(1), ..., X(n) with transition matrix Q, starting from a state eta0.
Sim.Markov.Chain(Q, n, eta0)
Sim.Markov.Chain(Q, n, eta0)
Q |
Transition probability matrix (r x r); |
n |
length of series; |
eta0 |
inital value in (1,...,r). |
x |
Simulated Markov chain |
Bouchra R Nasri and Bruno N Rémillard, January 31, 2019
Q <- matrix(c(0.8, 0.3, 0.2, 0.7),2,2) ; sim <- Sim.Markov.Chain(Q,eta0=1,n=100)
Q <- matrix(c(0.8, 0.3, 0.2, 0.7),2,2) ; sim <- Sim.Markov.Chain(Q,eta0=1,n=100)
Generates a univariate regime-switching random walk with Gaussian regimes starting from a given state eta0, using the inverse method from noise u.Can be useful when generating multiple time series.
SimHMMGaussianInv(u, mu, sigma, Q, eta0)
SimHMMGaussianInv(u, mu, sigma, Q, eta0)
u |
series of uniform i.i.d. series (n x 1); |
mu |
vector of means for each regime (r x 1); |
sigma |
vector of standard deviations for each regime (r x 1); |
Q |
Transition probality matrix (r x r); |
eta0 |
Initial value for the regime; |
x |
Simulated Data |
eta |
Probability of regimes |
Bouchra R Nasri and Bruno N Rémillard, January 31, 2019
Nasri & Remillard (2019). Copula-based dynamic models for multivariate time series. JMVA, vol. 172, 107–121.
Q <- matrix(c(0.8, 0.3, 0.2, 0.7),2,2) set.seed(1) u <-runif(250) mu <- c(-0.3 ,0.7) sigma <- c(0.15,0.05); eta0=1 x <- SimHMMGaussianInv(u,mu,sigma,Q,eta0)
Q <- matrix(c(0.8, 0.3, 0.2, 0.7),2,2) set.seed(1) u <-runif(250) mu <- c(-0.3 ,0.7) sigma <- c(0.15,0.05); eta0=1 x <- SimHMMGaussianInv(u,mu,sigma,Q,eta0)
This function computes the Cramér-von Mises statistic Sn for goodness-of-fit of the null hypothesis of a univariate uniform distrubtion over [0,1]
Sn(U)
Sn(U)
U |
vector of pseudos-observations (apprimating uniform variates) |
Sn |
Cramér-von Mises statistic |
Bouchra R Nasri and Bruno N Rémillard, January 31, 2019