Blair’s Science Desk

Using a power law to model the number of billionaires

Here’s a brief tutorial about using power laws to numerically model the number of billionaires in a society. We’ll build the model in R.

To get started, let’s install the poweRlaw package:

# install the powerRlaw package
install.packages("poweRlaw")

The poweRlaw package comes with many useful functions for working with power laws. Today, we’ll use the rplcon function, which generates random numbers drawn from a continuous power law.

The rplcon function takes 3 inputs:

  1. n: the size of our random sample of numbers
  2. xmin: the lower threshold for our power law
  3. alpha: the power-law exponent

For example, here’s how we’d sample 10 numbers from a power-law distribution that has a lower threshold of 100 and a power-law exponent of 3:

# example
rplcon(n = 10, xmin = 100, alpha = 3 )

Here’s the type of output we’d get

> rplcon(n = 10, xmin = 100, alpha = 3 )

[1] 108.6708 155.1333 112.7486 116.7339 137.0705 230.2303 127.7278 133.8075 122.2766 131.7206

Counting billionaires

Let’s suppose that our power law describes a distribution of wealth. Given this distribution, we want to know the portion of people who are billionaires. Here’s how we do it.

As, before, we’ll use the rplcon function to generate a power-law sample. Let’s set the minimum wealth to xmin = 1000 and use a power-law exponent of alpha = 2. And let’s suppose we have a population of 10 million people (n = 1e7). Here’s a sample of the distribution of wealth:

# simulate a power-law distribution of wealth
wealth = rplcon( n = 1e7, xmin = 1000, alpha = 2)

Now that we have simulated data for wealth, we can count the number of billionaires. To isolate our billionaires, we’ll use the which function. We’ll take our wealth vector and ask which elements have values greater than or equal to 1 billion. (In scientific notation, 1 billion is 1e9.)

# isolate billionaires
billionaires = which( wealth >= 1e9 )

Next, we’ll calculate the number of billionaires per capita. To do that, we’ll use the length function, which tells us how many elements are in a given vector.

To count the number of billionaires, we’ll use length(billionaires). And to count the whole population, we’ll use length(wealth). The number of billionaires per capita is then the ratio of the two numbers:

# count billionaires
n_billionaires = length( billionaires ) 

# count population
n_population = length( wealth )

# calculate the number of billionaires per capita
billionaires_per_capita = n_billionaires / n_population

Doing the math, I calculate a billionaire density of about 0.8 per million:

> billionaires_per_capita
[1] 8e-07

Since we’re working with random numbers, we’ll get a slightly different billionaire density each time we take the sample. Still, there should be a pattern between our power-law parameters and the number of billionaires. To explore the pattern, play around with different input values and see what you find.