Suppose we have a "magic" die $[1-6]$ that never rolls the same number consecutively.
That means you will never find the same number repeated in a row.
Now let's suppose that we roll this die $1000$ times.
How can I find the PDF, expected number of times and variance of getting a specific value?
On
[Assuming the dice is otherwise fair - that is, if you just rolled $x$, the next roll is discrete uniform on $\{1,2,3,4,5,6\}-x$]
Say the rolls are ordered $X_1,X_2,\cdots$ and $x_i\in\{1,2,3,4,5,6\}$ for all $i$. Then the joint pdf is:
$$\mathbb{P}(X_1=x_1, X_2=x_2,\cdots,X_n=x_n)=\frac{1}{6}\left(\frac{1}{5}\right)^{n-1}\prod_{i=1}^{n-1}\mathbb{I}_{\displaystyle\{x_i\neq x_{i+1}\}}$$
and each individual roll has the uniform distribution, i.e.
$$\mathbb{P}(X_n=x_n)=\frac{1}{6}$$
This preserved uniformity is easily proven by induction and conditioning, i.e. working through the relation $\mathbb{P}(X_n)=\mathbb{P}(X_n\vert X_{n-1})\mathbb{P}(X_{n-1})$.
On
Here is some information based upon the Goulden-Jackson Cluster Method.
We consider the set of words $\in \mathcal{V}^{\star}$ of length $n\geq 0$ built from an alphabet $$\mathcal{V}=\{1,2,3,4,5,6\}$$ and the set $B=\{11,22,33,44,55,66\}$ of bad words, which are not allowed to be part of the words we are looking for. We derive a function $f(s)$ with the coefficient of $s^n$ being the number of wanted words of length $n$.
According to the paper (p.7) the generating function $f(s)$ is \begin{align*} f(s)=\frac{1}{1-ds-\text{weight}(\mathcal{C})}\tag{1} \end{align*} with $d=|\mathcal{V}|=6$, the size of the alphabet and with the weight-numerator $\mathcal{C}$ with \begin{align*} \text{weight}(\mathcal{C})=\text{weight}(\mathcal{C}[11])+\cdots+\text{weight}(\mathcal{C}[66]) \end{align*}
We calculate according to the paper \begin{align*} \text{weight}(\mathcal{C}[11])&=-s^2-\text{weight}(\mathcal{C}[11])s\\ &\cdots\\ \text{weight}(\mathcal{C}[66])&=-s^2-\text{weight}(\mathcal{C}[66])s\\ \end{align*} and get \begin{align*} \text{weight}(\mathcal{C}[11])=\cdots=\text{weight}(\mathcal{C}[66])=-\frac{s^2}{1+s} \end{align*}
It follows \begin{align*} f(s)&=\frac{1}{1-ds-\text{weight}(\mathcal{C})}\\ &=\frac{1}{1-6s+\frac{6s^2}{1+s}}\\ &=\frac{1+s}{1-5s}\\ &=(1+s)\sum_{n\geq 0}5^ns^n\\ &=1+6\sum_{n\geq 1}5^{n-1}s^n\\ &=1+6s+30s^2+150s^3+750s^4+3750s^5+\mathcal{O}(s^6)\tag{2} \end{align*}
We conclude, the number of sequences of length $n\geq 1$ from $\{1,2,\ldots,6\}$ with no consecutive numbers is $6\cdot 5^{n-1}$.
Since the number of all sequences of length $n$ is $6^n$, the probability for each valid sequence of length $n\geq 1$ is assuming a fair die \begin{align*} \mathbb{P}(\text{valid sequence of length } n)= \left(\frac{5}{6}\right)^{n-1}\qquad\qquad n\geq 1 \end{align*}
$$ $$
[2016-05-10] Generating function $f(s,u)$ keeping track of a specific digit
With a slight generalisation we can mark all occurrences of a specific digit let's say $6$. We use the variable $u$ to count occurrences of $6$ and derive a generating function $f(s,u)$ \begin{align*} f(s,u)=\sum_{n=0}^{\infty}a_{n,k}u^ks^n \end{align*} with $a_{n,k}$ the number of admissible sequences, i.e. containing no bad words from $B$ of length $n$ with exactly $k$ $6$s.
We start similarly as above but additionally mark each occurrence of $6$ with $u$. \begin{align*} \text{weight}(\mathcal{C}[11])&=-s^2-\text{weight}(\mathcal{C}[11])s\\ &\cdots\\ \text{weight}(\mathcal{C}[66])&=-\color{blue}{u^2}s^2-\text{weight}(\mathcal{C}[66])\color{blue}{u}s\\ \end{align*} and get \begin{align*} \text{weight}(\mathcal{C}[11])&=\cdots=\text{weight}(\mathcal{C}[55])=-\frac{s^2}{1+s}\\ \text{weight}(\mathcal{C}[66])&=-\frac{\color{blue}{u^2}s^2}{1+\color{blue}{u}s} \end{align*}
It follows \begin{align*} f(s,u)&=\frac{1}{1-((d-1)+\color{blue}{u})s-\text{weight}(\mathcal{C},\color{blue}{u})}\\ &=\frac{1}{1-(5+\color{blue}{u})s+\frac{5s^2}{1+s}+\frac{\color{blue}{u^2}s^2}{1+\color{blue}{u}s}}\\ &=\frac{(1+s)(1+us)}{1-4s-5us^2}\tag{3}\\ &=1+(5+u)s+10(2+u)s^2+5(16+13u+u^2)s^3\\ &\qquad+10(32+36u+7u^2)s^4+5(256+368u+121u^2+5u^3)s^5+\mathcal{O}(s^6)\tag{4} \end{align*}
If we take a look e.g. at the coefficient of $s^4$ of the series expansions (2) and (4) we see the number of admissible sequences of length $4$ is $750$ with \begin{align*} &10\left.(32+36u+7u^2)\right|_{u=1}\\ &\qquad=\left.(320+360u+70u^2)\right|_{u=1}\\ &\qquad=750 \end{align*} The coefficient $70$ is the number of admissible sequences of length $4$ containing exactly two $6s$ and similar interpretations hold for the other summands.
Let $a_{k,n}$ denote the number of admissible sequences of length $n$ containing exactly $k$ $6$s. We use the coefficient of operator $[s^n]$ to denote the coefficient of $s^n$ in a series. We obtain from (3) \begin{align*} a_{k,n}=[u^ks^n]f(s,u)=[u^ks^n]\frac{1+(1+u)s+us^2}{1-4s-5us^2}\\ \end{align*}
Since \begin{align*} [u^ks^n]&\frac{1}{1-4s-5us^2}\\ &=[u^ks^n]\sum_{l=0}^\infty s^l(4+5us)^l\\ &=\sum_{l=0}^n[u^ks^{n-l}]\sum_{j=0}^l\binom{l}{j}(5us)^j4^{l-j}\\ &=\sum_{l=0}^n[u^k]\binom{l}{n-l}(5u)^{n-l}4^{2l-n}\\ &=\binom{n-k}{k}5^k4^{n-2k} \end{align*}
we obtain \begin{align*} a_{k,n}&=[u^ks^n]\frac{1+(1+u)s+us^2}{1-4s-5us^2}\\ &=[u^ks^n]\left(1+s+us+us^2\right)\frac{1}{1-4s-5us^2}\\ &=\left([u^ks^n]+[u^ks^{n-1}]+[u^{k-1}s^{n-1}]+[u^{k-1}s^{n-2}]\right)\frac{1}{1-4s-5us^2}\\ &=\binom{n-k}{k}5^k4^{n-2k}+\binom{n-k-1}{k}5^k4^{n-1-2k}\\ &\qquad+\binom{n-k}{k-1}5^{k-1}4^{n-2k+1}+\binom{n-k-1}{k-1}5^{k-1}4^{n-2k}\\ &=5^k4^{n-2k}\left(\binom{n-k}{k}+\frac{1}{4}\binom{n-k-1}{k} +\frac{4}{5}\binom{n-k}{k-1}+\frac{1}{5}\binom{n-k-1}{k-1}\right)\\ &=5^k4^{n-2k}\left(\frac{5}{4}\binom{n-k-1}{k} +2\binom{n-k-1}{k-1}+\frac{4}{5}\binom{n-k-1}{k-2}\right)\tag{5} \end{align*}
A small check of (5) gives e.g. \begin{align*} a_{5,2}&=5^2\cdot2^2\left(\frac{5}{4}\binom{2}{2}+2\binom{2}{1}+\frac{4}{5}\binom{2}{0}\right)\\ &=100\left(\frac{5}{4}+4+\frac{4}{5}\right)\\ &=125+400+80\\ &=605 \end{align*} in accordance with the coefficient of $u^2s^5$ in (4).
We conclude: The number $a_{1000,k}$ of admissible sequences of length $1000$ containing exactly $k$ $6$s is \begin{align*} 5^k4^{1000-2k} \left(\frac{5}{4}\binom{999-k}{k}+2\binom{999-k}{k-1}+\frac{4}{5}\binom{999-k}{k-2}\right) \end{align*}
in accordance with a result of @joriki.
On
let $P_{k,n}$ be the probability to get $k$ times the result "1" in $n$ throws, and let $$P_{k,n} = p_{k,n} + q_{k,n}$$ where $p,q$ are the corresponding probabilities where you add the constraint that you finish with a 1 (for $p$) or anything but 1 (for $q$).
It is easy to get the equations (for $n \ge 1$) : $$5p_{k,n} = q_{k-1,n-1}$$ and $$5q_{k,n} = 5p_{k,n-1} + 4 q_{k,n-1}$$ from which you get that for $n \ge 2$, $$5q_{k,n} = q_{k-1,n-2} + 4q_{k,n-1}$$
Obviously, the $(p_{k,n})$ family satisfies the same linear equation, and since $P$ is the sum of $p$ and $q$, we have for $n \ge2$, that $$5P_{k,n} = P_{k-1,n-2} + 4P_{k,n-1}$$
Together with the data that $P_{0,0} = 1, P_{0,1} = \frac 16, P_{1,1} = \frac 56$ this is enough to determine all those numbers.
Let $f(x,y)$ be the formal power series $$f(x,y)=\sum P_{k,n} x^k y^n$$ From the recurrence on $P$ we immediately have that in $(5-4y-xy^2)f(x,y)$, every coefficient vanish except the coefficients of $x^0y^0, x^0y^1, x^1y^1$, and so that $f(x,y)$ is the power series of the rational fraction $$\frac {30+y+5xy}{6(5-4y-xyy)}$$
We can check that doing the partial evaluation $x=1$ we get $$f(1,y) = \frac {30+6y}{6(5-4y-yy)} = \frac{6(5+y)}{6(5+y)(1-y)} = \frac 1 {1-y} = 1 + y + y^2 + \dots$$so for each $n$, all the probabilities $P_{k,n}$ sum to $1$.
If $E_n$ is the expectation of the number of occurences of 1 in $n$ throws then $$E_n = \sum k P_{k,n}$$so that $$\sum E_n y^n = \frac {df}{dx} (1,y)$$
We compute $$\frac{ df}{dx} = \frac {(5y)(5-4y-xyy)-(-yy)(30+y+5xy)}{6(5-4y-xyy)^2} = \frac {y(5+y)^2}{6(5-4y-xyy)^2}$$
Then evaluating this at $x=1$ we get $$\sum E_ny^n = \frac {y(5+y)^2}{6(5+y)^2(1-y)^2} = \frac y {6(1-y)^2} = \sum \frac n6 y^n$$
Thus $E_n = n/6$ which is as we expected.
To compute the variance we need to compute the expectancy of the square of the number of occurences, that is, $F_n = \sum k^2 P_{k,n}$. Then $$\sum F_n y^n = \frac {d^2f}{dx^2}(1,y) + \frac {df}{dx}(1,y) = \frac{2y^3(5+y)^2}{6(5+y)^3(1-y)^3} + \frac y{6(1-y)^2} = \frac{2y^3+y(5+y)(1-y)}{6(5+y)(1-y)^3} = \frac{5y-4y^2+y^3}{6(5+y)(1-y)^3}$$
$$\sum E_n^2 y^n = \sum \frac 1{36} n^2y^n = \frac {y(1+y)}{36(1-y)^3}$$
Computing their difference gives the variance :
$$\sum V_n y^n = \frac {30y-24y^2+6y^3-y(1+y)(5+y)} {36(5+y)(1-y)^3} = \frac {5y(5-6y+y^2)}{36(5+y)(1-y)^3} = \frac {5y(5-y)}{36(5+y)(1-y)^2} \sim_{y=1} \frac {20}{216(1-y)^2}$$
Doing the partial fraction decomposition gives a formula $$V_n = \frac {50}{1296} + \frac{20}{216}n - \frac {50}{1296}(-1/5)^n$$
I expect that after renormalization, your random variable (number of 1s obtained in $n$ throws) converges in law to a gaussian variable as $n \to \infty$.
On
An impressive array of heavy machinery has been brought to bear on this question. Here's a slightly more pedestrian solution.
The $6\cdot5^{999}$ admissible sequences of rolls are equiprobable, so we need to count how many of them contain $k$ occurrences of a given number, say, $6$.
Let $a_{jn}$ be the number of admissible sequences of length $n$ with exactly $j$ $6$s that begin and end with a non-$6$. Fix a non-$6$ at the beginning, glue a non-$6$ to the right of each of the $j$ $6$s and choose $j$ slots for the resulting glued blocks out of $n-j-1$ objects ($j$ glued blocks and $n-2j-1$ unglued non-$6$s). There are $5^{j+1}4^{n-2j-1}$ options for the non-$6$s, so that makes
$$ a_{jn}=\frac54\cdot4^n\left(\frac5{16}\right)^j\binom{n-j-1}j $$
sequences.
A sequence with $k$ $6$s can begin and end with non-$6$s, for a contribution of $a_{k,1000}$, begin with a $6$ and end with a non-$6$ or vice versa, for a contribution of $2a_{k-1,999}$, or begin and end with a $6$, for a contribution of $a_{k-2,998}$, for a total of
$$ a_{k,1000}+2a_{k-1,999}+a_{k-2,998}=4^{1000}\left(\frac5{16}\right)^k\left(\frac54\binom{999-k}k+2\binom{999-k}{k-1}+\frac45\binom{999-k}{k-2}\right)\;. $$
Dividing by the total number of admissible sequences yields the probability of rolling $k$ $6$s:
$$ \frac56\left(\frac45\right)^{1000}\left(\frac5{16}\right)^k\left(\frac54\binom{999-k}k+2\binom{999-k}{k-1}+\frac45\binom{999-k}{k-2}\right)\;. $$
We could use this to calculate the expectation and variance, but this is more easily done using the linearity of expectation.
By symmetry, the expected number of $6$s rolled is $\frac{1000}6=\frac{500}3=166.\overline6$.
With $X_i$ denoting the indicator variable for the $i$-th roll being a $6$, the variance is
\begin{align} \operatorname{Var}\left(\sum_iX_i\right) &= \mathbb E\left(\left(\sum_iX_i\right)^2\right)-\mathbb E\left(\sum_iX_i\right)^2 \\ &= 1000\left(\frac16-\frac1{36}\right)+2\sum_{i\lt j}\left(\textsf{Pr}(X_i=1\land X_j=1)-\textsf{Pr}(X_i=1)\textsf{Pr}(X_j=1)\right) \\ &= 1000\left(\frac16-\frac1{36}\right)+2\sum_{i\lt j}\textsf{Pr}(X_i=1)\left(\textsf{Pr}(X_j=1\mid X_i=1)-\textsf{Pr}(X_j=1)\right)\;. \end{align}
The conditional probability $p_{j-i}=\textsf{Pr}(X_j=1\mid X_i=1)$ depends only on $j-i$; it satisfies $p_{k+1}=\frac15(1-p_k)$ and has initial value $p_1=0$, with solution $p_k=\frac{1-\left(-\frac15\right)^{k-1}}6$. Thus
\begin{align} \operatorname{Var}\left(\sum_iX_i\right) &= 1000\left(\frac16-\frac1{36}\right)-\frac1{18}\sum_{k=1}^{999}(1000-k)\left(-\frac15\right)^{k-1} \\ &= 1000\left(\frac16-\frac1{36}\right)-\frac1{18}\cdot\frac{25}{36}\left(\frac65\cdot1000-1+\left(-\frac15\right)^{1000}\right) \\ &=\frac{60025-5^{-998}}{648} \\ &\approx\frac{60025}{648} \\ &\approx92.63\;, \end{align}
in agreement with mercio's result.
On
I just did a simple simulation for $1,000$ valid rolls $10,000$ times (to get a good average) and got the following average (absolute value) difference from the expected results of a fair die (which is $166.667$ of each of the $6$ numbers)
$1$: $~~7.7098$
$2$: $~~7.6851$
$3$: $~~7.7798$
$4$: $~~7.5926$
$5$: $~~7.6609$
$6$: $~~7.5286$
Interestingly, if I do the same simulation but for a fair die, I get:
$1$: $~~9.3069$
$2$: $~~9.4111$
$3$: $~~9.4203$
$4$: $~~9.4153$
$5$: $~~9.4366$
$6$: $~~9.3795$
I do not know why this is.
Getting a precise answer is possible, but not very enlightening (it would be very complicated, with a thousand terms to compute), so I'll focus on two approximate results: the law of large numbers, and the central limit theorem. With $1000$ throws on such a simple system, they are bound to be quite accurate.
To begin with, let's say that we throw a standard die a large number $N$ of times.
The standard die
Let $(X_n)_{n \geq 0}$ be the successive results of the die, which take their values in $\{1, \ldots, 6\}$. Then the $X_n$ are independent, and identically distributed (each has a probability of $1/6$ to yield any given number).
Let $f : \{1, \ldots, 6\} \to \mathbb{R}$ an observable. Then the sequence $(f(X_n))_{n \geq 0}$ is also a sequence of independent and identically distributed random variables. Put $S_N f := \sum_{k=0}^{N-1} f(X_k)$. For instance, if we take:
$f(x) = 1_{x=4}$: then $S_N f$ is the number of times we get the number $4$ in the first $N$ throws.
$f(x) = x$: then $S_N f$ is the sum of the first $N$ throws.
Since the expectation is linear,
$$\mathbb{E} (S_N f) = \mathbb{E} \left(\sum_{k=0}^{N-1} f(X_k)\right) = \sum_{k=0}^{N-1} \mathbb{E} \left(f(X_k)\right) = N\mathbb{E} \left(f(X_0)\right).$$
But that only deals with the average value of $S_N f$. Since $f$ is bounded, and thus integrable, the (weak) law of large numbers asserts that, in distribution,
$$\lim_{N \to + \infty} \frac{S_N f}{N} =\mathbb{E} \left(f(X_0)\right).$$
This gives you a first order result. This is not completely satisfactory; for instance, if $\mathbb{E} \left(f(X_0)\right) = 0$, this law tells us that the sum grows sub-linearly. We may want something more precise, which gives us a non degenerate limit distribution. Since $f$ is bounded, it is square-integrable, and we may apply the central limit theorem:
$$\lim_{N \to + \infty} \frac{S_N f-N\mathbb{E} \left(f(X_0)\right)}{\sqrt{N}} = \sigma(f) \mathcal{N},$$
where the convergence is in distribution, $\mathcal{N}$ is a standard centered Gaussian, $\overline{f} = f-\mathbb{E} (f)$, and:
$$\sigma(f)^2 = \mathbb{E} (\overline{f}^2) = \frac{1}{6} \sum_{k=1}^6 \overline{f}(k)^2,$$
For instance, if $f(x) = 1_{x=4}$, we get:
$$\lim_{N \to + \infty} \frac{S_N f}{N} = \frac{1}{6},$$
so roughly $1/6$th of the throws yield a $4$, and:
$$\lim_{N \to + \infty} \frac{S_N f-N/6}{\sqrt{N}} = \frac{\sqrt{5}}{6} \mathcal{N},$$
which tells us for instance that the number of $4$, with probability $0.95$, will be in the interval $[N/6-\sqrt{5N}/3, N/6+\sqrt{5N}/3]$. For $N = 1000$ throws, this is the interval $[143,190]$.
The magical die
I'll now show how these results apply to the magical die. However, I won't prove what I assert, because it tends to be quite high level (I do this in a graduate course, typically). I'll assume that the die is fair: for the first throw, the $6$ possible outcomes are equally likely, and for each throw past the first, the $5$ possible outcomes are also equally likely. The first remark is that, for any given throw,
$$\mathbb{P} (X_n = 1) = \ldots = \mathbb{P} (X_n = 6) = \frac{1}{6}.$$
This is because, if you relabel the faces of the die, it doesn't change the law of the sequence of throws, so at each throw, all the outcomes are equally likely.
The sequence $(X_n)_{n \geq 0}$ is then a stationary Markov chain: for all $n$, the throw $X_{n+1}$ depends only on $X_n$. I may encode it into the following matrix:
$$P:= \left(\begin{matrix} 0 & 1/5 & 1/5 & 1/5 & 1/5 & 1/5 \\ 1/5 & 0 & 1/5 & 1/5 & 1/5 & 1/5 \\ 1/5 & 1/5 & 0 & 1/5 & 1/5 & 1/5 \\ 1/5 & 1/5 & 1/5 & 0 & 1/5 & 1/5 \\ 1/5 & 1/5 & 1/5 & 1/5 & 0 & 1/5 \\ 1/5 & 1/5 & 1/5 & 1/5 & 1/5 & 0 \end{matrix}\right)$$
The coefficient $P_{ij}$ is the probability that, given that the die is currently on $i$, the next throw yields $j$. So $P_{ii} = 0$ for all $i$, and $P_{ij} = 1/5$ otherwise.
This Markov chain is finite-state (there are six possible states), transitive (starting from any state $i$, I can get to any state $j$ in two steps), and aperiodic (given any two states $i$ and $j$, I can go from $i$ to $j$ in $n$ steps, for all $n \geq 2$ - in other words, if you give me two possible outcomes of the die, I can see these possible outcomes $n$ throws apart, as long as $n \geq 2).
Let $f : \{1, \ldots, 6\} \to \mathbb{R}$ be an observation. Then the sequence $(f(X_n))_{n \geq 0}$ is still identically distributed (since the $X_n$ are identically distributed), but these random variables are no longer independent. Knowing $f(X_n)$ tends to constrain severely $f(X_{n+1})$: if $f$ is injective, then $f(X_{n+1}) \neq f(X_n)$. However, what will save us is that the dependence decays very fast: knowing $f(X_n)$ gives a lot of information on $f(X_{n+1})$, but very little on, say, $f(X_{n+10})$. And this influence decay exponentially fast. In other words, if the first throw is $1$, you can say reliably that the second throw won't be $1$, but there is not much you can say on the tenth or twentieth throw.
Given an observable $f$, note that we still have:
$$\mathbb{E} (S_N f) = \mathbb{E} \left(\sum_{k=0}^{N-1} f(X_k)\right) = N\mathbb{E} \left(f(X_0)\right),$$
because this argument only uses the stationarity of the process. We also still have a law of large numbers, in distribution,
$$\lim_{N \to + \infty} \frac{S_N f}{N} =\mathbb{E} \left(f(X_0)\right).$$
So the law of large number is still as simple, and is a consequence of Birkhoff's theorem and the transitivity of the Markov chain. Now, we also have a central limit theorem, but it is more complicated. Indeed,
$$\lim_{N \to + \infty} \frac{S_N f-N\mathbb{E} \left(f(X_0)\right)}{\sqrt{N}} = \sigma_{GK} (f) \mathcal{N},$$
where the convergence is in distribution, $\mathcal{N}$ is a standard centered Gaussian, $\overline{f} = f-\mathbb{E} (f)$ and:
$$\sigma_{GK} (f)^2 = \mathbb{E} (\overline{f}^2) + 2 \sum_{n=1}^{+ \infty} \mathbb{E} (\overline{f} (X_0) \overline{f}(X_n)).$$
the main change is the formula for the variance, which is given by Green-Kubo's formula.
Now, how do we compute this thing? An observable $f$ is essentially the same thing as an element of $\mathbb{R}^6$; for instance, $1_{x=4} = (0, 0, 0, 0, 1, 0)$. Now, one can prove that, for any functions $f$ and $g$,
$$\mathbb{E} (f (X_0) g(X_n)) = \mathbb{E} (f (X_0) (P^n) g (X_0)) = \frac{1}{6} \sum_{k=1}^6 f(k) (P^n) g (k),$$
and thanksfully, $P^n$ is easy to compute. If we write $A := P+I/5$ the $6 \times 6$ matrix whose coefficients are all $1/5$, then $A^k = (6/5)^{k-1} A$ for $k \geq 1$, so that:
$$P^n = \left( A-\frac{I}{5} \right)^n = \sum_{k=0}^n \binom{n}{k} A^k \left( - \frac{1}{5} \right)^{n-k} = \frac{5}{6} \sum_{k=0}^n \binom{n}{k} \left(\frac{6}{5} \right)^k \left( - \frac{1}{5} \right)^{n-k} A - \frac{5}{6}\left( - \frac{1}{5} \right)^n A+ \left( - \frac{1}{5} \right)^n I = A - \frac{5}{6}\left( - \frac{1}{5} \right)^n A+ \left( - \frac{1}{5} \right)^n I.$$
In addition, if $\mathbb{E} (f (X_0)) = 0$, then $Af = 0$. This is the case for $\overline{f}$. Hence,
$$P^n \overline{f} = \left( - \frac{1}{5} \right)^n \overline {f},$$
so that:
$$\sigma_{GK} (f)^2 = \mathbb{E} (\overline{f}^2) + 2 \sum_{n=1}^{+ \infty} \mathbb{E} (\overline{f} (X_0) (P^n \overline{f})(X_0)) = \mathbb{E} (\overline{f}^2) + 2 \sum_{n=1}^{+ \infty} \left( - \frac{1}{5} \right)^n \mathbb{E} (\overline{f}^2) = \left( 1+2 \sum_{n=1}^{+ \infty} \left( - \frac{1}{5} \right)^n\right) \mathbb{E} (\overline{f}^2) = \frac{2}{3} \sigma (f)^2.$$
For instance, if $f(x) = 1_{x=4}$, we get:
$$\lim_{N \to + \infty} \frac{S_N f}{N} = \frac{1}{6},$$
so roughly $1/6$th of the throws yield a $4$, and:
$$\lim_{N \to + \infty} \frac{S_N f-N/6}{\sqrt{N}} = \frac{\sqrt{10}}{6 \sqrt{3}} \mathcal{N},$$
which tells us for instance that the number of $4$, with probability $0.95$, will be in the interval $[N/6-\sqrt{10N}/(3\sqrt{3}), N/6+\sqrt{10N}/(3\sqrt{3})]$. For $N = 1000$ throws, this is the interval $[147,186]$.
What happens is that, if $f(X_n)$ is large, then since $X_{n+1} \neq X_n$, $f(X_{n+1})$ tends to be smaller. Large increments are immediately followed by smaller increments, and small increments are immediatly followed by larger increments. This phenomenon tends to reduce the spread of $S_N f$, which explains that the results are slighly more concentrated with the magical die than with a standard die.