07-normal-hmm.Rmd

# Gaussian HMM

Univariate Gaussian HMMs with `TMB` are very similar to the previous Poisson HMM.
The main changes are the distribution parameters being changed from $\bs{\lambda}$ to $\bs{\mu}$ and $\bs{\sigma}$.
In turn, these cause the density function to change from `dpois` to `dnorm` (and from `rpois` to `rnorm` for simulations).

Many functions available in *[functions/utils.R](#utils.R)* have been adapted for the Gaussian case in *[functions/norm_utils.R](#norm_utils.R)*.

We detail below an example of a Gaussian HMM with TMB.

## Dataset
The dataset we show in this example contains log-returns of the SP500 dataset pulled from Yahoo finance.
These log-returns are generated with the following code.
```{r 7sp500-generate, eval = FALSE}
library(tidyverse)
# Loading the sp500 dataset
SP500 <- read_csv("data/S&P500.csv") %>%
  mutate(lag = lag(`Adj Close`),
         lr = log(`Adj Close`/lag)) %>%
  select(lr) %>%
  drop_na()
SP500 <- 100 * SP500$lr
save(SP500, file = "data/SP500.RData")
```

## Likelihood function
The Gaussian negative log-likelihood in TMB can be coded as
```{Rcpp 7norm_hmm.cpp, code=readLines("code/norm_hmm.cpp"), eval=FALSE}
```

and requires the following utility functions.
```{Rcpp 7norm_utils.cpp, code=readLines("functions/norm_utils.cpp"), eval=FALSE}
```

## Estimation
The following R code illustrates how to compute the estimates and standard errors, using the objective function written above.
```{r 7example, cache=TRUE}
# Load packages and utility functions
source("code/packages.R")
source("functions/norm_utils.R")

# Compile and load the objective function for TMB
TMB::compile("code/norm_hmm.cpp")
dyn.load(dynlib("code/norm_hmm"))

# Parameters and covariates
load("data/SP500.RData")
m <- 2
# Means
mu <- c(-2 * mean(SP500), 2 * mean(SP500))
# Standard deviations
sigma <- c(0.5 * sd(SP500), 2 * sd(SP500))
# TPM
gamma <- matrix(c(0.9, 0.1,
                  0.1, 0.9), nrow = 2, byrow = TRUE)

# Parameters & covariates for TMB
working_params <- norm.HMM.pn2pw(m = m, mu = mu, sigma = sigma, gamma = gamma)
TMB_data <- list(x = SP500, m = m)
obj <- MakeADFun(TMB_data, working_params, DLL = "norm_hmm", silent = TRUE)

# Estimation
obj_tmb <- MakeADFun(data = TMB_data,
                     parameters = working_params,
                     DLL = "norm_hmm",
                     silent = TRUE)
mod_tmb <- nlminb(start = obj_tmb$par,
                  objective = obj_tmb$fn,
                  gradient = obj_tmb$gr,
                  hessian = obj_tmb$he)
```
Note that the script above makes a simple educated guess on the true parameters in order to choose initial parameters.

## Results
Estimates can be shown here along with their standard errors, as was done previously in the Poisson case.
As before, the matrix results are displayed in a column-wise format.
This can be seen by comparing estimates with a more readable view.
```{r 7results}
summary(sdreport(obj_tmb), "report")
# Readable estimates
obj_tmb$report()
```

An intuitive interpretation of these results is that on average, market prices either grow slowly (positive low mean) but surely (low variance) as people invest carefully, or the market prices decrease (negative larger mean) as people sell in panic for fear of losing too much (higher variance).
For more details on this, see e.g. \citet{bullaa} and get some information on bull and bear markets.

Note that the TPM is displayed column-wise.