How can I determine which of two sequences of coin flips is real and which is fake? This is an interesting problem I came across. I'm attempting to write a Python program to get a solution to it; however, I'm not sure how to proceed. So far, I know that I would expect the counts of heads to follow a binomial, and length of runs (of tails, heads, or both) to follow a geometric.

Below are two sequences of 300 “coin flips” (H for heads, T for
tails). One of these is a true sequence of 300 independent flips of a
fair coin. The other was generated by a person typing out H’s and T’s
and trying to seem random. Which sequence is truly composed of coin
flips?
Sequence 1:
TTHHTHTTHTTTHTTTHTTTHTTHTHHTHHTHTHHTTTHHTHTHTTHTHH
TTHTHHTHTTTHHTTHHTTHHHTHHTHTTHTHTTHHTHHHTTHTHTTTHH
TTHTHTHTHTHTTHTHTHHHTTHTHTHHTHHHTHTHTTHTTHHTHTHTHT
THHTTHTHTTHHHTHTHTHTTHTTHHTTHTHHTHHHTTHHTHTTHTHTHT
HTHTHTHHHTHTHTHTHHTHHTHTHTTHTTTHHTHTTTHTHHTHHHHTTT
HHTHTHTHTHHHTTHHTHTTTHTHHTHTHTHHTHTTHTTHTHHTHTHTTT
Sequence 2:
HTHHHTHTTHHTTTTTTTTHHHTTTHHTTTTHHTTHHHTTHTHTTTTTTH
THTTTTHHHHTHTHTTHTTTHTTHTTTTHTHHTHHHHTTTTTHHHHTHHH
TTTTHTHTTHHHHTHHHHHHHHTTHHTHHTHHHHHHHTTHTHTTTHHTTT
THTHHTTHTTHTHTHTTHHHHHTTHTTTHTHTHHTTTTHTTTTTHHTHTH
HHHTTTTHTHHHTHHTHTHTHTHHHTHTTHHHTHHHHHHTHHHTHTTTHH
HTTTHHTHTTHHTHHHTHTTHTTHTTTHHTHTHTTTTHTHTHTTHTHTHT
Both sequences have 148 heads, two
less than the expected number for a 0.5 probability of heads.


Looking at the HH, HT, TH, TT frequencies is probably the most straightforward way to approach the two series presented, given people's tendency to apply HT and TH more frequently when trying to appear random. More generally, however, that approach will fail to detect non-randomness even in sequences with obvious patterns. For instance, a repeating sequence of HHHTHTTT will produce balanced counts of pairs as well as triples (HHH, HHT, etc.).
Here's an idea for a more general solution that was initially inspired by answers discussing random walks. Start with the observation that for a random sequence, the number of heads in any given subsequence is binomially distributed. We can count the number of heads, $n_{i,j}$, between flip number $i$ and flip number $j>i$ for all values of $i,j$. Then compare these to the expected number of counts if all the $n_{i,j}$ were independent. Although they are obviously not independent, the comparison gives rise to a useful statistic: the maximum absolute difference between the observed counts and the expected counts.
Applying this approach to sequence 1 gives us 98.26 as our test statistic: there are 257 subsequences of length 44. If they were all composed of independent Bernoulli trials, the expected number of the 257 that contained exactly 22 heads is ~30.74, whereas sequence 1 contains 129 subsequences with exactly 22 heads (very underdispersed). 129 - 30.74 = 98.26, which is the maximum of these differences for sequence 1.
Performing the same calculations on sequence 2 gives a test statistic of 48.30: there are 197 subsequences of length 104. The expected number containing exactly 54 heads would be ~13.70. Sequence 2 contains 62, so the test statistic is 62 - 13.70 = 48.30.
The test statistics can be compared to those from a large number of random sequences of the same size. In this case, no samples are greater than the test statistic from sequence 1, and about 14% of samples are greater than the statistic from sequence 2.
Here it is all together with R code:
library(parallel)
set.seed(16)

seq1 <- utf8ToInt("TTHHTHTTHTTTHTTTHTTTHTTHTHHTHHTHTHHTTTHHTHTHTTHTHHTTHTHHTHTTTHHTTHHTTHHHTHHTHTTHTHTTHHTHHHTTHTHTTTHHTTHTHTHTHTHTTHTHTHHHTTHTHTHHTHHHTHTHTTHTTHHTHTHTHTTHHTTHTHTTHHHTHTHTHTTHTTHHTTHTHHTHHHTTHHTHTTHTHTHTHTHTHTHHHTHTHTHTHHTHHTHTHTTHTTTHHTHTTTHTHHTHHHHTTTHHTHTHTHTHHHTTHHTHTTTHTHHTHTHTHHTHTTHTTHTHHTHTHTTT") == utf8ToInt("H")
seq2 <- utf8ToInt("HTHHHTHTTHHTTTTTTTTHHHTTTHHTTTTHHTTHHHTTHTHTTTTTTHTHTTTTHHHHTHTHTTHTTTHTTHTTTTHTHHTHHHHTTTTTHHHHTHHHTTTTHTHTTHHHHTHHHHHHHHTTHHTHHTHHHHHHHTTHTHTTTHHTTTTHTHHTTHTTHTHTHTTHHHHHTTHTTTHTHTHHTTTTHTTTTTHHTHTHHHHTTTTHTHHHTHHTHTHTHTHHHTHTTHHHTHHHHHHTHHHTHTTTHHHTTTHHTHTTHHTHHHTHTTHTTHTTTHHTHTHTTTTHTHTHTTHTHTHT") == utf8ToInt("H")

# function to calculate test statistic S 
fS <- function(m, p = 0.5) {
  if (!is.matrix(m)) m <- matrix(m, ncol = 1)
  n <- nrow(m)
  steps <- rep.int(1:(n - 1L), 2:n)
  nEx <- dbinom(sequence(2:n, from = 0L), steps, p)*(n - steps)
  idx <- sequence((n - 1L):1)
  idx <- idx*(idx + 1L)/2L
  S <- numeric(ncol(m))
  for (i in 1:ncol(m)) S[i] <- max(abs(nEx - tabulate(idx + dist(cumsum(m[,i])), length(nEx))))
  S
}

(S1 <- fS(seq1))
#> [1] 98.26173
(S2 <- fS(seq2))
#> [1] 48.30014
# calculate S from 1e6 random sequences of length n (probably overkill)
n <- length(seq1)
cl <- makeCluster(detectCores() - 1L)
clusterExport(cl, list("fS", "n"))
system.time(simS <- unlist(parLapply(cl, 1:100, function(i) fS(matrix(sample(0:1, 1e4L*n, TRUE), n)))))
#>    user  system elapsed 
#>    0.00    0.03  339.42
stopCluster(cl)
nsim <- 1e6L
# calculate approximate p-values
sum(simS > S1)/nsim
#> [1] 0
sum(simS > S2)/nsim
#> [1] 0.139855
max(simS)
#> [1] 91.01078

