Monte Carlo Simulation and My Hayley Williams Crush

I've had a crush on Hayley Williams since 2012, my first semester of high school. The first Paramore song I heard was "crushcrushcrush", and I became a fan because of their pop punk sound.

I wanted to learn Monte Carlo simulation, so I gave it a question I like: in how many parallel worlds do Hayley and I live together?

Hayley Williams

Hayley Williams was born on December 27, 1988, in Meridian, Mississippi. She moved to Franklin, Tennessee in 2002, when she was 13, and Paramore formed there in 2004. "crushcrushcrush" is on Riot!, their second album, released on June 12, 2007. She wrote the song with Josh Farro, and it peaked at No. 54 on the Billboard Hot 100. She still sings it, and in February 2025 she sang it with Finneas at the Ryman Auditorium in Nashville.

She and Taylor York won Best Rock Song at the 2015 Grammys for "Ain't It Fun". In 2024, Paramore won Best Rock Album and Best Alternative Music Performance for This Is Why. In 2016 she started the hair color brand Good Dye Young with her longtime stylist Brian O'Connor, and their first product was Riot, her signature orange.

Her third solo album, Ego Death at a Bachelorette Party, came out in 2025 on her own imprint, Post Atlantic. It got four Grammy nominations in 2026. Paramore has been on a break since September 2025. In her words, "We always take huge breaks."

She lives in Nashville, about 16,400 km from Jakarta. Paramore has played Indonesia three times: Bali on August 17, 2011 and Jakarta on August 19, 2011, before I was a fan, and ICE BSD in Tangerang on August 25, 2018.

How Monte Carlo simulation started

In 1946, the mathematician Stanislaw Ulam was recovering from an illness and playing solitaire. He asked "what are the chances that a Canfield solitaire laid out with 52 cards will come out successfully?" He spent a lot of time on pure combinatorial calculations. Then he wondered whether a more practical method was "to lay it out say one hundred times and simply observe and count the number of successful plays."

That idea is a Monte Carlo simulation: you estimate a probability or an average by running many random trials and counting the results.

Ulam described the idea to John von Neumann later in 1946. On March 11, 1947, von Neumann sent a letter to Robert Richtmyer at Los Alamos with a statistical approach to neutron diffusion in fissionable material. Nicholas Metropolis wrote in the same article that he suggested the name, because "Stan had an uncle who would borrow money from relatives because he 'just had to go to Monte Carlo.'" The first Monte Carlo calculations on a computer ran on ENIAC in April 1948, and Metropolis and Ulam published The Monte Carlo Method in 1949.

Los Alamos still maintains MCNP, a Monte Carlo code for neutrons, photons, and other particles. Physically based renderers use Monte Carlo because the light transport equations generally have no analytic solution. AlphaGo combined Monte Carlo tree search with neural networks to beat the European Go champion 5 to 0.

Estimating pi with random points

How can random numbers give you pi?

Draw a square with sides of length 1. Inside it, draw a quarter circle with radius 1, centered on the bottom left corner. The square has an area of 1, and the quarter circle has an area of pi / 4. If you drop a point at a random place in the square, it lands inside the quarter circle with probability pi / 4.

A point (x, y) is inside the quarter circle when x² + y² ≤ 1. So you drop many points, count the ones inside, and multiply the fraction by 4.

JavaScript's Math.random() doesn't let you set the seed, so every run prints different numbers. A seed is the starting value of a random number generator, and the same seed always gives the same sequence. I want your output to match mine, so the code uses sfc32, a small generator by Chris Doty-Humphrey that's in the public domain. The JavaScript version comes from bryc's PRNG collection.

A generator with only 32 bits of state, such as mulberry32, is too small here. Art Owen's guideline is to use no more than about the square root of a generator's period. For a period of 2³², that's about 65,536 numbers, and a run of 100 million worlds draws more than a billion. sfc32 keeps 128 bits of state.

// random.ts
export function sfc32(seed: number) {
  let a = 0;
  let b = seed | 0;
  let c = 0;
  let d = 1;
  const next = (): number => {
    const t = (((a + b) | 0) + d) | 0;
    d = (d + 1) | 0;
    a = b ^ (b >>> 9);
    b = (c + (c << 3)) | 0;
    c = (c << 21) | (c >>> 11);
    c = (c + t) | 0;
    return (t >>> 0) / 4294967296;
  };
  // Mix the seed before the first use.
  for (let i = 0; i < 12; i++) next();
  return next;
}

next() returns a number from 0 up to, but not including, 1. PractRand's own seeding code also runs the generator 12 times before use.

// pi.ts
import { sfc32 } from "./random.ts";

