# When is ARMA(p,q) strictly (strongly) stationary and what is its unconditional distribution?

It seems that the discussions on ARMA are always focused on weak (second-order) stationary, but what about strong stationary? What are the conditions on the coefficients for it to be strictly stationary and what is the corresponding unconditional distribution?

Furthermore, is there a general procedure for deriving its unconditional distribution? I have tried deriving it in terms of the characteristic function, but not sure if it is correct. For example, for $AR(1)$ given by

$$x_t = ax_{t-1} + z_t$$

we get

$$\psi_X(s) := \mathbb{E}[\exp(i s x_t)] = \mathbb{E}[\exp(i a s x_{t-1})]\mathbb{E}[\exp(i s z_t)] = \psi_X(as)\exp(-0.5\sigma^2s^2)$$

assuming $z_t \sim N(0,\sigma^2)$, but it is not clear to me how to proceed from here. One can guess a solution and check it by substitution, so, using the knowledge of the variance of a weak stationary $AR(1)$, we may guess that

$$\psi_X(s) = \exp\left(-\frac{1}{2}\frac{\sigma^2}{1-a^2}s^2\right)$$

which gives

$$\psi_X(as) = \exp\left(-\frac{1}{2}\frac{\sigma^2 a^2}{1-a^2}s^2\right)$$

and we find that

$$\frac{\psi_X(s)}{\psi_X(as)} = \exp\left(-\frac{1}{2}\frac{\sigma^2}{1-a^2}s^2\right)\exp\left(\frac{1}{2}\frac{\sigma^2 a^2}{1-a^2}s^2\right) = \exp\left(-\frac{1}{2}\sigma^2s^2\right)$$

which satisfies the equation. So, it would seem that weak stationary is a sufficient condition for $AR(1)$ to be also strongly stationary (probably this only holds for this special case of Gaussian noise, and I guess there already must be a theorem that proves this). But how would one proceed in a general case, where noise is not restricted to Gaussian and without guessing? And what about more complex stationary processes, say, like $GARCH$? For example, $ARCH(1)$ is

$$x_t = \sigma_t z_t$$ $$\sigma_t^2 = a_0 + a_1 x_{t-1}^2$$

which we can re-write as

$$x_t^2 = (a_0 + a_1 x_{t-1}^2) z_t^2 .$$

which then gives

$$\psi(s) := \mathbb{E}[\exp(i s x_t^2)] = \mathbb{E}[\exp(i s a_0 z_t^2)] \mathbb{E}[\exp(i s a_1 x_{t-1}^2 z_t^2)] .$$

If we again assume that $z_t \sim N(0,1)$, then $z_t^2 \sim \chi^2(1)$ and $c z_t^2 \sim \Gamma(0.5,2c)$, then we have

$$\mathbb{E}[\exp(i s a_0 z_t^2)] = (1 - 2 a_0 i s)^{-0.5}$$

and, conditional on $x_{t-1}$,

$$\mathbb{E}_{t-1}[\exp(i s a_1 x_{t-1}^2 z_t^2) | x_{t-1}] = (1 - 2 a_1 x_{t-1}^2 i s)^{-0.5}$$

so that

$$\psi(s) = (1 - 2 a_0 i s)^{-0.5} \mathbb{E}[(1 - 2 a_1 x_{t-1}^2 i s)^{-0.5}]$$

and guessing the solution is not so easy now.

Add 1 It might be possible to use fractional integral (using Cauchy formula) to write the above as

$$\psi(s) = (1 - 2 a_0 i s)^{-0.5} (2 a_1 i s)^{-0.5} \mathbb{E}\left[\left(\frac{1}{2 a_1 i s} - x^2\right)^{-0.5}\right] = (1 - 2 a_0 i s)^{-0.5} (2 a_1 i s)^{-0.5}\Gamma(0.5)(J^{0.5}f)\left(\frac{1}{2 a_1 i s}\right)$$

where $f$ is the pdf of $x^2$. So we get a fractional integral equation (although I am not sure about the accuracy of the above).

You are making the problem much harder than needed. A necessary condition that the process ${X_t}$ be stationary is that $X_t$ must have the same distribution as $X_{t-1}$. In particular, $X_t$ must have the same mean and variance as $X_{t-1}$. Now, in order for $E[X_t] = E[X_{t-1}]$ to hold, it must be that $$E[X_t] = aE[X_{t-1}] + E[Z_t] = aE[X_{t-1}] \implies E[X_t] = 0 ~\forall t.$$ Similarly, in order for $\operatorname{var}(X_t) = \operatorname{var}(X_{t-1})$ to hold, it must be that $$\operatorname{var}(X_t) = a^2 \operatorname{var}(X_{t-1}) + \operatorname{var}(Z_t) = a^2 \operatorname{var}(X_t) + \sigma^2 \implies \operatorname{var}(X_t) = \frac{\sigma^2}{1-a^2}.$$ Now, the distribution of $X_t$ is the convolution of the distributions of $aX_{t-1}$ and $Z_t$, and if $X_{t-1}$ and $Z_t$ are (independent) normal random variables, so is $X_t$ a normal random variable.

With regard to strict stationarity when the $X$'s and $Z$'s are normal random variables, see the last two paragraphs of this answer of mine over on dsp.SE.

• Ok, I deleted my own comments accordingly. – Jarle Tufto Aug 15 '17 at 21:32
• It is only Gaussian for which a weak stationary implies a strong one or is it true for all elliptical distributions (or some other class) and why? – Confounded Aug 16 '17 at 8:19

Stationary distribution for the general ARMA model: The condition of strong stationarity is not sufficient to determine the distribution of the time-series process, even when restricted by the recursive ARMA equation. However, for practical purposes, when we want a strongly stationary ARMA process with normally distributed errors, we assume that the process has a multivariate normal distribution. For this process, any vector of observable values $$\mathbf{X} = (X_1,...,X_T)$$ has the distribution:

$$\mathbf{X} \sim \text{N} ( \mu \mathbf{1}, \mathbf{\Sigma} ) \quad \quad \quad \mathbf{\Sigma} = \begin{bmatrix} \gamma(0) & \gamma(1) & \cdots & \gamma(T-1) \\ \gamma(1) & \gamma(0) & \cdots & \gamma(T-2) \\ \vdots & \vdots & \ddots & \vdots \\ \gamma(T-1) & \gamma(T-2) & \cdots & \gamma(0) \\ \end{bmatrix},$$

where $$\mu \in \mathbb{R}$$ is the mean of the process, and the auto-covariance function $$\gamma (k) \equiv \mathbb{C}(X_t, X_{t+k})$$ is found from the parameters of the ARMA process via standard methods. (Strict stationarity also requires restrictions on the parameters of the model, ensuring that the roots of the autoregressive characteristic polynomial are outside the unit circle.) This joint distribution corresponds to the case where we have an ARMA process with normally distributed error terms, and strict stationarity. Note that the auto-covariance function will be a function of the auto-regressive parameter $$\boldsymbol{\phi}$$ and the moving-average parameter $$\boldsymbol{\theta}$$.

Stationary distribution for the AR(1) model: In the special case of an $$\text{AR}(1)$$ process with auto-regression parameter $$|\phi| < 1$$, the auto-covariance function for the process is:

$$\gamma (k) \equiv \mathbb{C}(X_t, X_{t+k}) = \frac{\phi^k}{1-\phi^2} \cdot \sigma^2.$$