One of the libraries in R has a runs.test procedure, which you can explore. My purpose here here is to give
an idea how looking at runs can help you decide whether
your observations are randomly sampled from the same population.
To begin we look specifically at sequences of Bernoulli
trials, as mentioned in your Question. (Randomness tests
for other distributions can be made 'Bernoulli' by looking
at numbers of observations above or below the sample mean or median.)
Simple example. Suppose we have a sample of size $N=10,$ purported to be from a population of Bernoulli trials with Success probability $p = 1/2,$ and that five of the observations are Successes (1s) and five are Failures (0). Then there are ${10 \choose 5} = 252$ possible arrangements
of the 0s and 1s.
choose(10,5)
[1] 252
Possible numbers of runs
are between $2$ (all five 0s first or all five 1s first) and $10$ (alternating 0s and 1s). One can show that the average number is 6 (halfway between).
There are only two ways out of 252 to get $2$ runs
and only two ways to get $10$ runs.
So, under the null
hypothesis that 0s and 1s occur at random, the probability of seeing one of these extreme numbers of runs is $4/252 \approx 0.016$, and we would reject the null hypothesis.
Perhaps we are seeing five observations from a Bernoulli process with $p = .1$ followed by five from a different Bernoulli process with $p=.9.$ Or output from a (non-independent) Markov process that alternates easily between states 0 and 1, but seldom stays in the same state on successive steps. (As here..)
Example with 100 Bernoulli observations. Suppose we have $N=100$ observations, 0 or 1, from a process purported to be a random sample 100, all randomly chosen from the same Bernoulli process.
Now we can use some help counting the runs.
The R procedure rle (for Run Length Encoding) shows
the number of runs, the value during each run, and the length of each run. For example, one sample of $N=100$ might have $m=29$ 1s, and $r=43$ runs.
rle(x)
Run Length Encoding
lengths: int [1:43] 1 1 1 2 2 1 2 1 4 1 ...
values : int [1:43] 0 1 0 1 0 1 0 1 0 1 ...
length(rle(x)$val)
[1] 43
table(x)
x
0 1
71 29
sum(x==1)
[1] 29
Given the number $m$ of Successes and the number of
Failures $n = N-m,$ there are formulas for the average number $\mu = E(R)$ of runs in a random sample and the variance $\sigma^2 = Var(R):$
$$ \mu = \frac{2mn}{N} + 1,\;\; \sigma^2 =\frac{(\mu-1)(\mu-2)}{N-1}$$
Moreover, for a sample as large as $N=100,$ the distribution of $R$ is nearly normal (especially in the tails, where it matters). So we can reject
the null hypothesis that the data are a random sample from a single population if $Z = \frac{R-\mu}{\sigma},$
has $|Z| \ge 1.96.$ [Perhaps see Wikipedia.]
The following simulation illustrates that such a test at the 5% level actually does reject about 5% of the time
for truly random data.
set.seed(2020)
B = 10^4; z = numeric(B)
for(i in 1:B){
x = rbinom(100, 1, .3)
m = sum(x==1)
n = sum(x==0)
r = length(rle(x)$val)
a = 2*m*n; N = m+n;
mu = a/N+1; vr = (mu-1)*(mu-2)/(N-1)
z[i] = (r-mu)/sqrt(vr) }
mean(abs(z) >=1.96)
[1] 0.049 # aprx P-reject = 0.05
mean(z); sd(z)
[1] 0.0003448186 # aprx E(Z) = 0
[1] 0.9963706 # aprx Var(Z) = 0
The following histogram shows the simulated distribution
of the approximate test statistic. [A histogram with half as many bars (without parameter br=30) looks a lot closer to normal near $0,$ but this one
gives a more honest view.]
hdr="Simulated Z with Standard Normal PDF"
hist(z, prob=T, br=30, col="skyblue2", main=hdr)
curve(dnorm(x), add=T, col="red", lwd=2)
abline(v = c(-1.96,1.96), lty="dotted")
Note: If the line of the program for generating
the Bernoulli sample is changed as shown below,
suggesting a non-random mixture of two Bernoulli processes
with $p = 0.05$ and $p = 0.55$ (giving 30% successes 'on average'), then the rejection rate
increases to about 80%.
...
x = c(rbinom(50, 1, .05),rbinom(50,1,.55))
...
mean(abs(z) >=1.96)
[1] 0.7954
By contrast, a test that the proportion of Successes
is $p = 0.3$ is not rejected for one such nonrandom sample.
set.seed(1234)
x = c(rbinom(50, 1, .05),rbinom(50,1,.55))
table(x)
x
0 1
63 37
prop.test(37,100, p=.3)
1-sample proportions test
with continuity correction
data: 37 out of 100, null probability 0.3
X-squared = 2.0119, df = 1, p-value = 0.1561
alternative hypothesis: true p is not equal to 0.3
95 percent confidence interval:
0.2772627 0.4728537
sample estimates:
p
0.37
