John D. Cook

Suppose you do N trials of something that can succeed or fail. After your experiment you want to present a point estimate and a confidence interval. Or if you’re a Bayesian, you want to present a posterior mean and a credible interval. The numerical results hardly differ, though the two interpretations differ. If you got […]The post Confidence interval widths first appeared on John D. Cook.

The previous post looked at a cellular automaton introduced by Edward Fredkin. It has only two states: a cell is on or off. At each step, each cell is set to the sum of the states of the neighboring cells mod 2. So a cell is on if it had an odd number neighbors turned […]The post Multicolor reproducing cellular automata first appeared on John D. Cook.

Edward Fredkin is best known these days for the Fredkin gate, a universal reversible circuit. I recently found out that Fredkin was one of the pioneers in cellular automata. His student Edwin Banks describes a cellular automaton introduced by Fredkin [1]. Edward Fredkin of MIT has described the following interesting cellular space. If the states […]The post Self-reproducing cellular automata first appeared on John D. Cook.

The other day I woke up with a song in my head I hadn’t heard in a long time, the hymn Beneath the Cross of Jesus. The name of the tune is St. Christopher. When I thought about the tune, I realized it has some fairly sophisticated harmony. My memory of the hymns I grew […]The post Humming St. Christopher first appeared on John D. Cook.

The previous post looked at repeatedly applying the function s(n) which is the sum of the divisors of n less than n. It is an open question whether the sequence s( s( s( … s(n) … ) ) ) always converges or enters a loop. In fact, it’s an open question of whether the sequence […]The post Aliquot ratio distribution first appeared on John D. Cook.

Let s(n) be the sum of the proper divisors of n, the divisors less than n itself. A number n is called excessive, deficient, or perfect depending on whether s(n) – n is positive, negative, or 0 respectively, as I wrote about a few days ago. The number s(n) is called the aliquot sum of […]The post The iterated aliquot problem first appeared on John D. Cook.

The term “cologarithm” was once commonly used but now has faded from memory. Here’s a plot of the frequency of the terms cololgarithm and colog from Google’s Ngram Viewer. The cologarithm base b is the logarithm base 1/b, or equivalently, the negative of the logarithm base b. cologb x = log1/b x = -logb x […]The post Cologarithms and Entropy first appeared on John D. Cook.

When students learn about decimals, they’re told that every fraction either has a terminating decimal expansion or a repeating decimal expansion. The previous post gives a constructive proof of the converse: given a repeating fraction (in any base) it shows how to find what rational number it corresponds to. Maybe you learned this so long […]The post Corollary of a well-known fact first appeared on John D. Cook.

My previous post ended with a discussion of repeating binary decimals such as 0.00110011…two = 1/5. For this post I’ll explain how calculations like that are done, how to convert a repeating decimal in any base to a fraction. First of all, we only need to consider repeating decimals of the form 0.1b, 0.01b, 0.001b, […]The post Repeating decimals in any base first appeared on John D. Cook.

Yesterday I wrote about the bits in powers of 3. That post had a low-resolution image, which has its advantages, but here’s a higher resolution image that also has its advantages. The image looks chaotic. I say this in the colloquial sense, not in the technical sense as in period three implies chaos. I just […]The post Chaotic image out of regular slices first appeared on John D. Cook.

The previous post showed that the bits of prime powers look kinda chaotic. When you plot them, they form a triangular shape because the size of the numbers is growing. The numbers are growing geometrically, so their binary representations are growing linearly. Here’s the plot for powers of 5: You can crop the triangle so […]The post Evolution of random number generators first appeared on John D. Cook.

I was thumbing through A New Kind of Science [1] and one of the examples that jumped out at me was looking at the bits in the binary representation of the powers of 3. I wanted to reproduce the image myself and here’s the result. Here a white square represents a 1 and a blue […]The post Powers of 3 in binary first appeared on John D. Cook.

I was curious about the orders of the sporadic groups, specifically their prime factors, so I made a bar graph to show how often each factor appears. To back up a little bit, simple groups are to groups roughly what prime numbers are to integers. Finite simple groups have been fully classified, and they fall […]The post Factors of orders of sporadic groups first appeared on John D. Cook.