const random = sfc32(2012);

function estimatePi(n: number): number {
  let inside = 0;
  for (let i = 0; i < n; i++) {
    const x = random();
    const y = random();
    if (x * x + y * y <= 1) inside++;
  }
  return (4 * inside) / n;
}

const sizes = [100, 10_000, 1_000_000];
for (const n of sizes) {
  const pi = estimatePi(n);
  const error = Math.abs(pi - Math.PI);
  console.log(n, pi, error.toFixed(5));
}

$ node pi.ts
100 3.04 0.10159
10000 3.1568 0.01521
1000000 3.142596 0.00100

I ran these files with Node 24.19, which runs TypeScript files without a build step. The seed 2012 is the year I heard "crushcrushcrush".

Why more points give a better estimate

The law of large numbers says that the average of many independent random trials gets closer to the true value as the number of trials grows. Each point in estimatePi() is one trial that scores 4 or 0, and the estimate is the average score.

The typical error of a Monte Carlo estimate is the standard deviation of one trial divided by the square root of the number of trials, so it shrinks slowly. One more correct decimal digit needs 100 times more trials. For the pi estimate, the standard deviation of one trial is 4 times the square root of (pi / 4)(1 − pi / 4), about 1.64.

Each row has 100 times more points than the row above it, and the typical error gets 10 times smaller. The actual error is random, so it can land above or below the typical error.

Simulating the worlds

Each world starts in 2012, when I hear "crushcrushcrush". Every year until 2026, the world can take the next step in a chain. The chances are my own guesses.

A world takes at most one step a year, and a step that fails one year can happen in a later year. So the answer depends on when each step happens, and multiplying the six chances doesn't answer the question.

// worlds.ts
import { sfc32 } from "./random.ts";

const random = sfc32(2012);

const steps: [string, number][] = [
  ["see Paramore live", 1 / 20],
  ["meet her after a show", 1 / 1000],
  ["she remembers me", 1 / 10],
  ["we become friends", 1 / 5],
  ["I move to Nashville", 1 / 50],
  ["we live together", 1 / 100],
];

const YEARS = 15; // 2012 to 2026

function liveOneWorld(): number {
  let step = 0;
  for (let y = 0; y < YEARS; y++) {
    if (step === steps.length) break;
    const [, odds] = steps[step];
    if (random() < odds) step++;
  }
  return step;
}

function simulate(n: number): number[] {
  const reached = steps.map(() => 0);
  for (let i = 0; i < n; i++) {
    const last = liveOneWorld();
    for (let s = 0; s < last; s++) {
      reached[s]++;
    }
  }
  return reached;
}

const n = Number(process.argv[2]);
const reached = simulate(n);
steps.forEach(([name], i) => {
  console.log(name, reached[i]);
});

liveOneWorld() returns how many steps one world reached. simulate() counts, for each step, how many worlds got that far.

$ node worlds.ts 1000000
see Paramore live 537100
meet her after a show 4209
she remembers me 1483
we become friends 646
I move to Nashville 38
we live together 1

Out of a million worlds, 537,100 have me at a Paramore show. One world has Hayley and me living together.

Rare events need more worlds

Does one world out of a million mean the chance is 1 in a million?

It might not. When the event is rare, the error depends on the number of hits you count, not on the number of worlds you run. The relative error is about 1 divided by the square root of the expected number of hits. With about one hit, the relative error is about 100 percent, so one hit can't tell 1 in 500,000 apart from 1 in 5 million.

To get the relative error down to about 10 percent, I need about 100 hits. For an event near 1 in a million, that means about 100 million worlds.

$ node worlds.ts 100000000
see Paramore live 53674077
meet her after a show 424208
she remembers me 146467
we become friends 64906
I move to Nashville 3212
we live together 66

That run took 14 seconds on my Apple M5 laptop. 66 worlds out of 100 million is about 1 in 1.5 million, with a relative error of about 12 percent.

This chain is small enough to calculate exactly, year by year, so I checked the simulation against it. The exact answer is 1 in 1,617,113. Over 100 million worlds, that predicts about 62 hits, and 66 is within the error.

What the number means

The 1 in 1.6 million comes from chances I made up. If I change "I move to Nashville" from 1 in 50 a year to 1 in 5, the exact answer becomes about 1 in 205,000, close to 8 times higher. The simulation counts carefully, but it can't make my guesses more accurate.

In the many worlds picture, the weight of a branch comes from quantum physics, not from a count of imagined worlds, so my count says nothing about real branches.

In this world, I still haven't seen Paramore live, and the band is on a break. I hope the first step of my chain happens in this one.