Expected Value of rolling mean of AR(1) process Consider a random variable following an AR(1) model:
$$x_t = \mu+\rho x_{t-1} + \epsilon_t$$
The unconditional expected value of this proces is:
$$E(x)=\mu/(1-\rho)$$
Now take a rolling mean of $x_t$ of n periods including the current:
$$\bar{x} = \frac{1}{n}\sum_{i=0}^{n-1} x_{t-i}$$
Based on Wiki: Sum of normally distributed random variables the expected value is still just the sum of the expected value for each  $x_t$ which is given above, hence the expected value of $\bar{x}_t$ is:
$$E(\bar{x}_t)=1/n \sum_{i=0}^{n-1} E(x_{t-i})=n/n  E(x_t)=E(x_t)$$
However, another approach is to roll back each observation in the summand to a common starting point using (this is with thanks to Aleksej at this post here):
$$x_t = \mu \sum_{i=0}^{k-1} \rho^i + \sum_{i=0}^{k-1} \rho^i \epsilon_{t-i} + \rho^k x_{t-k}$$
Inserting this into the expression for $\bar{x}_t$ yields:
$$\bar{x}_t=1/n\sum_{i=0}^{n-1}(\rho^{n-i}x_{t-n}+\mu\sum_{j=0}^{k-i-1}\rho^j+\sum_{j=0}^{k-i-1}\rho^{j}\epsilon_{t-i-j})$$
$$=1/n(\sum_{i=0}^{n-1}x_{t_n}\rho^{n-i}+ \mu\sum_{i=0}^{n-1}\sum_{j=0}^{n-i-1}\rho^{j}+\sum_{i=0}^{n-1}\sum_{j=0}^{n-i-1}\rho^j\epsilon_{t-i-j})$$
Taking the expection to this yields:
$$E(\bar{x})=1/n(E[x_{x-t}]\sum_{i=0}^{n-1}\rho^{n-i}+\mu\sum_{i=0}^{n-1}\sum_{j=0}^{n-i-1}\rho^{j})$$
$$=\frac{1}{n}(\frac{\mu}{1-\rho}\sum_{i=0}^{n-1}\rho^{n-i}+\mu\sum_{i=0}^{n-1}\sum_{j=0}^{n-i-1}\rho^{j})$$
But this does not coincide with the previous results using the sum of normal distributed variables.
Can someone point out, if I have made a mistake and which method is the correct one?
EDIT: Simulation study
Based on the comment from mlofton, I have conducted a simulation study, and I find that the expected difference and the simulated difference are indeed very close.
Using the formula from the model yields the same as average of the sum of the unconditional means for $$x_t$$
import matplotlib.pyplot as plt
import numpy as np


def simulate_ar_process(initial_value=None, intercept=-0.01, slope=0.5, noise_vol=0.005, obs=100):
    # If initial value is None use the unconditional mean
    if initial_value is None:
        initial_value = intercept / (1-slope)

    errors = np.random.normal(0, noise_vol, obs)
    pd_levels = []

    # Estimate PD-levels
    for err in errors:
        if len(pd_levels) == 0:
            pd_levels.append(initial_value + err)
        else:
            new_obs = intercept + pd_levels[-1] * slope + err
            pd_levels.append(new_obs)

    return pd_levels


def calculate_rolling_mean(pd_levels, look_back=20, burn_in=0):
    mas = []
    for i in range(look_back+burn_in, len(pd_levels)):
        current_range = pd_levels[i-look_back:i]
        mas.append(np.mean(current_range))

    return np.mean(mas)


intercept = -0.01
slope = 0.9
noise_vol = 0.005
obs = 1000
look_back = 20

average_of_rolling_means = []
average_of_x = []
for i in range(100):
    x = simulate_ar_process(intercept=intercept, slope=slope, noise_vol=noise_vol, obs=obs)
    average_of_rolling_means.append(calculate_rolling_mean(x, look_back=look_back))
    average_of_x.append(np.mean(x[20:]))

plt.hist(average_of_x, label="Average of x", alpha=0.3)
plt.hist(average_of_rolling_means, label="Average of rolling mean", alpha=0.3)
plt.legend()

# Calcualte difference
diffs = [x-y for x, y in zip(average_of_x, average_of_rolling_means)]
plt.hist(diffs)
np.average(diffs)

# Theoretical diff
unconditional_mean = intercept / (1-slope)
unconditional_mean_mva = unconditional_mean * np.sum([slope**(look_back-i) for i in range(0, look_back)])
unconditional_mean_mva += intercept * np.sum([slope ** j for i in range(0, look_back) for j in range(0, look_back-i)])
unconditional_mean_mva /= look_back

print("Expected difference: {}".format(unconditional_mean - unconditional_mean_mva))
print("Found difference: {}".format(np.average(diffs)))
#Expected difference: 0.0
#Found difference: 9.324487397407016e-06

However, I cannot see from the formula above, that $$E(\bar{x})=E(x)$$
Can someone point me in the right direction in the deviation of this equality?
 A: There seems to be an error in your calculations. Using a slightly different notation for clarity.
Let the AR(1) process be:
$$X_t = \mu + \phi X_{t-1} + e_t$$
Define mean of $n$ periods be:
$$\mu_n \equiv \frac{1}{n}\sum\limits_{t=1}^n X_t$$
The expectation operator is linear regardless of dependence structure between the random variables. So, your first expression is correct:
$$\mathbb E(\mu_n) = \mathbb E(X_t) = \frac{\mu}{1-\phi}$$
Second approach, for $t>1$:
\begin{align}
X_t &= \mu + \phi X_{t-1} + e_t \\
&= \phi^{t-1}X_1+ \sum\limits_{i=0}^{t-2} \Big(\phi^{i}\mu +\phi^{i}e_{t-i}\Big) \tag{1}
\end{align}
Using $(1)$:
\begin{align}
n\mu_n &= \sum\limits_{t=1}^n\phi^{t-1}X_1 +\sum\limits_{t=2}^n \sum\limits_{i=0}^{t-2} \Big(\phi^{i}\mu +\phi^{i}e_{t-i}\Big) \\
&= \sum\limits_{t=1}^n\phi^{t-1}X_1 +\sum\limits_{t=2}^n \sum\limits_{i=0}^{t-2} \phi^{i}\mu +\sum\limits_{t=2}^n \sum\limits_{i=0}^{t-2}\phi^{i}e_{t-i} \tag{2}
\end{align}
In equation (2) the middle term is non-random and the last term has expected value of $0$. So,
\begin{align}
n\mathbb E(\mu_n) &= \sum\limits_{t=1}^n\phi^{t-1}\mathbb E(X_1)  + \sum\limits_{t=2}^n \sum\limits_{i=0}^{t-2} \phi^{i}\mu \\
&= \frac{\mu}{1-\phi}\sum\limits_{t=1}^n\phi^{t-1}  + \sum\limits_{t=2}^n \sum\limits_{i=0}^{t-2} \phi^{i}\mu \\
&= \frac{\mu}{1-\phi}\sum\limits_{t=1}^n\phi^{t-1} + \sum\limits_{t=2}^n \bigg( \frac{\mu(1-\phi^{t-1})}{1-\phi}\bigg) \\
&= \frac{\mu}{1-\phi} \bigg( \sum\limits_{t=1}^n\phi^{t-1} + \sum\limits_{t=2}^n (1-\phi^{t-1})\bigg) \\
&=n \frac{\mu}{1-\phi}
\end{align}
