How likely is a random integer to be divisible by a group of small primes?

Oscar Olvera Astivia

The first part of this post first appeared in an earlier version “Some probability and number theory” at psychometroscar.com

Oscar Olvera Astivia is a post doctoral research fellow in the The Human Early Learning Partnership (HELP) at the University of British Columbia, Vancouver, BC, Canada.

Oscar completed his PhD in Measurement, Evaluation and Research Methodology (MERM) in 2017 and is in charge of psychometric and statistical analysis of EDI, MDI and other sources of data at HELP.

When I was doing my undergraduate studies, Number Theory was one of my very favourite courses.

I thought I’d go into it but then my interests switched to Probability & Statistics and kind of moved away from the area a little bit.

I still dabble in Number Theory from time to time though, which prompts the subject of this post.

A Computer Science friend of mine and I were talking about cryptography and he casually mentioned that, just as a quick-and-easy check, if you have an unknown number that is not divisible by any of the first 4 prime digits – 2, 3, 5, and 7 – then there’s a “pretty good chance” that said unknown number is prime.

Counting cases

Just to formalize things a little bit, I would like to know what proportion of natural numbers are divisible by 2,3,5 or 7.

Let us start by setting $n\in\mathbb{N}$ . We would like to know what is the probability that $n$ is divisible by 2, 3, 5 or 7.

Since 2 divides half of the set of $\mathbb{N}$ we know right from the beginning that the probability is > 0.5.

The least common multiple of 2, 3, 5 and 7 is 2x3x5x7=210.

We know that 210 has the following list of factors. We build this list by systematically multiplying the numbers 2,3,5,7:

1-element factors: 1, 2, 3, 5, 7

2-element factors: 2×3=6, 2×5=10, 2×7=14, 3×5=15, 3×7=21, 5×7=35

3-element factors: 2x3x5=30, 2x3x7=42, 2x5x7=70, 3x5x7=105

4-element factor: 2x3x5x7=210

Let $S(n)=\{210(n-1)+1, 210(n-1)+2, ..., 210n\}, n=1,2,...$

It is trivial to show $S\subset \mathbb{N}$ and in $S(n)$ there are 210 integers. So by the inclusion-exclusion principle:

$210-(105+70+42+30-35-21-15-14-10-6+7+5+3+2+1)=$ $210-162=48$

which implies there are 48 such integers not divisible by 2, 3, 5 or 7.

Probability of being divisible by 2, 3, 5, or 7

Overall, it is possible to say then that $\frac{48}{210}=\frac{8}{35}\approx 0.22857$ approximates the proportion of natural numbers that are not a multiple of 2, 3, 5 or 7 as $n$ becomes larger and larger.

By the Prime Number Theorem, we know that $\frac{\pi(x)}{x}\sim \frac{1}{ln(x)}$ where $\frac{\pi(x)}{x}$ is the prime-counting function.

So whereas the probability of finding a prime number becomes smaller as $x$ becomes larger, the probability of finding a natural number not divisible by 2, 3, 5 or 7 approaches 0.22857 as $n$ becomes larger and larger.

In general, around 77% of natural numbers are divisible by 2, 3, 5 or 7.

A simulation check

Here’s a quick simulation, written in the R programming language, to verify the reasoning above:

n <- 1:1000000 ##large n, the larger n the closer the true probability
n2 <- n%%2 ## %% gives me the division remainder
n3 <- n%%3
n5 <- n%%5
n7 <- n%%7
dd <- data.frame(n, n2, n3, n5, n7)<
roww <- apply(dd, 1, function(row) all(row !=0 ))
dd1<-dd[roww,]
dim(dd1)[1]/1000000
[1] 0.228571

An extension

(Gary Davis)

Simulation

Oscar’s nice exploration lends itself to an extension, where we ask:

“what’s the likelihood that a randomly chosen number is divisible by at least one of the first k primes (not just the first 4), and how does that probability vary with k?”

We’re going to ignore, for the moment, the niceties of probability measures on countably infinite sets and continue in a vein that Oscar described in the last part of his post: namely we will approximate these probabilities by relative frequencies in simulations.

