# Is the AR(1) process always Gaussian given Gaussian innovations?

I found that $AR(1)$ process $x_t=\phi x_{t-1}+\epsilon_t$ is not always Gaussian given Gaussian innovations $\epsilon_t$. This only happens when the $AR(1)$ model coefficient is very large. This goes against with the theory (If the white noise $\epsilon_t$ is a Gaussian process then $x_t$ is also a Gaussian process). How to explain this then? The following R codes gives that only $363$ out of the 1000 generations of $AR(1$) processes are gaussian process by the Shapiro-Wilk test.

decision=c()
for (k in 1:1000){ # 1000 runs
phi=0.9 # AR(1) coefficient
x=0.5 # value of x_0
X=c() # generated AR(1) process
for (i in 1:1000){  # length of the process is 1000
x=phi*x+rnorm(1) # Gaussian innovations
X=c(X,x)
}
decision=c(decision,shapiro.test(X)$p.value) # Gaussian test } sum(decision>0.05)/1000 # percentage of Gaussian processes  • How did you determine that it is "not Gaussian"? Sep 7, 2013 at 1:18 • @AlecosPapadopoulos I performed the Shapiro test to test for normal. Sep 7, 2013 at 1:22 • Yes, that's correct, the Shapiro Wilk, as stated at the link I gave, assumes independence. Because of the way it works, it is fairly robust to mild dependence, however. Note that up at$\phi=0.9$you're getting toward the region of nonstationarity and there you don't have unconditional normality. One piece of advice: run your simulation of the AR(1) for quite a few steps before recording the values for your sample, or it won't be sufficiently close to stationary. You might want to investigate the use of filter or rollapply (in package zoo) to generate your AR, or even easier arima.sim. Sep 7, 2013 at 1:56 • To clarify further (I hope); by 'toward the region of nonstationarity' I mean that when you get close to$\phi=1$that even though your process is stationary, the sample behaviour in samples can sometimes mimic a nonstationary process for a fair while (and increasingly so as you get closer); you need bigger and bigger samples to see it. Similarly, effects of initial values propagate further and further; with larger$\phi$it takes a longer warmup before your series behaves like a stationary AR. Sep 7, 2013 at 2:08 • You can simply use something like this to generate your AR(1): ts.sim <- arima.sim(list(order = c(1,1,0), ar = 0.9), n = 1000) – Stat Sep 7, 2013 at 2:28 ## 1 Answer By way of summary of the comments, so that this question has an answer (if Alecos would like to present a summary I'll happily delete this and upvote it instead, as long as the question ends up with an answer) -- As Alecos points out, the AR(1) process for$\phi=0.9\$ is definitely Gaussian.

Since it's stationary, observations from it should all have the same Gaussian distribution; the OP is right to expect that it should, since it's undeniably so. As Alecos says, the pen-and-paper result should suggest there's a problem with either the code or the test.

There are in fact several such problems:

(1) The Shapiro Wilk test assumes independence. We don't have it so the test doesn't apply. The OP notes that the test doesn't seem to have problems even with moderately large values of the parameter, and that's not surprising - because of the way the test works, it should be fairly robust to mild dependence. The Shapiro-Wilk doesn't look at the dependence in consecutive values. The most noticeable effect of the dependence on the distribution of the order statistics will be to increase their variance, but the Shapiro Wilk won't notice that at all. There will be a tendency for the tails to wander more from the straight line than with an independent series and eventually that sort of deviation will become detectable.

(2) The code doesn't have any warmup period. The mean and variance of the early values won't correspond to the mean and variance of the AR(1) process (so we don't have independence OR identically distributed values).