---
title: "STAT 436 / 536 - Lecture 15"
output: pdf_document
---

```{r setup, include=FALSE}
knitr::opts_chunk$set(echo = TRUE)
knitr::opts_chunk$set(warning = FALSE)
knitr::opts_chunk$set(message = FALSE)
knitr::opts_chunk$set(fig.align= 'center')
knitr::opts_chunk$set(fig.height= 6)
knitr::opts_chunk$set(fig.width = 6)
library(tidyverse)
library(gridExtra)
library(readr)
library(ggfortify)
```


### State Space Models

- Recall the innovations form of the state space model that we discussed last time has the following form:
\begin{eqnarray*}
  y_t &=& \widetilde{w}^T \widetilde{x}_{t} + \epsilon_t\\
  \widetilde{x}_t &=& F \widetilde{x}_{t-1} + g \epsilon_t,
\end{eqnarray*}
where only a single error term, $\epsilon_t,$ exists.
\vfill
- Now, contrast this with a state space model below

\vfill
\vfill
\vfill

#### ARMA Models
- The ARMA(p,q) model can be written compactly as
\vfill
\vfill
- Consider the following specification: $\widetilde{w}^T = (1 \; 0 )$, $\widetilde{x}_t = ,$ $v_t = 0$, $F = \begin{bmatrix}
    \alpha_{1}       & 1 \\
    \alpha_{2}       & 0 
\end{bmatrix}$, $R = (1 \; 0)^T$
\begin{eqnarray*}
y_t &=& \left(1 \; 0 \right) \begin{pmatrix}
    x_{1,t} \\
    x_{2,t}
\end{pmatrix}\\
\begin{pmatrix}
    x_{1,t} \\
    x_{2,t}
\end{pmatrix} &=& \begin{bmatrix}
    \alpha_{1}       & 1 \\
    \alpha_{2}       & 0 
\end{bmatrix} \begin{pmatrix}
    x_{1,t-1} \\
    x_{2,t-1}
\end{pmatrix} + \begin{pmatrix}
    1 \\
    0
\end{pmatrix} w_t
\end{eqnarray*}
\vfill
- The observation equation is:
\vfill
- For the state equation, we have:
\vfill
\vfill
\newpage

- The final step is to simplify $x_{1,t}$ by substituting $x_{2,t-1}$. Thus, 
\vfill

\vfill

- Similarly, an ARMA(2,1) model can be specified as: $\widetilde{w}^T = (1 \; 0 )$, $v_t = 0$, $F = \begin{bmatrix}
    \alpha_{1}       & 1 \\
    \alpha_{2}       & 0 
\end{bmatrix}$, $R = (1 \; \beta)^T$
\begin{eqnarray*}
y_t &=& \left(1 \; 0 \right) \begin{pmatrix}
    x_{1,t} \\
    x_{2,t}
\end{pmatrix}\\
\begin{pmatrix}
    x_{1,t} \\
    x_{2,t}
\end{pmatrix} &=& \begin{bmatrix}
    \alpha_{1}       & 1 \\
    \alpha_{2}       & 0 
\end{bmatrix} \begin{pmatrix}
    x_{1,t-1} \\
    x_{2,t-1}
\end{pmatrix} + \begin{pmatrix}
    1 \\
    \beta
\end{pmatrix} w_t
\end{eqnarray*}
Watch for a homework question about this!
\vfill
- Recall (536 students) there were two assumptions for state space models.
\vfill
1. The state vector 
and 
\vfill
2. The observed values are 
\vfill
\vfill

\newpage

- A major benefit of the state-space model framework is that we can easily integrate the stochastic modeling approach for serial correlation (ARMA) directly with additional covariates.

\vfill

- For instance consider the following model:
\begin{eqnarray*}
y_t &=& \left(1 \; 0 \;  \;\;  \right) \begin{pmatrix}
    x_{1,t} \\
    x_{2,t}\\
   ~
\end{pmatrix}\\
\begin{pmatrix}
    x_{1,t} \\
    x_{2,t} \\
    ~
\end{pmatrix} &=& \begin{bmatrix}
    \alpha_{1}       & 1  &0\\
    \alpha_{2}       & 0 &0\\
    0 & 0 & 1
\end{bmatrix} \begin{pmatrix}
    x_{1,t-1} \\
    x_{2,t-1} \\
   ~
\end{pmatrix} + \begin{pmatrix}
    1 \\
    0 \\
    0
\end{pmatrix} w_t
\end{eqnarray*}
\vfill
- What is this model? Specifically, what are $z_t$ and $\beta$? How could we add dynamics to $\beta$?
\vfill
- The `dlm` package contains functions for fitting general models of this type. Standard MCMC algorithms can also be constructed for these models, but we will also have a 536-specific lab that introduces Sequential Monte Carlo (SMC) methods that can be used for this type of model in general.
\vfill
- The `auto.arima()` function in the `forecast` package does have the option to include covariates.
\vfill

\newpage

#### Exercise
- Use the `auto.arima()` function to model the pastry count. Assume that we are looking at this through the lens of explanatory inference, in other words you need to explain to a baker what you see in your model.:
\vfill
1. Change the frequency of the ts object from 1 to 7, what do you see? why are there differences?
\vfill
2. Include the number of coffee drinks purchased on the same day as a covariate, how does the model change?
\vfill
3. How do the model using the ARMA structure compare to a standard linear regression model?
\vfill
4. As in Lab 7, assume our goal was to predict the number pastries consumed tomorrow, how would your model change?
\vfill
5. What if the baker was interested in both the number of pastries *and* the number of cups of coffee to prepare for tomorrow. How would this be modeled?
\vfill
```{r}
library(readr)
library(dplyr)
library(forecast)
bakery_sales <- read_csv('http://math.montana.edu/ahoegh/teaching/timeseries/data/BreadBasket.csv')

pastry_count <- bakery_sales %>% 
  filter(Item %in% c('Pastry','Scandinavian','Medialuna','Muffin','Scone')) %>% 
  group_by(Date) %>% tally() %>% rename(num_pastry = n) 

drink_count <- bakery_sales %>% filter(Item %in% c('Coffee','Tea')) %>% 
  group_by(Date) %>% tally() %>% rename(num_drink = n) 

ggtsdisplay(ts(pastry_count$num_pastry))

```
\vfill


\newpage