I learned recently that weekly excess deaths in the US have dipped into negative territory [1], and wondered whether we should start speaking of deficient deaths by analogy with excessive and deficient numbers. The ancient Greeks defined a number n to be excessive, deficient, or perfect according to whether the sum of the number’s proper […]The post Excessive, deficient, and perfect numbers first appeared on John D. Cook.

The grep utility searches text files for regular expressions, but it can search for ordinary strings since these strings are a special case of regular expressions. However, if your regular expressions are in fact simply text strings, fgrep may be much faster than grep. Or so I’ve heard. I did some benchmarks to see. Strictly […]The post Is fast grep faster? first appeared on John D. Cook.

This post will look at the “tower of powers” and ask what it means, when it converges, and how to compute it. Along the way I’ll tie in two recent posts, including one that should come as a surprise. First of all, the expression is right-associative. That is, it is the limit of x^(x^(x^…)) and […]The post Tower of powers and convergence first appeared on John D. Cook.

I was studying a statistics paper the other day in which the author said to solve t log( 1 + n/t ) = k for t as part of an algorithm. Assume 0 < k < n. Is this well posed? First of all, can this equation be solved for t? Second, if there is […]The post Lambert W strikes again first appeared on John D. Cook.

Suppose you need to compare two files with very long lines. Maybe each file is one line. The diff utility will show you which lines differ. But if the files are each one line, this tells you almost nothing. It confirms that there is a difference, but it doesn’t show you where the difference is. […]The post Comparing files with very long lines first appeared on John D. Cook.

In October 2020 there were 619 new codes added to the list of ICD-10 diagnosis codes. Some of these make sense, such as four codes added for COVID: U07.1 COVID-19 Z11.52 Encounter for screening for COVID-19 Z20.822 Contact with and (suspected) exposure to COVID-19 Z86.16 Personal history of COVID-19 There are also 24 new codes […]The post New ICD-10 diagnosis code weirdness first appeared on John D. Cook.

Allen William Johnson [1] discovered the following magic square whose entries are all squares. The following Python code verifies that this is a magic square. import numpy as np M = np.array( [[ 30**2, 246**2, 172**2, 45**2], [ 93**2, 116**2, 66**2, 258**2], [126**2, 138**2, 237**2, 44**2], [260**2, 3**2, 54**2, 150**2] ]) def verify(M): m, n […]The post Magic square of squares first appeared on John D. Cook.

There is a lot of talk about Mars right now, and understandably so. The flight of Ingenuity today was awesome. As Daniel Oberhaus pointed out on Twitter, … the atmosphere on the surface of Mars is so thin that it’s the equivalent of flying at ~100k feet on Earth. No rotorcraft, piloted or uncrewed, has […]The post Martian gravity first appeared on John D. Cook.

A cryptographic hash is also known as a one-way function because given an input x, one can quickly compute the hash h(x), but it is extremely time-consuming to try to recover x if you only know h(x). Even if the hashing algorithm is considered “broken,” it may take an enormous effort to break it. Google […]The post Hashing phone numbers first appeared on John D. Cook.

I ran across the word duodecimal recently, the Latin-derived term for base 12, and realized it’s been years, maybe decades, since I’d heard that term. I’ve used hexadecimal explicitly in 30 different blog posts, and implicitly in more, but this is the first time I’ve used duodecimal. Hexadecimal comes up constantly in application. I suppose […]The post Duodecimal vs Hexadecimal first appeared on John D. Cook.

Mathematics notation changes slowly over time, generally for the better. I can’t think of an instance that I think was a step backward. Gauss introduced the notation [x] for the greatest integer less than or equal to x in 1808. The notation was standard until relatively recently, though some authors used the same notation to […]The post Floor, ceiling, bracket first appeared on John D. Cook.

Start with a unit circle and draw the largest square you can inside the circle and the smallest square you can outside the circle. In geometry lingo these are the inscribed and circumscribed squares. The circle fits inside the square better than the square fits inside the circle. That is, the ratio of the area […]The post Box in ball in box in high dimension first appeared on John D. Cook.

