Compute share moving between deciles of a stationary AR(1) process

Question

I want to compute the probability $P_{ij}$ to move from decile $i$ in one period to decile $j$ in the next period in the distribution of a stationary AR(1) process

$$Y_t = \rho Y_{t - 1} + \upsilon_t,$$

where $\upsilon_t$ is an i.i.d. shock with distribution $N(0, \sigma^2)$ and $|\rho| < 1$. The stationary distribution of $Y$ is given by $N(0, \tau^2)$, where $\tau \equiv \sigma / \sqrt{1 - \rho^2}$.

Is the following expression what I should be trying to compute? Or is there some more straightforward way to do it?

\begin{align} P_{ij} &= \int_{(i - 1) / 10}^{i / 10} \underbrace{\Pr\biggl(\Phi^{-1}\left(\frac{j - 1}{10}\right) \leq \frac{Y_t}{\tau} < \Phi^{-1}\left(\frac{j}{10}\right)\biggr)}_{\text{probability to move to decile $j$ in time $t$}} \underbrace{\frac{\phi\left(\Phi^{-1}(x)\right)}{0.1}}_{\substack{\text{share of decile $i$} \\[0.5ex] \text{at time $t - 1$}}} d x, \end{align} where $\Phi^{-1}$ is the inverse cumulative distribution function of a standard normal variable and $\phi$ its probability density function.

Answer 1

Letting $Z_t=Y_t/\tau$ where $\tau^2=\sigma^2/(1-\rho^2)$ and $d_i=\Phi^{-1}(\frac i{10})$, $i=1,2,\dots,9$ denote the corresponding deciles, we have, from the definition of conditional probability, $$ \begin{align} P_{ij}&=\Pr(d_{j-1}<Z_t\le d_j|d_{i-1}<Z_{t-1}\le d_i) \\&=\frac{\Pr(d_{j-1}<Z_t\le d_j \cap d_{i-1}<Z_{t-1}\le d_i)} {\Pr(d_{i-1}<Z_{t-1}\le d_i)} \\&=\frac{\Phi_2(d_i,d_j)-\Phi_2(d_i,d_{j-1})-\Phi_2(d_{i-1},d_j)+\Phi_2(d_{i-1},d_{j-1})}{1/10} \end{align} $$ where $\Phi_2$ is the cdf of the bivariate standard normal distribution with correlation $\rho$ (available as the function mvtnorm::pmvnorm in R).

It does not follow that the process is a first order Markov chain with the $P_{ij}$'s as the transition probabilities. This is demonstrated numerically below by computing higher order transisition probabilities based on the joint multivariate normal distribution of $Z_t,Z_{t-1},\dots,Z_{t-n}$ using a straightforward generalisation of the above expression. Indeed, the resulting Markov chain does not appear to be of any finite order, analogous to how an for example an AR(1) process plus white noise becomes a non-Markovian ARMA(1,1) process with partial autocorrelation function only tailing off.

library(mvtnorm)
# Compute the probability of transitioning from state i to j
# If i is a vector, the probabability of transitioning to state j
# given that the previous states are i[1], i[2], ..., i[n] is
# computed.
p <- function(i, j, rho) {
  order <- length(i)
  num <- pmvnorm(lower = qnorm(c(j - 1, i - 1)/10), 
          upper = qnorm(c(j, i)/10), 
          sigma = toeplitz(rho^(0:order)), 
          keepAttr = FALSE)
  den <- pmvnorm(lower = qnorm((i - 1)/10), 
          upper = qnorm(i/10), 
          sigma = toeplitz(rho^(0:(order-1))), 
          keepAttr = FALSE)
  names(den) <- NULL
  num/den
}
# Probability of remaining in state 1
p(1, 1, .5)
#> [1] 0.3240152
# Probability of remaining in state 1 given that the processes
# has been in state 1 for up to five previous time points
p(c(1, 1), 1, .5)
#> [1] 0.3476454
p(c(1, 1, 1), 1, .5)
#> [1] 0.3520791
p(c(1, 1, 1, 1), 1, .5)
#> [1] 0.3533796
p(c(1, 1, 1, 1, 1), 1, .5)
#> [1] 0.353544

^{Created on 2023-12-30 with reprex v2.0.2}