The idea is to take the first k primes and approximate the probability that a random integer is divisible by at least one of these primes, and then plot that probability as a function of k.

Oscar carried out his simulation in R, and I will do this simulation in Mathematica®.

First we define a function that tells us which, if any, of the first k primes divides a number n:

f[k _, n _] := Or[Table[Mod[n, Prime[i]] == 0, {i, 1, k}]]

For a given k and n this function will return a list of k True or False values, which will be all False exactly when is not divisible by any of the first k primes.

So, given k, we pick out from the first n positive integers those that are divisible by at least one of the first k primes, and then divide the number of these by n to give us our approximate probability, when n is sufficiently large:

prob[k_, n_] := Length[Cases[Range[n], x_ /; MemberQ[f[k, x], True] ]]/n

Then we calculate these probability approximations for small values of k, for a suitably large choice of n, and plot the result:

n = 10^7;
T = Table[{k, prob[k, n]}, {k, 1, 5}]
ListPlot[T, Joined -> True, Mesh -> All]

with the result:

So by the time we’ve got up to the first 20 primes: 2, 3, 5, 7, 11, 13, 17, 19, 23, 29, 31, 37, 41, 43, 47, 53, 59, 61, 67, 71, the likelihood that a random integer is divisible by at least one of these 20 primes is approximately 87.2%

Hmm!! What is the rate at which these (approximate) probabilities are approaching 1?

Inclusion-exclusion

Oscar observes that the probability of a positive integer being divisible by 2 is 1/2.

What’s the probability of an integer being divisible by 2 or3?

In counting those n divisible 2 and those divisible by 3 we’ve counted those divisible by 6 twice, so the probability that n is divisible by 2 or by 3 is $1/2+1/3-1/2\times 1/3 = 2/3 \approx 0.667$ .

This is an instance of the “measure” of a union of two sets set (the probabilities of two infinite sets of positive integers) being the sum of their measures minus the measure of their intersection.

Similarly, the probability that n is divisible by 2 or 3 or 5 is $1/2+1/3+1/5 -(1/2\times 1/3+1/2\times 1/5+1/3\times 1/5)+1/2\times 1/3\times 1/5 = 11/15 \approx 0.733$ .

When we extend this to the 4 primes 2, 3, 5, 7 we see that the probability that a positive integer n is divisible by at least one of these 4 primes is:

$1/2 + 1/3 + 1/5 + 1/7 - (1/2 \times 1/3 + 1/2\times 1/5 + 1/2\times 1/7 + 1/3\times 1/5 + 1/3\times 1/7 + 1/5\times 1/7) + (1/2\times 1/3\times 1/5 + 1/2\times 1/3\times 1/7 + 1/2\times 1/5\times 1/7 + 1/3\times 1/5\times 1/7) - (1/2\times 1/3\times 1/5\times 1/7) = 27/35 \approx 0.771$

To 6 decimal places the result is 0.771429 in excellent agreement with our simulation.

This idea can be extended to the first k primes.

Let’s denote the probability that a random positive integer n is divisible by at least one of the first k primes as prob(k).

Then by using the formula for the probability of a union of events we can calculate prob(k) recursively as follows:

prob[1] = 1/2;
prob[k_] := prob[k] = prob[k – 1] + 1/Prime[k] – prob[k – 1]/Prime[k]

Using this recursive formula we can calculate the exact probabilities and their approximations for k from 1 through 20, for example:

A = Table[{k, prob[k], N[prob[k]]}, {k, 1, 20}]
TableForm[A, TableHeadings -> {{}, {“k”, “prob(k)”, “Approximation”}}

We can imagine asking Oscar’s computer science friend if he would rather have a 95% likelihood that a random number was divisible by the first k primes: would he like to know k?

k = 1;
While[prob[k] < 0.95, k++]
k

with the result: 7398.

So to be 95% confident that a number n is divisible by the first k primes we need to check up to the $7398^{th}$ prime, which is 75037.