Variance of number of runs of consecutive heads in $n$ biased coin flips Suppose we have a sequence of $n$ coin flips of a biased coin with probability $p$ of being heads. Let $X_n$ be the random variable that records the number of runs of consecutive heads, e.g. THHTHTHH has 3 runs of heads. By linearity of expectation, $E(X_n)$ is the sum of $n$ probabilities of getting a run of heads starting at coin flip $j$ for $1 \leq j \leq n$. For $j = 1$ this probability is $p$, and for $j > 1$ we must have flip $j-1$ is tails and flip $j$ is heads, which has probability $p(1-p)$. So $E(X_n) = p + (n-1)p(1-p)$. What about Var$(X_n)$? Is there any way to compute a formula for it? Or if not, at least get asymptotics for it?
 A: Following the idea of your derivation: let $t_i$ be $1$ if a run of heads starts at coin number $i$  ($i=1 \cdots n$), $0$ otherwise. Then $X_n = \sum_{i=1}^n t_i$ and
$$ 
P(t_i=1)=E(t_i)=
\begin{cases}
p  & \text{if $i=1$} \\
p(1-p) & \text{otherwise}  \\
\end{cases}
$$
The vars $t_i$ are not independent, but $E(X_n)=\sum E(t_i) = p + (n-1)p (1-p)$ as you already found. Now let's compute
$$E(X^2_n)= E\left[\left(\sum t_i \right)^2\right]=\sum_{i=1}^n\sum_{j=1}^n E(t_i t_j)$$
The terms of the sum can be grouped in:
A) $i=j$: (with $n$ terms)  $E(t_i t_j)=E(t_i^2)=E(t_i)$. The sum over these terms is $p + (n-1)p (1-p)$.
B) $|i-j|=1$: (with $2 \times (n-1)$ terms). Here one of the two must be zero,so  $E(t_i t_j)=0$
In the remainder, the variables $t_i,t_j$ are independent. Hence
C) $|i-j|>1$ and $i=1$ or $j=1$ : (with $2\times(n-2)$ terms) $E(t_i t_j) = p^2(1-p)$
D) $|i-j|>1$ and $i>1$ and $j>1$ : (with $n^2-5n +6$ terms) $E(t_i t_j) = p^2(1-p)^2$
Finally
$$E(X_n^2)=p+\left( n-1\right) \,\left( p-1\right) \,p+ 2\,\left( n-2\right) \,\left( 1-p\right) \,{p}^{2}+\left( {n}^{2}-5\,n+6\right) \,{\left( 1-p\right) }^{2}\,{p}^{2}$$
$$Var(X_n) = E(X_n^2) - E(X_n)^2 = \left( 1-p\right) \,p\,\left( 3\,n\,{p}^{2}-5\,{p}^{2}- 3 \, n\,p+3 p+n\right)=\\
=\left( 1-p\right) \,p\,\left[  n\,\left( 3\,{p}^{2}-3\,p+1\right) -p\,\left( 5\,p-3\right)  \right]$$
As pointed out by Did in a comment, this is valid for $n\ge 4$; if $n<4$ the group D vanishes; if $n\le 2$, only case A survives. Then
$$Var(X_1)=p(1-p) \hspace{9pt} Var(X_2)=p(1-p)^2(2-p)\hspace{9pt} Var(X_3)=p(1-p)^2(4p^2-6p+3)$$
A: Here is a modest  contribution. Classify binary sequences according to
the difference  in the  number of  runs of heads  and tails.   This is
either  zero, one, or  minus one.  We will  represent these  with four
bivariate probability generating functions with $z$ counting the total
sequence length and $u$ counting the number of runs of heads.
For the same number of  runs we get the following generating functions
(depending on whether we start with heads or tails):
$$G_1(z, u) = 
\sum_{k=0}^\infty u^k
\left(\frac{pz}{1-pz}\frac{(1-p)z}{1-(1-p)z}\right)^k
= {\frac {1+ \left( -{p}^{2}+p \right) {z}^{2}-z}{1+p \left( u-1
 \right)  \left( -1+p \right) {z}^{2}-z}}$$
and
$$G_2(z, u) = 
\sum_{k=0}^\infty u^k
\left(\frac{(1-p)z}{1-(1-p)z}\frac{pz}{1-pz}\right)^k
= {\frac {1+ \left( -{p}^{2}+p \right) {z}^{2}-z}{1+p \left( u-1
 \right)  \left( -1+p \right) {z}^{2}-z}}.$$
The symmetry here is obvious. If there is one more run of heads than tails we get
$$G_3(z, u) =
u\frac{pz}{1-pz}
\sum_{k=0}^\infty u^k
\left(\frac{pz}{1-pz}\frac{(1-p)z}{1-(1-p)z}\right)^k
\\=
{\frac {pz \left( -up{z}^{2} \left( -1+p \right) -uz+u \right) }{
 \left( -pz+1 \right)  \left( 1+p \left( u-1 \right)  \left( -1+p
 \right) {z}^{2}-z \right) }}.$$
Finally if there is one more run of tails than heads we get
$$G_4(z, u) = 
\frac{(1-p)z}{1-(1-p)z}
\sum_{k=0}^\infty u^k
\left(\frac{(1-p)z}{1-(1-p)z}\frac{pz}{1-pz}\right)^k
\\ =
{\frac { \left( 1-p \right) z \left( 1+ \left( -{p}^{2}+p \right) {z}^
{2}-z \right) }{ \left( 1- \left( 1-p \right) z \right)  \left( 1+p
 \left( u-1 \right)  \left( -1+p \right) {z}^{2}-z \right) }}.$$
