Benford's law: Difference between revisions

Content added Content deleted

Inline

Revision as of 00:15, 4 May 2013

This page uses content from Wikipedia. The original article was at Benford's_law. The list of authors can be seen in the page history. As with Rosetta Code, the text of Wikipedia is available under the GNU FDL. (See links for details on variance)

Benford's law, also called the first-digit law, refers to the frequency distribution of digits in many (but not all) real-life sources of data. In this distribution, the number 1 occurs as the first digit about 30% of the time, while larger numbers occur in that position less frequently: 9 as the first digit less than 5% of the time. This distribution of first digits is the same as the widths of gridlines on a logarithmic scale. Benford's law also concerns the expected distribution for digits beyond the first, which approach a uniform distribution.

This result has been found to apply to a wide variety of data sets, including electricity bills, street addresses, stock prices, population numbers, death rates, lengths of rivers, physical and mathematical constants, and processes described by power laws (which are very common in nature). It tends to be most accurate when values are distributed across multiple orders of magnitude.

A set of numbers is said to satisfy Benford's law if the leading digit d (d ∈ {1, ..., 9}) occurs with probability

For this task, write (a) routine(s) to calculate the distribution of first significant (non-zero) digits in a collection of numbers, then display the actual vs. expected distribution in the way most convenient for your language (table / graph / histogram / whatever).

Use the first 1000 numbers from the Fibonacci sequence as your data set. No need to show how the Fibonacci numbers are obtained. You can generate them or load them from a file; whichever is easiest. Display your actual vs expected distribution.

For extra credit: Show the distribution for one other set of numbers from a page on Wikipedia. State which Wikipedia page it can be obtained from and what the set enumerates. Again, no need to display the actual list of numbers or the code to load them.

See also:

numberphile.com.

Perl 6

<lang perl6>sub benford (@a) { @a.grep(/<[1..9]>/)».match(/<[1..9]>/).bag }

sub dump (%distribution, $base = 10) {

   printf "%9s %9s  %s\n", <Actual Expected Deviation>;
   for 1 .. 9 -> $digit {
       my $actual = %distribution{$digit} * 100 / [+] %distribution.values;
       my $expected = (1 + 1 / $digit).log($base) * 100;
       printf "%d: %5.2f%% | %5.2f%% | %.2f%%\n",
         $digit, $actual, $expected, abs($expected - $actual);
   }

}

( 1, 1, 2, *+* ... * )[^1000].&benford.&dump;</lang>

Output: First 1000 Fibonaccis

   Actual  Expected  Deviation
1: 30.10% | 30.10% | 0.00%
2: 17.70% | 17.61% | 0.09%
3: 12.50% | 12.49% | 0.01%
4:  9.60% |  9.69% | 0.09%
5:  8.00% |  7.92% | 0.08%
6:  6.70% |  6.69% | 0.01%
7:  5.60% |  5.80% | 0.20%
8:  5.30% |  5.12% | 0.18%
9:  4.50% |  4.58% | 0.08%

Extra credit: Square Kilometers of land under cultivation, by country / territory. First column from Wikipedia: Land use statistics by country.

   Actual  Expected  Deviation
1: 33.33% | 30.10% | 3.23%
2: 18.31% | 17.61% | 0.70%
3: 13.15% | 12.49% | 0.65%
4:  8.45% |  9.69% | 1.24%
5:  9.39% |  7.92% | 1.47%
6:  5.63% |  6.69% | 1.06%
7:  4.69% |  5.80% | 1.10%
8:  5.16% |  5.12% | 0.05%
9:  1.88% |  4.58% | 2.70%

Python

Works with Python 3.X & 2.7 <lang python>from __future__ import division from itertools import islice from math import log10 from random import randint

expected = [log10(1+1/d) for d in range(1,10)]

def fib():

   a,b = 1,1
   while True:
       yield a
       a,b = b,a+b

powers of 3 as a test sequence

def power_of_threes():

   s = 1
   while True:
       yield s
       s *= 3

def heads(s):

   for a in s: yield int(str(a)[0])