A few weeks ago I wrote several blog posts about very simple approximations that are surprisingly accurate. These approximations are valid over a limited range, but with range reduction they can be used over the full range of the functions. In this post I want to look again at and Padé approximation It turns out […]The post More on why simple approximations work first appeared on John D. Cook.

Is 0 a stable fixed point of f(x) = 4x (1-x)? If you set x = 0 and compute f(x) you will get exactly 0. Apply f a thousand times and you’ll never get anything but zero. But this does not mean 0 is a stable attractor, and in fact it is not stable. It’s […]The post Mathematical stability vs numerical stability first appeared on John D. Cook.

Some programming languages, such as Perl, have an infix operator <=> that returns a three-state comparison. The expression a <=> b evaluates to -1, 0, or 1 depending on whether a < b, a = b, or a > b. You could think of <=> as a concatenation of <, =, and >. The <=> operator […]The post Spaceship operator in Python first appeared on John D. Cook.

The previous post Does chaos imply period 3? ended with looking at a couple cubic polynomials whose roots have period 3 under the mapping f(x) = 4x(1-x). These are 64 x³ – 112 x² + 56 x – 7 and 64 x³ – 96 x² + 36 x – 3. I ended the post by […]The post Real radical roots first appeared on John D. Cook.

Sharkovskii’s theorem says that if a continuous map f from an interval I to itself has a point with period 3, then it has a point with period 5. And if it has a point with period 5, then it has points with order 7, etc. The theorem has a chain of implications that all […]The post Does chaos imply period 3? first appeared on John D. Cook.

The previous post explained what is meant by period three implies chaos. This post is a follow-on that looks at Sarkovsky’s theorem, which is mostly a generalization of that theorem, but not entirely [1]. First of all, Mr. Sarkovsky is variously known Sharkovsky, Sharkovskii, etc. As with many Slavic names, his name can be anglicized […]The post Sarkovsky’s theorem first appeared on John D. Cook.

One of the most famous theorems in chaos theory, maybe the most famous, is that “period three implies chaos.” This compact statement comes from the title of a paper [1] by the same name. But what does it mean? This post will look at what the statement means, and the next post will look at […]The post Period three implies chaos first appeared on John D. Cook.

