# Autocorrelation of a MA(1) process with correlated noise sources

I have the following problem. I want to find the expression for the autocorrelation of the following random variable or MA(1) process.

$X_n = \epsilon_n - \epsilon_{n-1}$

Where $\epsilon$ is a stationary correlated noise with finite mean $\mu_\epsilon$ and standard deviation $\sigma_\epsilon$.

So the mean and variance of $X$ is $\mu_X = 0$ and $\sigma_X^2 = 2\sigma^2_\epsilon$ (this was the main mistake in the derivation)

I start from this expression for the autocorrelation.

$\displaystyle R^k_X = \frac{E[(X_n-\mu_X) (X_{n+k}-\mu_X)]}{\sigma_X^2}$

For simplicity I worked with the lag-one, so this is:

$\displaystyle R^1_X = \frac{E[(\epsilon_n-\epsilon_{n-1})(\epsilon_{n+1}-\epsilon_{n})]}{2\sigma^2_\epsilon}$

Expanding and taking the expectancy it gives:

$\displaystyle R^1_X = \frac{E[\epsilon_n\epsilon_{n+1}] + E[\epsilon_n\epsilon_{n-1}]-E[\epsilon_{n+1}\epsilon_{n-1}]-E[\epsilon^2_n]}{2\sigma^2_\epsilon}$

The first two terms are equivalent (correct me if I am wrong), and they are related to the lag one autocovariance of the noise, but they don't have their means subtracted, and the third term is related to the lag two. Using $E[\epsilon_n \epsilon_{n+k}] = Cov[\epsilon_n, \epsilon_{n+k}] + \mu^2_\epsilon$ and adding and subtracting the squared means allow rewriting to:

$\displaystyle R^1_X = R^1_\epsilon -\frac{R^2_\epsilon}{2} + \frac{\mu_\epsilon^2 - E[\epsilon^2_n]}{2\sigma^2_\epsilon}$

Then for the last term, intead of covariance we have variance $E[\epsilon_n^2] = \sigma^2_\epsilon + \mu^2_\epsilon$ and this gives:

$\displaystyle R^1_X = R^1_\epsilon - \frac{R^2_\epsilon}{2} - 0.5$

The problem arises when checking this expression with simulated processes. For the correlated $\epsilon$ noise I used another MA(1) process with uncorrelated normal noise source $\xi$, thus having something like $\epsilon_n = \xi_n+\alpha \xi_{n-1}$, so the lag-two term is theoretically zero.

The following picture shows the deviations from the value of -0.5 for different $\alpha$, it is only valid when $\alpha$ is 0 or 1. I hope there is a problem in the deviation of the formula.

This is the python/numpy code for the simulation, R1 measures the lag-one crosscorrelation:

def R1(x,y):
return mean((x[1:]-x.mean())*(y[:-1]-y.mean()))/std(x)/std(y)

xi = randn(100000)

eps = xi[1:] + alpha*xi[:-1]
X = eps[1:] - eps[:-1]

out = R1(X,X) - R1(eps,eps)

• If it is "only valid when $\alpha$ is 0 or 1," then what sense does it make to plot values for all $\alpha$ in the range $[-2,2]$?? To understand your plot we would need to know how you generated its values.
– whuber
Commented Jun 28, 2017 at 17:29
• What I meant is that the analytic result doesn't agree with the simulation because the curve should be constant. The simulation is straightforward as is explained at the end. Uncorrelated noise -> MA(1) with alpha -> MA(1) with alpha=-1 (original random variable analized) Commented Jun 28, 2017 at 17:49
• "Straightforward," yes. But I believe the problem lies both in how you present the simulation results and in how you implement the simulation. Unfortunately you haven't documented either one adequately so far. You need to do that in order to resolve the discrepancy between simulation and calculation.
– whuber
Commented Jun 28, 2017 at 17:52
• How should I present the simulation results? Commented Jun 28, 2017 at 17:58
• You should explain how the simulations were done, how many iterations were used in each simulation, and summarize each one with the estimated correlation and a standard error for it. That will reveal what data your plot is based on as well is indicating how precise the results are.
– whuber
Commented Jun 28, 2017 at 18:03

### Theory

Adding a constant (say $\lambda$) to $\epsilon_n$ makes no difference, because $$(\epsilon_n + \lambda) - (\epsilon_{n-1} + \lambda) = \epsilon_n - \epsilon_{n-1} = x_n.$$ By choosing $\lambda=-\mu_\epsilon$ you may assume $E[\epsilon_n]=0$.

This makes covariances easy to compute in terms of the lagged correlations of $(\epsilon_n)$ given by $E[\epsilon_{n+i}\epsilon_n]=\rho_i \sigma_\epsilon^2$ for any $i$:

\eqalign{\operatorname{Cov}(x_{n+j}, x_{n}) &= E[x_{n+j}x_{n}]= E[(\epsilon_{n+j}-\epsilon_{n+j-1})(\epsilon_{n}-\epsilon_{n-1})]\\&=E[\epsilon_{n+j}\epsilon_n] - E[\epsilon_{n+j}\epsilon_{n-1}] - E[\epsilon_{n+j-1}\epsilon_n] + E[\epsilon_{n+j-1}\epsilon_{n-1}] \\&= \sigma_\epsilon^2(\rho_j - \rho_{j+1} - \rho_{j-1} + \rho_j) = -\sigma_\epsilon^2(\rho_{j-1} - 2\rho_{j} + \rho_{j+1}).}

A particular case of this is the variance, whose expression can be simplified from the facts $\rho_{-i}=\rho_i$ and $\rho_0=1$:

$$\operatorname{Var}(x_n) = \operatorname{Cov}(x_n,x_n) = -\sigma_\epsilon^2(\rho_{-1} - 2\rho_0 + \rho_1)=2\sigma_\epsilon^2(1-\rho_{1}).$$

Dividing the covariances by the variance yields the lag-$j$ correlation of $(x_n)$, written $\rho_j(x)$, in terms of the correlations of $(\epsilon_n)$,

$$\rho_j(x) = -\frac{\rho_{j-1} - 2\rho_{j} + \rho_{j+1}}{2(1-\rho_1)}.\tag{1}$$

Similar calculations from the definition $\epsilon_n = \xi_n + \alpha \xi_{n-1}$ give

$$\rho_1(\epsilon) = \frac{\alpha}{1+\alpha^2}\tag{2}$$

and $\rho_j(\epsilon) =0$ for $j\gt 1$. Plugging $(2)$ into $1$ and simplifying gives

$$\rho_1(x) = -\frac{1}{2}\frac{1-2\alpha+\alpha^2}{1-\alpha+\alpha^2}.$$

Note that this implies $-2/3 \le \rho_1(x) \le 0$, a range that is fully covered as $\alpha$ ranges from $-1$ through $1$.

### Simulation

$50$ separate simulations of $(\xi_n)$ were used to produce $21$ simulated series $(x_n)$, $n=1,\ldots, 10,000$, for an equally spaced range of values of $\rho_1$ from $-2/3$ through $0$. (The R code is at the end of this message.) The plots of the observed lag-1 autocorrelation coefficients against the values of $\rho_1$ (shown below in red) form a tight envelope around the line of equality (shown below in gray), visibly demonstrating the correctness of these calculations.

We must conclude that the code used to create the plot in the question is erroneous.

n <- 1e4                            # Length of (x) series
r <- seq(-2/3, 0, length.out=21)    # Array of rho_1(x) values to study
alpha <- ifelse(r==-1/2, 0, (1 + r - sqrt(-2*r - 3*r^2)) / (1 + 2*r))
#
# Simulate (x) many times.
#
X <- lapply(1:50, function(i) {
xi <- rnorm(n+2)
sapply(alpha, function(a) {
epsilon <- xi[-1] + xi[-(n+2)] * a
x <- epsilon[-(n+1)] - epsilon[-1]
acf(ts(x), lag.max=1, plot=FALSE)$acf[2,1,1] # acf(ts(epsilon), lag.max=2, plot=FALSE)$acf[2,1,1]
})
})
#
# Prepare to plot the simulated results.
#
rho.1 <- alpha / (1 + alpha^2)    # Lag-1 acf of epsilon
rho.2 <- 0                        # Lag-2 acf of epsilon
rho <- (2*rho.1 - (1 + rho.2)) / (2*(1-rho.1))    # Lag-1 acf of x
#
# Draw a plot window and reference line.
#
plot(range(rho), range(X), type="n", asp=1,
xlab=expression(rho[1]), ylab="Observed",
main="Lag-1 Correlations of X")
abline(c(0,1), col="gray", lwd=2)
#
# Plot the results.
# (Plotting them as points is commented out.)
# invisible(lapply(X, function(y) points(rho, y, pch=16, col="#ff000010")))
invisible(lapply(X, function(y) lines(rho, y, pch=16, col="#ff000020")))

• Thanks for showing a clearer reasoning. I can't prove equation 1. Do you have a reference for its derivation? Commented Jun 28, 2017 at 20:10
• It's the standard formula for correlation! Divide the lag-$j$ covariance by the variance.
– whuber
Commented Jun 28, 2017 at 20:22
• If I had that formula I wouldn't have asked the question! I started from the other standard formula en.wikipedia.org/wiki/Autocorrelation. How is the variance affected by the lag-1? Isn't it the lag-0? Any standard reference to see its derivation would be greatful? Commented Jun 28, 2017 at 20:26
• It's even more basic than that: please refer to the usual definition of correlation, such as at Wikipedia.
– whuber
Commented Jun 28, 2017 at 20:34
• Excellent, thanks. My mistake was when I said that $\sigma_X^2 = 2\sigma^2_\epsilon$, this holds for independent variables! Commented Jun 28, 2017 at 20:45