ベンフォードの法則を確かめてみる

Python 30分プログラム

30分プログラム、その584。素数の分布はベンフォードの法則に従っているらしいので試してみる。

「ベンフォードの法則」とは、ある数値群をみたとき、最高桁が「1」である数値は（15や189や1088など）は全体の約30%、「2」であるものは約18%、「3」であるもの約12%･･･「9」であるものは約5%という割合になっているという法則である。この法則は物理学者フランク・ベンフォードが1938年に発見した分布法則であり、市場分析や不正検出アルゴリズムなどにも応用されているものである。

スペインの数学者がBartolo LuqueとLucas Lacasaは素数にもこの法則が当てはまることをこの度解明したそうだ（論文要旨）。

というわけで、これを最高桁の分布を調べるプログラムを書いてみた。

とりあえず5000までの素数で確認したけど、特にベンフォードの法則には従っていなかった。まあ、そんなに簡単に分かることだったら論文にもならないだろうし、当然な気もする。

使い方

$ python benford.py
1: 160 23.92 %
2: 146 21.82 %
3: 139 20.78 %
4: 139 20.78 %
5: 17 2.54 %
6: 18 2.69 %
7: 18 2.69 %
8: 17 2.54 %
9: 15 2.24 %

ソースコード

#! /usr/bin/python
# -*- mode:python; coding:utf-8 -*-
#
# benford.py -
#
# Copyright(C) 2009 by mzp
# Author: MIZUNO Hiroki / mzpppp at gmail dot com
# http://howdyworld.org
#
# Timestamp: 2009/05/15 21:55:28
#
# This program is free software; you can redistribute it and/or
# modify it under MIT Lincence.
#

def sieve(xs):
    if xs == []:
        return []
    else:
        x  = xs[0]
        return [x] + sieve([y for y in xs if y % x != 0])

def primes(n):
    return sieve(xrange(2,n))

def left_most(n):
    return int(str(n)[0])

def count(xs):
    statics = [0] * 9
    for x in xs:
        statics[left_most(x)-1] += 1
    return statics

def show(xs):
    total = sum(xs)
    for (i,x) in zip(xrange(0,9),xs):
        print "%d: %d %.2f %%" % (i+1,x,100* x / float(total))

def benford(xs):
    show(count(xs))

benford(primes(5000))

参考