A while back I presented a very simple algorithm for computing natural logs: log(x) ≈ (2x – 2)(x + 1) for x between exp(-0.5) and exp(0.5). It’s accurate enough for quick mental estimates. I recently found an approximation by Ronald Doerfler that is a little more complicated but much more accurate: log(x) ≈ 6(x – […]The post Better approximation for ln, still doable by hand first appeared on John D. Cook.

It occurred to me recently that a problem I solved numerically years ago could be solved analytically, the problem of determining beta distribution parameters so that the distribution has a specified mean and variance. The calculation turns out to be fairly simple. Maybe someone has done it before. Problem statement The beta distribution has two […]The post Beta distribution with given mean and variance first appeared on John D. Cook.

The following equation is almost true. And by almost true, I mean correct to well over 200 decimal places. This sum comes from [1]. Here I will show why the two sides are very nearly equal and why they’re not exactly equal. Let’s explore the numerator of the sum with a little code. >>> from […]The post Close but no cigar first appeared on John D. Cook.

The previous post made use of both the arithmetic and geometric means. It also showed how both of these means correspond to different points along a continuum of means. This post combines those ideas. Let a and b be two positive numbers. Then the arithmetic and geometric means are defined by A(a, b) = (a […]The post Arithmetic-geometric mean first appeared on John D. Cook.

The previous post looked at a simple method of finding square roots that amounts to a special case of Newton’s method, though it is much older than Newton’s method. We can extend Newton’s method to find cube roots and nth roots in general. And when we do, we begin to see a connection to r-means. […]The post Higher roots and r-means first appeared on John D. Cook.

Here’s a simple way to estimate the square root of a number x. Take a guess g at the root and compute the average of g and x/g. If you want to compute square roots mentally or with pencil and paper, how accurate can you get with this method? Could you, for example, get within […]The post Calculating square roots first appeared on John D. Cook.

A couple years ago I wrote a blog post on Kepler’s equation M + e sin E = E. Given mean anomaly M and eccentricity e, you want to solve for eccentric anomaly E. There is a simple way to solve this equation. Define f(E) = M + e sin E and take an initial guess at the solution and […]The post Efficiently solving Kepler’s equation first appeared on John D. Cook.

The exponential sum page on this site draws lines between the consecutive partial sums of where m is the month, d is the day, and y is the last two digits of the year. The sum for today is unusually round: By contrast, the sum from yesterday is nowhere near round: Out of curiosity, I […]The post Unusually round exponential sum first appeared on John D. Cook.

Ashley Kanter left a comment on Tuesday’s post Within one percent with an approximation I’d never seen. One that I find handy is the hypotenuse of a right-triangle with other sides a and b (where a<b) can be approximated to within 1% by 5(a+b)/7 when 1.04 ≤ b/a ≤1.50. That sounds crazy, but it’s right. […]The post Hypotenuse approximation first appeared on John D. Cook.

Richard Feynman said nearly everything is really interesting if you go into it deeply enough. In that spirit I’m going to dig into the units on Coulomb’s constant. This turns out to be an interesting rabbit trail. Coulomb’s law says that the force between two charged particles is proportional to the product of their charges […]The post Coulomb’s constant first appeared on John D. Cook.

This post looks at some common approximations and determines the range over which they have an error of less than 1 percent. So everywhere in this post “≈” means “with relative error less than 1%.” Whether 1% relative error is good enough completely depends on context. Constants The familiar approximations for π and e are […]The post Within one percent first appeared on John D. Cook.

The regular expression search utility grep has a recursive switch -R, but it may not work like you’d expect. Suppose want to find the names of all .org files in your current directory and below that contain the text “cheese.” You have four files, two in the working directory and two below, that all contain […]The post Recursive grep first appeared on John D. Cook.

Dogleg severity (DLS) is essentially what oilfield engineers call curvature. Oil wells are not simply vertical holes in the ground. The boreholes curve around underground, even if the intent is to drill a perfectly straight hole. With horizontal drilling the curvature can be substantial. Dogleg severity is calculated by measuring inclination and azimuth every 100 […]The post Calculating dogleg severity first appeared on John D. Cook.

Möbius transformations are functions of the form where ad – bc ≠ 0. A Möbius transformation is uniquely determined by its values at three points. Last year I wrote a post that mentioned how to determine the coefficients of a Möbius transformation. There I said The unique bilinear transform sending z1, z2, and z3 to […]The post Solving for Möbius transformation coefficients first appeared on John D. Cook.

The following table gives the best rational approximations to π, e, and φ (golden ratio) for a given accuracy goal. Here “best” means the fraction with the smallest denominator that meets the accuracy requirement. I found these fractions using Mathematica’s Convergents function. For any irrational number, the “convergents” of its continued fraction representation give a […]The post Smallest denominator for given accuracy first appeared on John D. Cook.

I find simple approximations more interesting than highly accurate approximations. Highly accurate approximations can be interesting too, in a different way. Somebody has to write the code to compute special functions to many decimal places, and sometimes that person has been me. But somewhat ironically, these complicated approximations are better known than simple approximations. One […]The post Simple approximation for Gamma function first appeared on John D. Cook.

A few days ago I wrote about how to estimate 10x. This is an analogous post for exp(x) = ex. We will assume -0.5 ≤ x ≤ 0.5. You can bootstrap your way from there to other values of x. For example, exp(1.3) = exp(1 + 0.3) = e exp(0.3) and exp(0.8) = exp(1 – […]The post Mentally computing e^x first appeared on John D. Cook.

In my post on mentally calculating logarithms, I showed that log10 x ≈ (x – 1)/(x + 1) for 1/√10 ≤ x ≤ √10. You could convert this into an approximation for logs in any base by multiplying by the right scaling factor, but why does it work out so simply for base 10? Define […]The post What makes the log10 trick work? first appeared on John D. Cook.