The generating function for all cases is given by
$$G(z, u) =
G_1(z,u) + G_2(z,u) + G_3(z, u) + G_4(z, u)\\
= {\frac {{z}^{2}{p}^{2}u-{z}^{2}{p}^{2}-{z}^{2}pu+upz+{z}^{2}p-pz-z+2}{
{z}^{2}{p}^{2}u-{z}^{2}{p}^{2}-{z}^{2}pu+{z}^{2}p-z+1}}.$$
Observe that
$$G(z, 1) = \frac{2-z}{1-z} = \frac{2}{1-z} - \frac{z}{1-z}$$
so that for $n>1$ we have
$$[z^n] G(z, 1) = 2 - 1 = 1,$$
which confirms that $[z^n] G(z,u)$ is indeed a probability generating function (probabilities sum to one).
Now we have for the expectation that
$$\mathrm{E}[X_n] = [z^n] \left. \frac{d}{du} G(z, u) \right|_{u=1}.$$
This works out to
$$[z^n] \frac{pz (1-pz)}{(1-z)^2}
= p [z^{n-1}] \frac{1}{(1-z)^2}
- p^2 [z^{n-2}] \frac{1}{(1-z)^2}
\\= p n - p^2 (n-1)
= p + (n-1)\times p(1-p),$$
confirming the result quoted in the question text.
For the variance compute the required factorial moment, which is
$$\mathrm{E}[X_n (X_n-1)] = [z^n] 
\left. \left(\frac{d}{du}\right)^2 G(z, u) \right|_{u=1}.$$
This works out to (formula is valid for $n\ge 4$ as explained in the definitive answer)
$$[z^n] \frac{2p^2 z^3 (1-p) (1-pz)}{(1-z)^3}\\
= 2p^2 (1-p) [z^{n-3}] \frac{1}{(1-z)^3}
- 2p^3 (1-p) [z^{n-4}] \frac{1}{(1-z)^3} \\=
2p^2 (1-p) {n-3+2\choose 2} 
- 2p^3 (1-p) {n-4+2\choose 2}
\\=  \left( -1+p \right) {p}^{2} \left( n-2 \right)  \left( pn-n-3\,p+1
 \right)$$
To conclude recall that
$$\mathrm{Var}(X) = \mathrm{E}(X^2) - \mathrm{E}(X)^2
= \mathrm{E}(X(X-1)) + \mathrm{E}(X) - \mathrm{E}(X)^2$$
to get
$$\mathrm{Var}(X_n)
=- \left( -1+p \right) p \left( 3\,{p}^{2}n-3\,pn-5\,{p}^{2}+n+3\,p
 \right)$$
which is
$$(1-p)p \left((1-3p(1-p))\times n + 3p-5p^2\right).$$
Sequence A001792 is relevant to this computation.
Addendum. The reader may observe that $G_{1,2,3,4}(z,u)$ are all multiples of $G_1(z,u)$ although this does not appear to simplify the calculation.
A: We seem to have two different results for the variance. To investigate this I wrote some Maple code to compute $[z^n] G(z,u)$ by direct enumeration. As far as what I posted earlier the direct generating function seems to match the algebraic one and the same goes for expectation and variance.
Here is $[z^3] G(z, u)$
$$ \left( 1-p \right) ^{3}+3\,p \left( 1-p \right) ^{2}u+2\,{p}^{2} \left( 1-p \right) u+{p}^{2} \left( 1-p \right) {u}^{2}+{p}^{3}u $$
and this is $[z^4] G(z, u)$
$$  \left( 1-p \right) ^{4}+4\, \left( 1-p \right) ^{3}pu+3\,{p}^{2} \left( 1-p \right) ^{2}u+3\,{p}^{2} \left( 1
-p \right) ^{2}{u}^{2}\\+2\,{p}^{3} \left( 1-p \right) u+2\,{p}^{3} \left( 1-p \right) {u}^{2}+{p}^{4}u.$$
These two generating functions can be verified with pen and paper. To conclude, here is how I computed them. The two functions mom_ex and mom_var compute the moments i.e. expectation and variance. Note that $[z^{11}] G(z,u)$ has only $37$ terms, so this calculation is quite reasonable.


ev :=
proc(n)
    option remember;
    local q, d, pos, pb, pf, r, gf;


    gf := 0;

        for q from 2^n to 2^(n+1)-1 do
            d := convert(q, base, 2);

            if d[1] = 1 then
                r := 1;
                pf := p;
                pb := p;
            else
                r := 0;
                pf := 1-p;
                pb := 1-p;
            fi;

            for pos to n-1 do
                if d[pos] <> d[pos+1] then
                    if d[pos+1] = 1 then r := r+1; fi;
                    pf := 1-pf;
                fi;
                pb := pb*pf;
            od;

            gf := gf+pb*u^r;
        od;

    gf;
end;

mom_ex :=
proc(n)
        expand(subs(u=1, diff(ev(n), u)));
end;

mom_var :=
proc(n)

        expand(subs(u=1, diff(ev(n), u$2)))
        + mom_ex(n)-mom_ex(n)^2;
end;

