# (Lack of) Statistical Thinking in the NYT

The [New York Times yesterday](https://www.nytimes.com/2023/11/02/podcasts/run-up-trump-biden-election.html) (Nov 2, 2023) had a story and podcast about the "**only \[US county\] that has voted for the winner of the \[US\] presidential race every year since 1980**."

They write:

![](https://cdn.hashnode.com/res/hashnode/image/upload/v1698956571667/a43ad636-67aa-4fa8-b08a-81b7691c66f7.png align="center")

**"Prediction," you say?** That will awaken any quant!

Also, **it's about Clallam County, Washington — here in my home state.** To paraphrase a famous political quote, I can see it from my house (at least its mountains). It's a great place, and one of my favorite bookstores is there in Port Angeles, WA. You should check out the superb photography in the NYT story!

OK, it's a human interest story — complete with the required diner! — and shouldn't be taken seriously. Still, it is in the *Times* and it makes a somewhat explicit statistical claim. **So it's fair game to examine whether voters in Clallam County, Washington (in a great diner) are "resembling a prediction."**

In this post, I examine the following question: **given 11 elections, how many of the 3000 counties would be on the "winning" side in all of them, simply due to random chance?**

Is the NYT reporting something interesting? Or are they making a story from random data? Here are some simple ways I answer that.

---

### Approach 1: Mental Arithmetic

For simplicity, let's assume that every election is exactly 50/50. If we make a single guess, we have p=0.5 odds of getting it correct, or 1 chance out of 2 outcomes. With a second election, there are 4 possible outcomes (2^2), and the random odds of getting both correct are 1 out of 4.

In 11 elections there are 2^11 outcomes. If you know the powers of 2, as many programmers do, that's 2048 outcomes.

Or you might approximate it as follows: `2^8 = 256` and `2^3 = 8`, so `2^11 is 256 * 8`, or about `250 * 8`, which is 2000. Or `250 * 10` less 20%, which is `2500 - 500` = 2000.

From that, **we might expect about 1/2000 counties to get 11 elections correct, randomly**. An observation of 1/3000 is roughly what we'd expect randomly.

---

### Approach 2: Binomial Confidence Intervals

We can formalize that with a binomial test. Suppose we have 11 independent elections, with 0.5 odds in each one, 3000 counties, and observe 1 success case.

Here's the R code:

```r
# number of US presidential elections 1980 - 2020 
n.elec <- (2020-1980) / 4 + 1
# supposed odds of voting with the "winning" side in one election
p.elec <- 0.50
# number of counties in the US
n.cty  <- 3000   # roughly

# odds of observing 1 success in 3000 counties, assuming 11 elections at 0.5 odds
binom.test(x=1, n=n.cty, p = p.elec ^ n.elec)
```

The result is:

![](https://cdn.hashnode.com/res/hashnode/image/upload/v1698957538972/aa487560-b120-46a6-8c3a-4baf1bd392d0.png align="center")

To expand on that, let's *multiply the confidence intervals by the sample size* to get **the expected number of successes** rather than the *proportion*:

```r
# confidence intervals around that observation 
binom.test(x=1, n=n.cty, p = p.elec ^ n.elec)$conf.int * n.cty  # (0-6 expected "successes")
```

That is:

![](https://cdn.hashnode.com/res/hashnode/image/upload/v1698957633185/9f9e2c35-1381-4bb1-95fc-defad776077a.png align="center")

Under our assumptions, **we expect 0 - 6** counties to get it right. **So observing 1 county is not interesting statistically.** (It's still a well-photographed story!)

R has other binomial estimation methods that are robust to extremely small proportions, in the `binom` package. Here's the code and result:

```r
# estimates using multiple methods
library(binom)              # install if needed
binom.confint(x=1, n=n.cty)
# see the confidence intervals (0-7 successes expected)
binom.confint(x=1, n=n.cty)[ , c("lower", "upper")] * n.cty
```

![](https://cdn.hashnode.com/res/hashnode/image/upload/v1698957790516/23785b37-383a-4f36-8250-98e9ae5c8e0d.png align="center")

**Using every binomial estimation method, we expect somewhere between 0 and 3-7 counties to get it right** (again, under the assumptions we made).

---

### Approach 3: Binomial Distribution

Another binomial approach is to examine the density — i.e., the likelihood — of particular outcomes. In R, `dbinom()` is the density function for the exact binomial test.

**What's the most likely observation for the number of counties?** Let's plot the likelihood of observing 0:10 successes (i.e., counties) that correctly call all 11 elections in a row under our assumptions. First the R code:

```r
# plot of the density function shows the most likely number of successes
plot(x=0:10, y=dbinom(0:10, n.cty, p = p.elec ^ n.elec), 
     main="Expected using coin flip odds", xlab="Number of counties, assuming independent random events")
```

Here's the plot:

![](https://cdn.hashnode.com/res/hashnode/image/upload/v1698961611519/ac8534bb-c0be-4c75-98b5-a6fb965175ce.png align="center")

If the county results are random and independent, **the single most likely result is that exactly 1 county would align with 11 elections**. Maybe 0 or 2, or 3, or 4. Almost certainly fewer than 8.

**Although I love Clallam County, WA, and enjoy reading about the people there, the result highlighted in the New York Times appears to be random.**

---

### But It's Not Random

If you follow elections closely, you'll be aware of the limitations of such simple analyses. Each variable in my simple binomial assumptions is imprecise:

1. **The odds of aligning with the outcome are not 50/50**. When there is a landslide election (like Reagan in 1984) more counties will have a higher likelihood of aligning with the vote.
    
2. **Votes are not equally distributed across counties**. Most counties have small populations, and it's possible to win in a relative landslide while winning few counties ([In 2012, Obama](https://www.nbcnews.com/id/wbna50073771) won 51% of the vote and 62% of the electoral college, but only 22% of counties.)  
    *Side question*: are there articles about counties that are *wrong* in every election? Statistically, that would be just as interesting!
    
3. Counties are not predictors of, or otherwise related to, national elections. Instead, *for national elections*, **counties are just a way to aggregate and report votes**. The nests are `Electoral College <== States <== Voters`, and separately `Counties <== Voters`. (Although counties are nested inside states geographically, they are not meaningful units in an electoral sense.)
    
4. **Results from one election to the next — the 11 wins in a row — are not independent probabilities**. Many counties never vary, always voting the same way in presidential elections.
    
5. **Counties are correlated and not independent of one another**, so there are not 3000 independent actors.
    

One might model each of those things under various assumptions (hello, [FiveThirtyEight](https://projects.fivethirtyeight.com/polls/)). Overall I would *expect* that such modeling to have modest effects — highly significant in terms of *predicting elections* but unlikely to have a vastly different answer to the question we examined.

However, the primary point of this post remains: **simple stats and mental arithmetic suggest that although the New York Times story might be fun, it is not serious**. They should know better.

But did I mention that the photos are excellent? I recommend visiting Clallam County!

---

### All the Code

Here's the complete R code — plus **a brief, bonus snippet** for the Electoral College vote shares. (For more on binomial models, see Chapter 4 in the [R book](https://r-marketing.r-forge.r-project.org) or [Python book](https://link.springer.com/book/10.1007/978-3-030-49720-0).)

```r
# probability estimates
# re New York Times, https://www.nytimes.com/2023/11/02/podcasts/run-up-trump-biden-election.html

# number of US presidential elections 1980 - 2020 
n.elec <- (2020-1980) / 4 + 1
# supposed odds of voting with the "winning" side in one election
p.elec <- 0.50
# number of counties in the US
n.cty  <- 3000   # roughly

# odds of observing 1 success in 3000 counties, assuming 11 elections at 0.5 odds
binom.test(x=1, n=n.cty, p = p.elec ^ n.elec)
# confidence intervals around that observation 
binom.test(x=1, n=n.cty, p = p.elec ^ n.elec)$conf.int * n.cty  # (0-6 expected "successes")

# estimates using multiple methods
library(binom)              # install if needed
binom.confint(x=1, n=n.cty)
# see the confidence intervals (0-7 successes expected)
binom.confint(x=1, n=n.cty)[ , c("lower", "upper")] * n.cty

# plot of the density function shows the most likely number of successes
plot(x=0:10, y=dbinom(0:10, n.cty, p = p.elec ^ n.elec), 
     main="Expected using coin flip odds", xlab="Number of counties, assuming independent random events")

# history of electoral college proportions
margin.elec <- c(0.909, 0.976, 0.792, 0.688, 0.704, 0.504, 0.532, 0.678, 0.617, 0.565, 0.569)
(mean.elec   <- exp(mean(log(margin.elec)))) # geometric mean

# expected number of "correct calls", using electoral college probabilities
plot(x=0:60, y=dbinom(0:60, n.cty, p = mean.elec ^ n.elec), 
     main="Expected using electoral college", xlab="Number of 'wins' out of 3000")
```

![](https://cdn.hashnode.com/res/hashnode/image/upload/v1748018891345/c59ac669-83a5-44db-a75b-f0e772b875f6.png align="center")