def show_dist(title, s):

   f,size = [0] * 9, 0
   for x in s:
       f[x-1] += 1
       size += 1
   res = [c/size for c in f]

   print("\n%s Benfords deviation" % title)
   for r, e in zip(res, expected):
       print("%5.1f%% %5.1f%%  %5.1f%%" % (r*100., e*100., abs(r - e)*100.))

def rand1000():

   while True: yield randint(1,9999)

if __name__ == '__main__':

   show_dist("fibbed", islice(heads(fib()), 1000))
   show_dist("threes", islice(heads(power_of_threes()), 1000))

   # just to show that not all kind-of-random sets behave like that
   show_dist("random", islice(heads(rand1000()), 10000))</lang>

Output:

fibbed Benfords deviation
 30.1%  30.1%    0.0%
 17.7%  17.6%    0.1%
 12.5%  12.5%    0.0%
  9.6%   9.7%    0.1%
  8.0%   7.9%    0.1%
  6.7%   6.7%    0.0%
  5.6%   5.8%    0.2%
  5.3%   5.1%    0.2%
  4.5%   4.6%    0.1%

threes Benfords deviation
 30.0%  30.1%    0.1%
 17.7%  17.6%    0.1%
 12.3%  12.5%    0.2%
  9.8%   9.7%    0.1%
  7.9%   7.9%    0.0%
  6.6%   6.7%    0.1%
  5.9%   5.8%    0.1%
  5.2%   5.1%    0.1%
  4.6%   4.6%    0.0%

random Benfords deviation
 11.2%  30.1%   18.9%
 10.9%  17.6%    6.7%
 11.6%  12.5%    0.9%
 11.1%   9.7%    1.4%
 11.6%   7.9%    3.7%
 11.4%   6.7%    4.7%
 10.3%   5.8%    4.5%
 11.0%   5.1%    5.9%
 10.9%   4.6%    6.3%

Scala

<lang scala>// Fibonacci Sequence (begining with 1,1): 1 1 2 3 5 8 13 21 34 55 ... val fibs : Stream[BigInt] = { def series(i:BigInt,j:BigInt):Stream[BigInt] = i #:: series(j, i+j); series(1,0).tail.tail }

/**

* Given a numeric sequence, return the distribution of the most-signicant-digit 
* as expected by Benford's Law and then by actual distribution.
*/

def benford[N:Numeric]( data:Seq[N] ) : Map[Int,(Double,Double)] = {

 import scala.math._
 
 val maxSize = 10000000  // An arbitrary size to avoid problems with endless streams
 
 val size = (data.take(maxSize)).size.toDouble
 
 val distribution = data.take(maxSize).groupBy(_.toString.head.toString.toInt).map{ case (d,l) => (d -> l.size) }
 
 (for( i <- (1 to 9) ) yield { (i -> (log10(1D + 1D / i), (distribution(i) / size))) }).toMap

}

{

 println( "Fibonacci Sequence (size=1000): 1 1 2 3 5 8 13 21 34 55 ...\n" )
 println( "%9s %9s %9s".format( "Actual", "Expected", "Deviation" ) )

 benford( fibs.take(1000) ).toList.sorted foreach { 
   case (k, v) => println( "%d: %5.2f%% | %5.2f%% | %5.4f%%".format(k,v._2*100,v._1*100,math.abs(v._2-v._1)*100) ) 
 }

}</lang>

Output:

Fibonacci Sequence (size=1000): 1 1 2 3 5 8 13 21 34 55 ...

   Actual  Expected Deviation
1: 30.10% | 30.10% | 0.0030%
2: 17.70% | 17.61% | 0.0909%
3: 12.50% | 12.49% | 0.0061%
4:  9.60% |  9.69% | 0.0910%
5:  8.00% |  7.92% | 0.0819%
6:  6.70% |  6.69% | 0.0053%
7:  5.60% |  5.80% | 0.1992%
8:  5.30% |  5.12% | 0.1847%
9:  4.50% |  4.58% | 0.0757%