Benford's law: Difference between revisions
Add ABC
(Add F# version) |
Not a robot (talk | contribs) (Add ABC) |
||
(45 intermediate revisions by 23 users not shown) | |||
Line 31:
* A starting page on Wolfram Mathworld is {{Wolfram|Benfords|Law}}.
<br><br>
=={{header|8th}}==
<
: n:log10e ` 1 10 ln / ` ;
Line 142 ⟶ 101:
benford-test
bye
</syntaxhighlight>
{{out}}
Line 157 ⟶ 116:
9 0.0458 0.045
</pre>
=={{header|11l}}==
{{trans|D}}
<syntaxhighlight lang="11l">F get_fibs()
V a = 1.0
V b = 1.0
[Float] r
L 1000
r [+]= a
(a, b) = (b, a + b)
R r
F benford(seq)
V freqs = [(0.0, 0.0)] * 9
V seq_len = 0
L(d) seq
I d != 0
freqs[String(d)[0].code - ‘1’.code][1]++
seq_len++
L(&f) freqs
f = (log10(1.0 + 1.0 / (L.index + 1)), f[1] / seq_len)
R freqs
print(‘#9 #9 #9’.format(‘Actual’, ‘Expected’, ‘Deviation’))
L(p) benford(get_fibs())
print(‘#.: #2.2% | #2.2% | #.4%’.format(L.index + 1, p[1] * 100, p[0] * 100, abs(p[1] - p[0]) * 100))</syntaxhighlight>
{{out}}
<pre>
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%
</pre>
=={{header|ABC}}==
<syntaxhighlight lang="abc">HOW TO RETURN fibonacci.numbers n:
PUT 1, 1 IN a, b
PUT {} IN fibo
FOR i IN {1..n}:
INSERT a IN fibo
PUT b, a+b IN a, b
RETURN fibo
HOW TO RETURN digit.distribution nums:
PUT {} IN digits
FOR i IN {1..9}: PUT i IN digits["`i`"]
PUT {} IN dist
FOR i IN {1..9}: PUT 0 IN dist[i]
FOR n IN nums:
PUT digits["`n`"|1] IN digit
PUT dist[digit] + 1 IN dist[digit]
FOR i IN {1..9}:
PUT dist[i] / #nums IN dist[i]
RETURN dist
PUT digit.distribution fibonacci.numbers 1000 IN observations
WRITE "Digit"<<6, "Expected">>10, "Observed">>10/
FOR d IN {1..9}:
WRITE d<<6, ((10 log (1 + 1/d))>>10)|10, observations[d]>>10/</syntaxhighlight>
{{out}}
<pre>Digit Expected Observed
1 0.30102999 0.301
2 0.17609125 0.177
3 0.12493873 0.125
4 0.09691001 0.096
5 0.07918124 0.08
6 0.06694678 0.067
7 0.05799194 0.056
8 0.05115252 0.053
9 0.04575749 0.045</pre>
=={{header|Ada}}==
Line 162 ⟶ 198:
The program reads the Fibonacci-Numbers from the standard input. Each input line is supposed to hold N, followed by Fib(N).
<
procedure Benford is
Line 210 ⟶ 246:
Ada.Text_IO.New_Line;
end loop;
end Benford;</
{{out}}
Line 230 ⟶ 266:
Input is the list of primes below 100,000 from [http://www.mathsisfun.com/numbers/prime-number-lists.html]. Since each line in that file holds prime and only a prime, but no ongoing counter, we must slightly modify the program by commenting out a single line:
<
We can also edit out the declaration of the variable "Counter" ...or live with a compiler warning about never reading or assigning that variable.
Line 251 ⟶ 287:
=={{header|Aime}}==
<
sum(text a, text b)
{
Line 338 ⟶ 374:
0;
}</
{{out}}
<pre> expected found
Line 350 ⟶ 386:
8 5.115 5.300
9 4.575 4.5 </pre>
=={{header|ALGOL 68}}==
{{works with|ALGOL 68G|Any - tested with release 2.8.3.win32}}
Uses Algol 68G's LONG LONG INT which has programmer specifiable precision.
<
# set the number of digits for LONG LONG INT values #
PR precision 256 PR
Line 380 ⟶ 415:
PROC compare to benford = ( []REAL actual )VOID:
FOR i TO 9 DO
print( ( "Benford: ", fixed( log( 1 + ( 1 / i ) ), -7, 3
, " actual: ", fixed( actual[ i ], -7, 3 )
, newline
)
)
OD # compare to benford # ;
# generate 1000 fibonacci numbers #
Line 390 ⟶ 429:
# compare to the probabilities expected by Benford's law #
compare to benford( digit probability( fn ) )
END</
{{out}}
<pre>
Line 403 ⟶ 442:
Benford: 0.046 actual: 0.045
</pre>
=={{header|APL}}==
{{works with|Dyalog APL}}
<syntaxhighlight lang="apl">task←{
benf ← ≢÷⍨(⍳9)(+/∘.=)(⍎⊃∘⍕)¨
fibs ← (⊢,(+/¯2↑⊢))⍣998⊢1 1
exp ← 10⍟1+÷⍳9
obs ← benf fibs
⎕←'Expected Actual'⍪5⍕exp,[1.5]obs
}</syntaxhighlight>
{{out}}
<pre>Expected Actual
0.30103 0.30100
0.17609 0.17700
0.12494 0.12500
0.09691 0.09600
0.07918 0.08000
0.06695 0.06700
0.05799 0.05600
0.05115 0.05300
0.04576 0.04500</pre>
=={{header|Arturo}}==
<syntaxhighlight lang="arturo">fib: [a b]: <= [0 1]
do.times:998 -> [a b]: @[b 'fib ++ <= a+b]
leading: fib | map 'x -> first ~"|x|"
| tally
print "digit actual expected"
loop 1..9 'd ->
print [
pad.right ~"|d|" 8
pad.right ~"|leading\[d] // 1000|" 9
log 1 + 1//d 10
]</syntaxhighlight>
{{out}}
<pre>digit actual expected
1 0.301 0.3010299956639811
2 0.177 0.1760912590556812
3 0.125 0.1249387366082999
4 0.095 0.09691001300805641
5 0.08 0.0791812460476248
6 0.067 0.06694678963061322
7 0.056 0.05799194697768673
8 0.053 0.05115252244738128
9 0.045 0.04575749056067514</pre>
=={{header|AutoHotkey}}==
{{works with|AutoHotkey_L}}(AutoHotkey1.1+)
<
fib := NStepSequence(1, 1, 2, 1000)
Out := "Digit`tExpected`tObserved`tDeviation`n"
Line 448 ⟶ 540:
}
return, (Carry ? Carry : "") . Result
}</
'''Output:'''
<pre>Digit Expected Observed Deviation
Line 462 ⟶ 554:
=={{header|AWK}}==
<syntaxhighlight lang="awk">
# syntax: GAWK -f BENFORDS_LAW.AWK
BEGIN {
Line 490 ⟶ 582:
function abs(x) { if (x >= 0) { return x } else { return -x } }
function log10(x) { return log(x)/log(10) }
</syntaxhighlight>
{{out}}
<pre>
Line 504 ⟶ 596:
9 4.5757 4.5000 0.0757
</pre>
=={{header|BCPL}}==
BCPL doesn't do floating point well, so I use integer math to compute the most significant digits of the Fibonacci sequence and use a table that has the values of log10(d + 1/d)
<syntaxhighlight lang="bcpl">
GET "libhdr"
Line 549 ⟶ 640:
RESULTIS 0
}
</syntaxhighlight>
{{Out}}
<pre>
Line 566 ⟶ 657:
9 0.045 0.046
</pre>
=={{header|BASIC256}}==
<syntaxhighlight lang="vb">n = 1000
dim actual(n) fill 0
for nr = 1 to n
num$ = string(fibonacci(nr))
j = int(left(num$,1))
actual[j] += 1
next
print "First 1000 Fibonacci numbers"
print "Digit ", "Actual freq ", "Expected freq"
for i = 1 to 9
freq = frequency(i)*100
print " "; ljust(i,4), rjust(actual[i]/10,5), rjust(freq,5)
next
end
function frequency(n)
return (log10(n+1) - log10(n))
end function
function fibonacci(f)
f = int(f)
a = 0 : b = 1 : c = 0 : n = 0
while n < f
a = b
b = c
c = a + b
n += 1
end while
return c
end function</syntaxhighlight>
{{Out}}
<pre>First 1000 Fibonacci numbers
Digit Actual freq Expected freq
1 30.1 30.1029995664
2 17.7 17.6091259056
3 12.5 12.4938736608
4 9.6 9.69100130081
5 8.0 7.91812460476
6 6.7 6.694679
7 5.6 5.79919469777
8 5.3 5.11525224474
9 4.5 4.57574905607</pre>
=={{header|C}}==
<
#include <stdlib.h>
#include <math.h>
Line 631 ⟶ 770:
return EXIT_SUCCESS;
}</
{{Out}}
Line 646 ⟶ 785:
8 0.053 0.051
9 0.045 0.046</pre>
=={{header|C++}}==
<
//if used prepackaged you can compile writing "g++ -std=c++11 -lcln yourprogram.cpp -o yourprogram"
#include <cln/integer.h>
Line 710 ⟶ 848:
return 0 ;
}
</syntaxhighlight>
{{out}}
<pre> found expected
Line 723 ⟶ 861:
9 : 4.5 % 4.58 %
</pre>
=={{header|Chipmunk Basic}}==
{{trans|Yabasic}}
{{works with|Chipmunk Basic|3.6.4}}
<syntaxhighlight lang="vbnet">100 cls
110 n = 1000
120 dim actual(n)
130 for nr = 1 to n
140 num$ = str$(fibonacci(nr))
150 j = val(left$(num$,1))
160 actual(j) = actual(j)+1
170 next
180 print "First 1000 Fibonacci numbers"
190 print "Digit Actual freq Expected freq"
200 for i = 1 to 9
210 freq = frequency(i)*100
220 print format$(i,"###");
230 print using " ##.###";actual(i)/10;
240 print using " ##.###";freq
250 next
260 end
270 sub frequency(n)
280 frequency = (log10(n+1)-log10(n))
290 end sub
300 sub log10(n)
310 log10 = log(n)/log(10)
320 end sub
330 sub fibonacci(n)
335 rem https://rosettacode.org/wiki/Fibonacci_sequence#Chipmunk_Basic
340 n1 = 0
350 n2 = 1
360 for k = 1 to abs(n)
370 sum = n1+n2
380 n1 = n2
390 n2 = sum
400 next k
410 if n < 0 then
420 fibonacci = n1*((-1)^((-n)+1))
430 else
440 fibonacci = n1
450 endif
460 end sub</syntaxhighlight>
=={{header|Clojure}}==
<
(:gen-class))
Line 799 ⟶ 979:
(show-benford-stats y))
</syntaxhighlight>
{{Output}}
<pre>
Line 836 ⟶ 1,016:
9 5.76 4.58 1.18
</pre>
=={{header|CoffeeScript}}==
<
a = 1; b = 0
return () ->
Line 858 ⟶ 1,037:
console.log "Digit\tActual\tExpected"
for i in [1..9]
console.log i + "\t" + actual[i - 1].toFixed(3) + '\t' + expected[i - 1].toFixed(3)</
{{out}}
<pre>Leading digital distribution of the first 1,000 Fibonacci numbers
Line 871 ⟶ 1,050:
8 0.053 0.051
9 0.045 0.046</pre>
=={{header|Common Lisp}}==
<
"Return the frequency distribution of the most significant nonzero
digits in the given list of numbers. The first element of the list
Line 915 ⟶ 1,093:
(map 'list #'list '(1 2 3 4 5 6 7 8 9)
actual-distribution
expected-distribution))))</
<pre>; *fib1000* is a list containing the first 1000 numbers in the Fibonnaci sequence
Line 929 ⟶ 1,107:
8 0.053 0.051
9 0.045 0.046</pre>
=={{header|Crystal}}==
{{trans|Ruby}}
<
EXPECTED = (1..9).map{ |d| Math.log10(1 + 1.0 / d) }
Line 970 ⟶ 1,147:
# just to show that not all kind-of-random sets behave like that
show_dist("random", random(10000))</
{{out}}
Line 1,007 ⟶ 1,184:
9: 12.3% 4.6% 7.7%
</pre>
=={{header|D}}==
{{trans|Scala}}
<
double[2][9] benford(R)(R seq) if (isForwardRange!R && !isInfinite!R) {
Line 1,033 ⟶ 1,209:
writefln("%d: %5.2f%% | %5.2f%% | %5.4f%%",
i+1, p[1] * 100, p[0] * 100, abs(p[1] - p[0]) * 100);
}</
{{out}}
<pre> Actual Expected Deviation
Line 1,049 ⟶ 1,225:
===Alternative Version===
The output is the same.
<
std.algorithm, std.array;
Line 1,066 ⟶ 1,242:
f * 100.0 / N, expected[i] * 100,
abs((f / double(N)) - expected[i]) * 100);
}</
=={{header|Delphi}}==
See [
=={{header|EasyLang}}==
<syntaxhighlight lang=easylang>
func$ add a$ b$ .
for i to higher len a$ len b$
a = number substr a$ i 1
b = number substr b$ i 1
r = a + b + c
c = r div 10
r$ &= r mod 10
.
if c > 0
r$ &= c
.
return r$
.
#
len fibdist[] 9
proc mkfibdist . .
# generate 1000 fibonacci numbers as
# (reversed) strings, because 53 bit
# integers are too small
#
n = 1000
prev$ = 0
val$ = 1
fibdist[1] = 1
for i = 2 to n
h$ = add prev$ val$
prev$ = val$
val$ = h$
ind = number substr val$ len val$ 1
fibdist[ind] += 1
.
for i to len fibdist[]
fibdist[i] = fibdist[i] / n
.
.
mkfibdist
#
len benfdist[] 9
proc mkbenfdist . .
for i to 9
benfdist[i] = log10 (1 + 1.0 / i)
.
.
mkbenfdist
#
numfmt 3 0
print "Actual Expected"
for i to 9
print fibdist[i] & " " & benfdist[i]
.
</syntaxhighlight>
{{out}}
<pre>
Actual Expected
0.301 0.301
0.177 0.176
0.125 0.125
0.096 0.097
0.080 0.079
0.067 0.067
0.056 0.058
0.053 0.051
0.045 0.046
</pre>
=={{header|Elixir}}==
<
def distribution(n), do: :math.log10( 1 + (1 / n) )
Line 1,087 ⟶ 1,331:
end
Benfords_law.task</
{{out}}
Line 1,104 ⟶ 1,348:
=={{header|Erlang}}==
<syntaxhighlight lang="erlang">
-module( benfords_law ).
-export( [actual_distribution/1, distribution/1, task/0] ).
Line 1,130 ⟶ 1,374:
[Key | _] = erlang:integer_to_list( N ),
dict:update_counter( Key - 48, 1, Dict ).
</syntaxhighlight>
{{out}}
<pre>
Line 1,145 ⟶ 1,389:
9 0.045 0.04575749056067514
</pre>
=={{header|F sharp|F#}}==
For Fibonacci code, see https://rosettacode.org/wiki/Fibonacci_sequence#F.89
<
let fibonacci = Seq.unfold (fun (x, y) -> Some(x, (y, x + y))) (0I,1I)
Line 1,170 ⟶ 1,413:
printfn "\nBenford's law for 1 through 9:"
benfordLawFigures |> List.iter (fun f -> printf $"{f:N5} ")
</syntaxhighlight>
{{out}}
Line 1,180 ⟶ 1,423:
</pre>
=={{header|Factor}}==
<
kernel math math.functions math.statistics math.text.utils
sequences ;
Line 1,210 ⟶ 1,452:
.header leading histogram [ .digit-report ] assoc-each ;
MAIN: main</
{{out}}
<pre>
Line 1,224 ⟶ 1,466:
9 0.0458 0.0450
</pre>
=={{header|Forth}}==
<
: f2drop fdrop fdrop ;
Line 1,269 ⟶ 1,510:
set-precision ;
: compute-benford tally report ;</
{{Out}}
<pre>Gforth 0.7.2, Copyright (C) 1995-2008 Free Software Foundation, Inc.
Line 1,286 ⟶ 1,527:
8 0.053 0.0512
9 0.045 0.0458 ok</pre>
=={{header|Fortran}}==
FORTRAN 90. Compilation and output of this program using emacs compile command and a fairly obvious Makefile entry:
<
Compilation started at Sat May 18 01:13:00
Line 1,297 ⟶ 1,537:
0.300999999 0.177000001 0.125000000 9.60000008E-02 7.99999982E-02 6.70000017E-02 5.70000000E-02 5.29999994E-02 4.50000018E-02 LEADING FIBONACCI DIGIT
Compilation finished at Sat May 18 01:13:00</
<
! compute most significant digits, Fibonacci like.
implicit none
Line 1,367 ⟶ 1,607:
write(6,*) (benfordsLaw(i),i=1,9),'THE LAW'
write(6,*) (count(i)/1000.0 ,i=1,9),'LEADING FIBONACCI DIGIT'
end program benford</
=={{header|FreeBASIC}}==
{{libheader|GMP}}
<
' compile with: fbc -s console
Line 1,460 ⟶ 1,699:
Print : Print "hit any key to end program"
Sleep
End</
{{out}}
<pre>First 1000 Fibonacci numbers
Line 1,488 ⟶ 1,727:
8 71038 7.10 % 5.12 % -1.989 %
9 70320 7.03 % 4.58 % -2.456 %</pre>
=={{header|Fōrmulæ}}==
{{FormulaeEntry|page=https://formulae.org/?script=examples/Benford%27s_law}}
'''Solution:'''
The following function calculates the distribution of the first digit of a list of integer, positive numbers:
[[File:Fōrmulæ - Benford's law 01.png]]
'''Example:'''
[[File:Fōrmulæ - Benford's law 02.png]]
[[File:Fōrmulæ - Benford's law 03.png]]
'''Benford distribution.''' Benford distribution is, by definition (using a precision of 10 digits):
[[File:Fōrmulæ - Benford's law 04.png]]
[[File:Fōrmulæ - Benford's law 05.png]]
'''Comparison chart.''' The following function creates a chart, in order to compare the distribution of the first digis in a givel list or numbers, and the Benford distribution:
[[File:Fōrmulæ - Benford's law 06.png]]
'''Testing.''' Testing with the previous list of numbers:
[[File:Fōrmulæ - Benford's law 07.png]]
[[File:Fōrmulæ - Benford's law 08.png]]
'''Case 1.''' Use the first 1,000 numbers from the Fibonacci sequence as your data set
[[File:Fōrmulæ - Benford's law 09.png]]
[[File:Fōrmulæ - Benford's law 10.png]]
'''Case 2.''' The sequence of the fisrt 1,000 natural numbers does not follow the Benford distribution:
[[File:Fōrmulæ - Benford's law 11.png]]
[[File:Fōrmulæ - Benford's law 12.png]]
'''Case 3.''' The following example is for the list of the fist 1,000 factorial numbers.
[[File:Fōrmulæ - Benford's law 13.png]]
[[File:Fōrmulæ - Benford's law 14.png]]
=={{header|FutureBasic}}==
<syntaxhighlight lang="futurebasic">
Short t, i, j, k, m
Double a(9), z
Double phi, psi
CFStringRef s
print @"Benford:"
for i = 1 to 9
a(i) = log10( 1 + 1 / i )
print fn StringWithFormat( @"%.3f ", a(i) ),
next
// Fibonacci according to DeMoivre and Binet
for t = 1 to 9 : a(t) = 0 : next // Clean array
phi = ( 1 + sqr(5) ) / 2
psi = ( 1 - sqr(5) ) / 2
for i = 1 to 1000
z = ( phi^i - psi^i ) / sqr( 5 )
s = fn StringWithFormat( @"%e", z) // Get first digit
t = fn StringIntegerValue( left( s, 1 ) )
a(t) = a(t) + 1
next
print @"\n\nFibonacci:"
for i = 1 to 9
print fn StringWithFormat( @"%.3f ", a(i) / 1000 ),
next
// Multiplication tables
for t = 1 to 9 : a(t) = 0 : next // Clean array
for i = 1 to 10
for j = 1 to 10
for k = 1 to 10
for m = 1 to 10
z = i * j * k * m
s = fn StringWithFormat( @"%e", z )
t = fn StringIntegerValue( left( s, 1 ) )
a(t) = a(t) + 1
next
next
next
next
print @"\n\nMultiplication:"
for i = 1 to 9
print fn StringWithFormat( @"%.3f ", a(i) / 1e4 ),
next
// Factorials according to DeMoivre and Stirling
for t = 1 to 9 : a(t) = 0 : next // Clean array
for i = 10 to 110
z = sqr(2 * pi * i ) * (i / exp(1) )^i
s = fn StringWithFormat( @"%e", z )
t = fn StringIntegerValue( left( s, 1 ) )
a(t) = a(t) + 1
next
print @"\n\nFactorials:"
for i = 1 to 9
print fn StringWithFormat( @"%.2f ", a(i) / 100 ),
next
handleevents
}</syntaxhighlight>
{{Out}}
<pre>
Benford:
0.301 0.176 0.125 0.097 0.079 0.067 0.058 0.051 0.046
Fibonacci:
0.301 0.177 0.125 0.096 0.080 0.067 0.056 0.053 0.045
Multiplication:
0.302 0.181 0.124 0.103 0.071 0.061 0.055 0.055 0.047
Factorials:
0.35 0.16 0.12 0.06 0.06 0.06 0.02 0.10 0.08
</pre>
=={{header|Go}}==
<
import (
Line 1,528 ⟶ 1,891:
math.Log10(1+1/float64(i+1)))
}
}</
{{Out}}
<pre>
Line 1,543 ⟶ 1,906:
9 0.045 0.046
</pre>
=={{header|Groovy}}==
'''Solution:'''<br>
Uses [[Fibonacci_sequence#Analytic_8|Fibonacci sequence analytic formula]]
{{trans|Java}}
<
def population = (0..<size).collect { generator(it) }
def firstDigits = [0]*10
Line 1,555 ⟶ 1,917:
}
firstDigits
}</
'''Test:'''
<
println "d actual predicted"
(1..<10).each {
printf ("%d %10.6f %10.6f\n", it, digitCounts[it]/1000, Math.log10(1.0 + 1.0/it))
}</
'''Output:'''
Line 1,575 ⟶ 1,937:
8 0.053000 0.051153
9 0.045000 0.045757</pre>
=={{header|GW-BASIC}}==
{{improve|GW-BASIC|Does not show actual vs expected values for Fibonaccis as task specifies.}}
<syntaxhighlight lang="GW-BASIC">
10 DEF FNBENFORD(N)=LOG(1+1/N)/LOG(10)
20 CLS
30 PRINT "One digit Benford's Law"
40 FOR i = 1 TO 9
50 PRINT i,:PRINT USING "##.######";FNBENFORD(i)
60 NEXT i
70 END
</syntaxhighlight>
{{out}}
<pre>
1 0.301030
2 0.176091
3 0.124939
4 0.096910
5 0.079181
6 0.066947
7 0.057992
8 0.051153
9 0.045758
</pre>
=={{header|Haskell}}==
<
import Data.Char (digitToInt)
Line 1,598 ⟶ 1,984:
| d <- [1 .. 9] ]
main = print tab</
{{out}}
<pre>[(1,0.301,0.301029995663981),
Line 1,610 ⟶ 1,996:
(9,0.045,0.0457574905606751)]
</pre>
=={{header|Icon}} and {{header|Unicon}}==
The following solution works in both languages.
<
procedure main()
Line 1,648 ⟶ 2,033:
if n%2 = 1 then return [c+d, d]
else return [d, c]
end</
Sample run:
Line 1,665 ⟶ 2,050:
->
</pre>
=={{header|J}}==
'''Solution'''
<syntaxhighlight lang="j">benford=: 10&^.@(1 + %) NB. expected frequency of first digit y
Digits=: '123456789'
firstSigDigits=: {.@(-. -.&Digits)@":"0 NB. extract first significant digit from numbers
freq=: (] % +/)@:<:@(#/.~)@, NB. calc frequency of values (x) in y</syntaxhighlight>
'''Required Example'''
<syntaxhighlight lang="j"> First1000Fib=: (, +/@:(_2&{.)) ^: (1000-#) 1 1
NB. Expected vs Actual frequencies for Digits 1-9
Digits ((] ,. benford)@"."0@[ ,. (freq firstSigDigits)) First1000Fib
1 0.30103 0.301
2 0.176091 0.177
3 0.124939 0.125
4 0.09691 0.096
5 0.0791812 0.08
6 0.0669468 0.067
7 0.0579919 0.056
8 0.0511525 0.053
9 0.0457575 0.045</syntaxhighlight>
'''Alternatively'''
We show the correlation coefficient of Benford's law with the leading digits of the first 1000 Fibonacci numbers is almost unity.
<
benford =: log10@:(1+%)
assert '0.30 0.18 0.12 0.10 0.08 0.07 0.06 0.05 0.05' -: 5j2 ": benford >: i. 9
Line 1,703 ⟶ 2,109:
assert '0.9999' -: 6j4 ": (normalize TALLY_BY_KEY) r benford >: i.9
assert '0.9999' -: 6j4 ": TALLY_BY_KEY r benford >: i.9 NB. Of course we don't need normalization</
=={{header|Java}}==
<
import java.util.Locale;
Line 1,734 ⟶ 2,140:
}
}
}</
The output is:
<pre>1 0.301000 0.301030
Line 1,746 ⟶ 2,152:
9 0.045000 0.045757</pre>
To use other number sequences, implement a suitable <tt>NumberGenerator</tt>, construct a <tt>Benford</tt> instance with it and print it.
=={{header|JavaScript}}==
<
.reduce(
(fib, _, i) => i < 2 ? (
Line 1,766 ⟶ 2,171:
]);
console.log(benford(fibseries(1000)))</
{{output}}
<pre>0: (3) [1, 0.301, 0.3010299956639812]
Line 1,777 ⟶ 2,182:
7: (3) [8, 0.053, 0.05115252244738129]
8: (3) [9, 0.045, 0.04575749056067514]</pre>
=={{header|jq}}==
{{works with|jq|1.4}}
This implementation shows the observed and expected number of occurrences together with the χ² statistic.
For the sake of being self-contained, the following includes a generator for Fibonacci numbers, and a prime number generator that is inefficient but brief and can generate numbers within an arbitrary range.<
# Numerical accuracy is insufficient beyond about 1450.
def fibonacci(n):
Line 1,888 ⟶ 2,292:
;
task</
{{out}}
<pre>First 100 fibonacci numbers:
Line 1,959 ⟶ 2,363:
χ² = 3204.8072</pre>
=={{header|Julia}}==
<syntaxhighlight lang="julia"># Benford's Law
P(d) = log10(1+1/d)
function benford(numbers)
firstdigit(n) = last(digits(n))
counts = zeros(9)
foreach(n -> counts[firstdigit(n)] += 1, numbers)
counts ./ sum(counts)
end
struct Fibonacci end
Base.iterate(::Fibonacci, (a, b) = big.((0, 1))) = b, (b, a + b)
sample = Iterators.take(Fibonacci(), 1000)
observed = benford(sample) .* 100
expected = P.(1:9) .* 100
table = Real[1:9 observed expected]
using Plots
plot([observed expected]; title = "Benford's Law",
seriestype = [:bar :line], linewidth = [0 5],
xticks = 1:9, xlabel = "first digit", ylabel = "frequency %",
label = ["1000 Fibonacci numbers" "P(d)=log10(1+1/d)"])
using Printf
println("Benford's Law\nFrequency of first digit\nin 1000 Fibonacci numbers")
println("digit observed expected")
foreach(i -> @printf("%3d%9.2f%%%9.2f%%\n", table[i,:]...), 1:9)</syntaxhighlight>
{{Out}}
<pre>Benford's Law
Frequency of first digit
in 1000 Fibonacci numbers
digit observed expected
1 30.10% 30.10%
8 5.30% 5.12%
9 4.50% 4.58%</pre>
=={{header|Kotlin}}==
<
interface NumberGenerator {
Line 2,020 ⟶ 2,448:
}
fun main(a: Array<String>) = println(Benford(FibonacciGenerator))</
=={{header|Liberty BASIC}}==
Using function from
http://rosettacode.org/wiki/Fibonacci_sequence#Liberty_BASIC
<syntaxhighlight lang="lb">
dim bin(9)
Line 2,061 ⟶ 2,488:
fiboI = a
end function
</syntaxhighlight>
{{out}}
Line 2,076 ⟶ 2,503:
9 0.045 0.046
</pre>
=={{header|Lua}}==
<
expected = {}
for i = 1, 9 do
Line 2,097 ⟶ 2,523:
for i = 1, 9 do
print(i, actual[i] / n, expected[i])
end</
{{out}}
<pre>digit actual expected
Line 2,109 ⟶ 2,535:
8 0.053 0.051152522447381
9 0.045 0.045757490560675</pre>
=={{header|Mathematica}} / {{header|Wolfram Language}}==
<
Table[{d, N@Count[fibdata, d]/Length@fibdata, Log10[1. + 1/d]}, {d, 1,
9}] // Grid</
{{out}}
<pre>1 0.301 0.30103
Line 2,124 ⟶ 2,549:
8 0.053 0.0511525
9 0.045 0.0457575</pre>
=={{header|MATLAB}}==
{{trans|Julia}}
<syntaxhighlight lang="MATLAB">
benfords_law();
function benfords_law
% Benford's Law
P = @(d) log10(1 + 1./d);
% Benford function
function counts = benford(numbers)
firstdigit = @(n) floor(mod(n / 10^floor(log10(n)), 10));
counts = zeros(1, 9);
for i = 1:length(numbers)
digit = firstdigit(numbers(i));
if digit ~= 0
counts(digit) = counts(digit) + 1;
end
end
counts = counts ./ sum(counts);
end
% Generate Fibonacci numbers
function fibNums = fibonacci(n)
fibNums = zeros(1, n);
a = 0;
b = 1;
for i = 1:n
c = b;
b = a + b;
a = c;
fibNums(i) = b;
end
end
% Sample
sample = fibonacci(1000);
% Observed and expected frequencies
observed = benford(sample) * 100;
expected = arrayfun(P, 1:9) * 100;
% Table
mytable = [1:9; observed; expected]';
% Plotting
bar(1:9, observed);
hold on;
plot(1:9, expected, 'LineWidth', 2);
hold off;
title("Benford's Law");
xlabel("First Digit");
ylabel("Frequency %");
legend("1000 Fibonacci Numbers", "P(d) = log10(1 + 1/d)");
xticks(1:9);
% Displaying the results
fprintf("Benford's Law\nFrequency of first digit\nin 1000 Fibonacci numbers\n");
disp(table(mytable(:,1),mytable(:,2),mytable(:,3),'VariableNames',{'digit', 'observed(%)', 'expected(%)'}))
end
</syntaxhighlight>
{{out}}
<pre>
Benford's Law
Frequency of first digit
in 1000 Fibonacci numbers
digit observed(%) expected(%)
_____ ___________ ___________
1 30 30.103
2 17.7 17.609
3 12.5 12.494
4 9.6 9.691
5 8 7.9181
6 6.7 6.6947
7 5.7 5.7992
8 5.3 5.1153
9 4.5 4.5757
</pre>
=={{header|NetRexx}}==
<
options replace format comments java crossref symbols nobinary
Line 2,167 ⟶ 2,673:
brenfordDeveation(fibList)
return
</syntaxhighlight>
{{out}}
<pre>
Line 2,182 ⟶ 2,688:
9: 4.500000% 4.575749% 0.075749%
</pre>
=={{header|Nim}}==
<
import strformat
Line 2,260 ⟶ 2,765:
let distrib = actualDistrib(fibSeq)
echo "Fibonacci numbers first digit distribution:\n"
distrib.display()</
{{out}}
Line 2,276 ⟶ 2,781:
8 0.0530 0.0512
9 0.0450 0.0458</pre>
=={{header|Oberon-2}}==
{{Works with|oo2c version 2}}
<
MODULE BenfordLaw;
IMPORT
Line 2,323 ⟶ 2,827:
END
END BenfordLaw.
</syntaxhighlight>
{{out}}
<pre>
Line 2,338 ⟶ 2,842:
9 0.045 0.046
</pre>
=={{header|OCaml}}==
For the Fibonacci sequence, we use the function from
https://rosettacode.org/wiki/Fibonacci_sequence#Arbitrary_Precision<br>
Note the remark about the compilation of the program there.
<
open Num
Line 2,374 ⟶ 2,877:
List.iter (Printf.printf "%f ") (List.map benfords_law xvalues) ;
Printf.printf "\n" ;;
</syntaxhighlight>
{{out}}
<pre>
Line 2,382 ⟶ 2,885:
0.301030 0.176091 0.124939 0.096910 0.079181 0.066947 0.057992 0.051153 0.045757
</pre>
=={{header|PARI/GP}}==
<
my(t=vector(9,n,sum(i=1,#v,v[i]==n)));
print("Digit\tActual\tExpected");
Line 2,392 ⟶ 2,894:
lucas(n)=fibonacci(n-1)+fibonacci(n+1);
dist(fibonacci)
dist(lucas)</
{{out}}
<pre>Digit Actual Expected
Line 2,415 ⟶ 2,917:
8 51 51
9 46 46</pre>
=={{header|Pascal}}==
<
{$IFDEF FPC}{$MODE Delphi}{$ELSE}{$APPTYPE CONSOLE}{$ENDIF}
uses
Line 2,469 ⟶ 2,970:
writeln(i:5,dgtCnt[i]:7,expectedCnt[i]:10,reldiff:10:5,' %');
end;
end.</
<pre>Digit Count Expected rel Diff
Line 2,481 ⟶ 2,982:
8 53 51 -3.92157 %
9 45 45 0.00000 %</pre>
=={{header|Perl}}==
<
use strict ;
use warnings ;
Line 2,507 ⟶ 3,007:
$result = sprintf ( "%.2f" , 100 * $expected[ $i - 1 ] ) ;
printf "%15s %%\n" , $result ;
}</
{{Out}}
<pre>
Line 2,521 ⟶ 3,021:
9 : 4.50 % 4.58 %
</pre>
=={{header|Phix}}==
{{trans|Go}}
<!--<syntaxhighlight lang="phix">(phixonline)-->
<span style="color: #008080;">with</span> <span style="color: #008080;">javascript_semantics</span>
<span style="color: #008080;">function</span> <span style="color: #000000;">benford</span><span style="color: #0000FF;">(</span><span style="color: #004080;">sequence</span> <span style="color: #000000;">s</span><span style="color: #0000FF;">,</span> <span style="color: #004080;">string</span> <span style="color: #000000;">title</span><span style="color: #0000FF;">)</span>
<span style="color: #004080;">sequence</span> <span style="color: #000000;">f</span> <span style="color: #0000FF;">=</span> <span style="color: #7060A8;">repeat</span><span style="color: #0000FF;">(</span><span style="color: #000000;">0</span><span style="color: #0000FF;">,</span><span style="color: #000000;">9</span><span style="color: #0000FF;">)</span>
<span style="color: #008080;">for</span> <span style="color: #000000;">i</span><span style="color: #0000FF;">=</span><span style="color: #000000;">1</span> <span style="color: #008080;">to</span> <span style="color: #7060A8;">length</span><span style="color: #0000FF;">(</span><span style="color: #000000;">s</span><span style="color: #0000FF;">)</span> <span style="color: #008080;">do</span>
<span style="color: #004080;">integer</span> <span style="color: #000000;">fdx</span> <span style="color: #0000FF;">=</span> <span style="color: #7060A8;">sprint</span><span style="color: #0000FF;">(</span><span style="color: #000000;">s</span><span style="color: #0000FF;">[</span><span style="color: #000000;">i</span><span style="color: #0000FF;">])[</span><span style="color: #000000;">1</span><span style="color: #0000FF;">]-</span><span style="color: #008000;">'0'</span>
<span style="color: #000000;">f</span><span style="color: #0000FF;">[</span><span style="color: #000000;">fdx</span><span style="color: #0000FF;">]</span> <span style="color: #0000FF;">+=</span> <span style="color: #000000;">1</span>
<span style="color: #008080;">end</span> <span style="color: #008080;">for</span>
<span style="color: #004080;">sequence</span> <span style="color: #000000;">res</span> <span style="color: #0000FF;">=</span> <span style="color: #0000FF;">{</span><span style="color: #000000;">title</span><span style="color: #0000FF;">,</span><span style="color: #008000;">"Digit Observed% Predicted%"</span><span style="color: #0000FF;">}</span>
<span style="color: #008080;">for</span> <span style="color: #000000;">i</span><span style="color: #0000FF;">=</span><span style="color: #000000;">1</span> <span style="color: #008080;">to</span> <span style="color: #7060A8;">length</span><span style="color: #0000FF;">(</span><span style="color: #000000;">f</span><span style="color: #0000FF;">)</span> <span style="color: #008080;">do</span>
<span style="color: #004080;">atom</span> <span style="color: #000000;">o</span> <span style="color: #0000FF;">=</span> <span style="color: #000000;">f</span><span style="color: #0000FF;">[</span><span style="color: #000000;">i</span><span style="color: #0000FF;">]/</span><span style="color: #7060A8;">length</span><span style="color: #0000FF;">(</span><span style="color: #000000;">s</span><span style="color: #0000FF;">)*</span><span style="color: #000000;">100</span><span style="color: #0000FF;">,</span>
<span style="color: #000000;">e</span> <span style="color: #0000FF;">=</span> <span style="color: #7060A8;">log10</span><span style="color: #0000FF;">(</span><span style="color: #000000;">1</span><span style="color: #0000FF;">+</span><span style="color: #000000;">1</span><span style="color: #0000FF;">/</span><span style="color: #000000;">i</span><span style="color: #0000FF;">)*</span><span style="color: #000000;">100</span>
<span style="color: #000000;">res</span> <span style="color: #0000FF;">=</span> <span style="color: #7060A8;">append</span><span style="color: #0000FF;">(</span><span style="color: #000000;">res</span><span style="color: #0000FF;">,</span><span style="color: #7060A8;">sprintf</span><span style="color: #0000FF;">(</span><span style="color: #008000;">" %d %9.3f %8.3f"</span><span style="color: #0000FF;">,</span> <span style="color: #0000FF;">{</span><span style="color: #000000;">i</span><span style="color: #0000FF;">,</span> <span style="color: #000000;">o</span><span style="color: #0000FF;">,</span> <span style="color: #000000;">e</span><span style="color: #0000FF;">}))</span>
<span style="color: #008080;">end</span> <span style="color: #008080;">for</span>
<span style="color: #008080;">return</span> <span style="color: #000000;">res</span>
<span style="color: #008080;">end</span> <span style="color: #008080;">function</span>
<span style="color: #008080;">function</span> <span style="color: #000000;">fib</span><span style="color: #0000FF;">(</span><span style="color: #004080;">integer</span> <span style="color: #000000;">lim</span><span style="color: #0000FF;">)</span>
<span style="color: #004080;">atom</span> <span style="color: #000000;">a</span><span style="color: #0000FF;">=</span><span style="color: #000000;">0</span><span style="color: #0000FF;">,</span> <span style="color: #000000;">b</span><span style="color: #0000FF;">=</span><span style="color: #000000;">1</span>
<span style="color: #004080;">sequence</span> <span style="color: #000000;">res</span> <span style="color: #0000FF;">=</span> <span style="color: #7060A8;">repeat</span><span style="color: #0000FF;">(</span><span style="color: #000000;">0</span><span style="color: #0000FF;">,</span><span style="color: #000000;">lim</span><span style="color: #0000FF;">)</span>
<span style="color: #008080;">for</span> <span style="color: #000000;">i</span><span style="color: #0000FF;">=</span><span style="color: #000000;">1</span> <span style="color: #008080;">to</span> <span style="color: #000000;">lim</span> <span style="color: #008080;">do</span>
<span style="color: #0000FF;">{</span><span style="color: #000000;">res</span><span style="color: #0000FF;">[</span><span style="color: #000000;">i</span><span style="color: #0000FF;">],</span> <span style="color: #000000;">a</span><span style="color: #0000FF;">,</span> <span style="color: #000000;">b</span><span style="color: #0000FF;">}</span> <span style="color: #0000FF;">=</span> <span style="color: #0000FF;">{</span><span style="color: #000000;">b</span><span style="color: #0000FF;">,</span> <span style="color: #000000;">b</span><span style="color: #0000FF;">,</span> <span style="color: #000000;">b</span><span style="color: #0000FF;">+</span><span style="color: #000000;">a</span><span style="color: #0000FF;">}</span>
<span style="color: #008080;">end</span> <span style="color: #008080;">for</span>
<span style="color: #008080;">return</span> <span style="color: #000000;">res</span>
<span style="color: #008080;">end</span> <span style="color: #008080;">function</span>
<span style="color: #004080;">sequence</span> <span style="color: #000000;">res</span> <span style="color: #0000FF;">=</span> <span style="color: #0000FF;">{</span><span style="color: #000000;">benford</span><span style="color: #0000FF;">(</span><span style="color: #000000;">fib</span><span style="color: #0000FF;">(</span><span style="color: #000000;">1000</span><span style="color: #0000FF;">),</span><span style="color: #008000;">"First 1000 Fibonacci numbers"</span><span style="color: #0000FF;">),</span>
<span style="color: #000000;">benford</span><span style="color: #0000FF;">(</span><span style="color: #7060A8;">get_primes</span><span style="color: #0000FF;">(-</span><span style="color: #000000;">10000</span><span style="color: #0000FF;">),</span><span style="color: #008000;">"First 10000 Prime numbers"</span><span style="color: #0000FF;">),</span>
<span style="color: #000000;">benford</span><span style="color: #0000FF;">(</span><span style="color: #7060A8;">sq_power</span><span style="color: #0000FF;">(</span><span style="color: #000000;">3</span><span style="color: #0000FF;">,</span><span style="color: #7060A8;">tagset</span><span style="color: #0000FF;">(</span><span style="color: #000000;">500</span><span style="color: #0000FF;">)),</span><span style="color: #008000;">"First 500 powers of three"</span><span style="color: #0000FF;">)}</span>
<span style="color: #7060A8;">papply</span><span style="color: #0000FF;">(</span><span style="color: #004600;">true</span><span style="color: #0000FF;">,</span><span style="color: #7060A8;">printf</span><span style="color: #0000FF;">,{</span><span style="color: #000000;">1</span><span style="color: #0000FF;">,{</span><span style="color: #008000;">"%-40s%-40s%-40s\n"</span><span style="color: #0000FF;">},</span><span style="color: #7060A8;">columnize</span><span style="color: #0000FF;">(</span><span style="color: #000000;">res</span><span style="color: #0000FF;">)})</span>
<!--</syntaxhighlight>-->
{{out}}
<pre>
First 1000 Fibonacci numbers First 10000 Prime numbers First 500 powers of three
Line 2,592 ⟶ 3,068:
9 4.500 4.576 9 10.060 4.576 9 4.600 4.576
</pre>
=={{header|Picat}}==
<syntaxhighlight lang="picat">go =>
N = 1000,
Fib = [fib(I) : I in 1..N],
check_benford(Fib),
nl.
% Check a list of numbers for Benford's law
check_benford(L) =>
Len = L.len,
println(len=Len),
Count = [F[1].to_integer() : Num in L, F=Num.to_string()].occurrences(),
P = new_map([I=D/Len : I=D in Count]),
println("Benford (percent):"),
foreach(D in 1..9)
B = benford(D)*100,
PI = P.get(D,0)*100,
Diff = abs(PI - B),
printf("%d: count=%5d observed: %0.2f%% Benford: %0.2f%% diff=%0.3f\n", D,Count.get(D,0),PI,B,Diff)
end,
nl.
benford(D) = log10(1+1/D).
% create an occurrences map of a list
occurrences(List) = Map =>
Map = new_map(),
foreach(E in List)
Map.put(E, cond(Map.has_key(E),Map.get(E)+1,1))
end.</syntaxhighlight>
{{out}}
<pre>Benford (percent):
1: count= 301 observed: 30.10% Benford: 30.10% diff=0.003
2: count= 177 observed: 17.70% Benford: 17.61% diff=0.091
3: count= 125 observed: 12.50% Benford: 12.49% diff=0.006
4: count= 96 observed: 9.60% Benford: 9.69% diff=0.091
5: count= 80 observed: 8.00% Benford: 7.92% diff=0.082
6: count= 67 observed: 6.70% Benford: 6.69% diff=0.005
7: count= 56 observed: 5.60% Benford: 5.80% diff=0.199
8: count= 53 observed: 5.30% Benford: 5.12% diff=0.185
9: count= 45 observed: 4.50% Benford: 4.58% diff=0.076
</pre>
'''Extra credit:'''
Using data from https://en.wikipedia.org/wiki/Land_use_statistics_by_country "Cultivated land (km^2)"
{{out}}
<pre>Benford (percent):
1: count= 72 observed: 33.64% Benford: 30.10% diff=3.542
2: count= 38 observed: 17.76% Benford: 17.61% diff=0.148
3: count= 29 observed: 13.55% Benford: 12.49% diff=1.058
4: count= 19 observed: 8.88% Benford: 9.69% diff=0.812
5: count= 20 observed: 9.35% Benford: 7.92% diff=1.428
6: count= 12 observed: 5.61% Benford: 6.69% diff=1.087
7: count= 9 observed: 4.21% Benford: 5.80% diff=1.594
8: count= 11 observed: 5.14% Benford: 5.12% diff=0.025
9: count= 4 observed: 1.87% Benford: 4.58% diff=2.707</pre>
=={{header|PicoLisp}}==
Picolisp does not support floating point math, but it can call libc math functions and convert the results to a fixed point number for e.g. natural logarithm.
<syntaxhighlight lang="picolisp">
(scl 4)
(load "@lib/misc.l")
Line 2,633 ⟶ 3,168:
(bye)
</syntaxhighlight>
{{Out}}
<pre>
Line 2,647 ⟶ 3,182:
9 0.045 0.046
</pre>
=={{header|PL/I}}==
<syntaxhighlight lang="pl/i">
(fofl, size, subrg):
Benford: procedure options(main); /* 20 October 2013 */
Line 2,689 ⟶ 3,223:
end tally;
end Benford;
</syntaxhighlight>
Results:
<pre>
Line 2,703 ⟶ 3,237:
9 0.04575749 0.04499817
</pre>
=={{header|PL/pgSQL}}==
<syntaxhighlight lang="sql">
WITH recursive
constant(val) AS
Line 2,740 ⟶ 3,273:
(select cast(corr(probability_theoretical,probability_real) as numeric(5,4)) correlation
from benford) c
</syntaxhighlight>
=={{header|PowerShell}}==
The sample file was not found. I selected another that contained the first two-thousand in the Fibonacci sequence, so there is a small amount of extra filtering.
<syntaxhighlight lang="powershell">
$url = "https://oeis.org/A000045/b000045.txt"
$file = "$env:TEMP\FibonacciNumbers.txt"
Line 2,761 ⟶ 3,293:
Remove-Item -Path $file -Force -ErrorAction SilentlyContinue
</syntaxhighlight>
{{Out}}
<pre>
Line 2,776 ⟶ 3,308:
9 45 0.045 0.04576
</pre>
=={{header|Prolog}}==
{{works with|SWI Prolog|6.2.6 by Jan Wielemaker, University of Amsterdam}}
Note: SWI Prolog implements arbitrary precision integer arithmetic through use of the GNU MP library
<
% Does the Fibonacci sequence follow Benford's law?
%~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
Line 2,823 ⟶ 3,354:
findall(B, (between(1,9,N), benford(N,B)), Benford),
findall(C, firstchar(C), Fc), freq(Fc, Freq),
writeHdr, maplist(writeData, Benford, Freq).</
{{out}}
<pre>?- go.
Line 2,836 ⟶ 3,367:
5.12% - 5.30%
4.58% - 4.50%</pre>
=={{header|PureBasic}}==
<
NewMap d1.i()
Dim fi.s(#MAX_N)
Line 2,877 ⟶ 3,407:
PrintN(~"\nPress Enter...")
Input()</
{{out}}
<pre>Dig. Cnt. Exp. Dif.
Line 2,892 ⟶ 3,422:
Press Enter...</pre>
=={{header|Python}}==
Works with Python 3.X & 2.7
<
from itertools import islice, count
from collections import Counter
Line 2,933 ⟶ 3,462:
# just to show that not all kind-of-random sets behave like that
show_dist("random", islice(heads(rand1000()), 10000))</
{{out}}
<pre>fibbed Benfords deviation
Line 2,967 ⟶ 3,496:
11.0% 5.1% 5.9%
10.9% 4.6% 6.3%</pre>
=={{header|R}}==
<syntaxhighlight lang="r">
pbenford <- function(d){
return(log10(1+(1/d)))
Line 3,001 ⟶ 3,529:
print(data)
</syntaxhighlight>
{{out}}
digit obs.frequency exp.frequency dev_percentage
Line 3,013 ⟶ 3,541:
8 0.053 0.05115 0.184748
9 0.045 0.04576 0.075749
=={{header|Racket}}==
<
(define (log10 n) (/ (log n) (log 10)))
Line 3,050 ⟶ 3,577:
;; 7: 5.8% 5.8%
;; 8: 5.1% 5.1%
;; 9: 4.6% 4.6%</
=={{header|Raku}}==
(formerly Perl 6)
{{Works with|rakudo|2016-10-24}}
<syntaxhighlight lang="raku"
sub show(%distribution) {
Line 3,068 ⟶ 3,594:
multi MAIN($file) { show benford $file.IO.lines }
multi MAIN() { show benford ( 1, 1, 2, *+* ... * )[^1000] }</
'''Output:''' First 1000 Fibonaccis
Line 3,099 ⟶ 3,625:
For the extra credit stuff, it was chosen to generate Fibonacci and factorials rather than find a
web─page with them listed, as each list is very easy to generate.
<
numeric digits length( e() ) - length(.) /*width of (e) for LN & LOG precision.*/
parse arg N .; if N=='' | N=="," then N= 1000 /*allow sample size to be specified. */
Line 3,131 ⟶ 3,657:
do f=1 for 9; say pad center(f,7) pad center(format(!.f/N,,length(N-2)),w1) #.f
end /*k*/
return</
{{out|output|text= when using the default (1000 numbers) for the input:}}
<pre>
Line 3,184 ⟶ 3,710:
9 0.0426 0.0457575
</pre>
=={{header|Ring}}==
<
# Project : Benford's law
Line 3,220 ⟶ 3,745:
if y = 1 return 1 ok
if y > 1 return fibonacci(y-1) + fibonacci(y-2) ok
</syntaxhighlight>
Output:
<pre>
Line 3,233 ⟶ 3,758:
8 5.300 5.115
9 4.500 4.576
</pre>
=={{header|RPL}}==
Here we use an interesting RPL feature that allows to pass an arithmetic expression or a function name as an argument. Thus, the code below can check the 'benfordance' of various series without needing to edit it. The resulting array being a matrix, some additional operations allow to compute the maximum difference at digit level between the tested function and Benford's law, as a distance indicator.
{{works with|Halcyon Calc|4.2.7}}
≪ → sn min max
≪ { 9 2 } 0 CON
min max FOR j
{ 0 1 } 1 j sn EVAL
IF DUP THEN
MANT FLOOR PUT
DUP2 GET 1 + PUT
ELSE
3 DROPN
END
NEXT
max min - 1 + /
1 9 FOR j
{ 0 2 } 1 j PUT
j INV 1 + LOG PUT
NEXT
DUP [ 1 0 ] *
OVER [ 0 1 ] * -
RNRM
≫
≫
'BENFD' STO
'→FIB' 1 100 BENFD
≪ LOG ≫ 1 10000 BENFD
Output:
<pre>
4: [[ 0.301 0.301029995664 ]
[ 0.177 0.176091259056 ]
[ 0.125 0.124938736608 ]
[ 0.096 9.69100130081E-02 ]
[ 0.08 7.91812460476E-02 ]
[ 6.7E-02 6.69467896306E-02 ]
[ 0.056 5.79919469777E-02 ]
[ 0.053 5.11525224474E-02 ]
[ 0.045 4.57574905607E-02 ]]
3: 1.8847477552E-02
2: [[ 9E-04 0.301029995664 ]
[ 9E-03 0.176091259056 ]
[ 0.09001 0.124938736608 ]
[ 0.90001 9.69100130081E-02 ]
[ 0.00001 7.91812460476E-02 ]
[ 0.00002 6.69467896306E-02 ]
[ 0.00001 5.79919469777E-02 ]
[ 0.00001 5.11525224474E-02 ]
[ 0.00002 4.57574905607E-02 ]]
1: 0.775161263392
</pre>
=={{header|Ruby}}==
{{trans|Python}}
<
def fib(n)
Line 3,273 ⟶ 3,850:
# just to show that not all kind-of-random sets behave like that
show_dist("random", random(10000))</
{{out}}
Line 3,310 ⟶ 3,887:
9: 10.8% 4.6% 6.2%
</pre>
=={{header|Run BASIC}}==
<
N = 1000
for i = 0 to N - 1
Line 3,343 ⟶ 3,919:
next i
end function
</syntaxhighlight>
<table border=1><TR bgcolor=wheat><TD>Digit<td>Actual<td>Expected</td><tr><tr align=right><td>1</td><td>30.100</td><td>30.103</td></tr><tr align=right><td>2</td><td>17.700</td><td>17.609</td></tr><tr align=right><td>3</td><td>12.500</td><td>12.494</td></tr><tr align=right><td>4</td><td> 9.500</td><td> 9.691</td></tr><tr align=right><td>5</td><td> 8.000</td><td> 7.918</td></tr><tr align=right><td>6</td><td> 6.700</td><td> 6.695</td></tr><tr align=right><td>7</td><td> 5.600</td><td> 5.799</td></tr><tr align=right><td>8</td><td> 5.300</td><td> 5.115</td></tr><tr align=right><td>9</td><td> 4.500</td><td> 4.576</td></tr></table>
=={{header|Rust}}==
{{works with|rustc|1.12 stable}}
Line 3,351 ⟶ 3,926:
This solution uses the ''num'' create for arbitrary-precision integers and the ''num_traits'' create for the ''zero'' and ''one'' implementations. It computes the Fibonacci numbers from scratch via the ''fib'' function.
<
extern crate num_traits;
extern crate num;
Line 3,403 ⟶ 3,978:
}
</syntaxhighlight>
{{out}}
<pre>
Line 3,417 ⟶ 3,992:
9: 0.045 v. 0.046
</pre>
=={{header|Scala}}==
<
val fibs : Stream[BigInt] = { def series(i:BigInt,j:BigInt):Stream[BigInt] = i #:: series(j, i+j); series(1,0).tail.tail }
Line 3,447 ⟶ 4,021:
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) )
}
}</
{{out}}
<pre>
Line 3,463 ⟶ 4,037:
9: 4.50% | 4.58% | 0.0757%
</pre>
=={{header|Scheme}}==
{{works with|Chez Scheme}}
<syntaxhighlight lang="scheme">; Compute the probability of leading digit d (an integer [1,9]) according to Benford's law.
(define benford-probability
(lambda (d)
(- (log (1+ d) 10) (log d 10))))
; Determine the distribution of the leading digit of the sequence provided by the given
; generator. Return as a vector of 10 elements -- the 0th element will always be 0.
(define leading-digit-distribution
(lambda (seqgen count)
(let ((digcounts (make-vector 10 0)))
(do ((index 0 (1+ index)))
((>= index count))
(let* ((value (seqgen))
(string (format "~a" value))
(leadchr (string-ref string 0))
(leaddig (- (char->integer leadchr) (char->integer #\0))))
(vector-set! digcounts leaddig (1+ (vector-ref digcounts leaddig)))))
(vector-map (lambda (digcnt) (/ digcnt count)) digcounts))))
; Create a Fibonacci sequence generator.
(define make-fibgen
(lambda ()
(let ((fn-2 0) (fn-1 1))
(lambda ()
(let ((fn fn-1))
(set! fn-1 (+ fn-2 fn-1))
(set! fn-2 fn)
fn)))))
; Create a sequence generator that returns elements of the given vector.
(define make-vecgen
(lambda (vector)
(let ((index 0))
(lambda ()
(set! index (1+ index))
(vector-ref vector (1- index))))))
; Read all the values in the given file into a list.
(define list-read-file
(lambda (filenm)
(call-with-input-file filenm
(lambda (ip)
(let accrue ((value (read ip)))
(if (eof-object? value)
'()
(cons value (accrue (read ip)))))))))
; Display a table of digit, Benford's law, sequence distribution, and difference.
(define display-table
(lambda (seqnam seqgen count)
(printf "~%~3@a ~11@a ~11@a ~11@a~%" "dig" "Benford's" seqnam "difference")
(let ((dist (leading-digit-distribution seqgen count)))
(do ((digit 1 (1+ digit)))
((> digit 9))
(let* ((fraction (vector-ref dist digit))
(benford (benford-probability digit))
(diff (- fraction benford)))
(printf "~3d ~11,5f ~11,5f ~11,5f~%" digit benford fraction diff))))))
; Emit tables of various sequence distributions.
(display-table "Fib/1000" (make-fibgen) 1000)
(display-table "Rnd/1T/1M" (lambda () (1+ (random 1000000000000))) 1000000)
(let ((craters (list->vector (list-read-file "moon_craters.lst"))))
(display-table "Craters/D" (make-vecgen craters) (vector-length craters)))</syntaxhighlight>
{{out}}
The first one thousand Fibonnaci numbers.
<pre>dig Benford's Fib/1000 difference
1 0.30103 0.30100 -0.00003
2 0.17609 0.17700 0.00091
3 0.12494 0.12500 0.00006
4 0.09691 0.09600 -0.00091
5 0.07918 0.08000 0.00082
6 0.06695 0.06700 0.00005
7 0.05799 0.05600 -0.00199
8 0.05115 0.05300 0.00185
9 0.04576 0.04500 -0.00076</pre>
{{out}}
One million pseudo-random numbers from one to one trillion.
<pre>dig Benford's Rnd/1T/1M difference
1 0.30103 0.11055 -0.19048
2 0.17609 0.11100 -0.06509
3 0.12494 0.11126 -0.01368
4 0.09691 0.11129 0.01438
5 0.07918 0.11141 0.03223
6 0.06695 0.11088 0.04394
7 0.05799 0.11120 0.05321
8 0.05115 0.11112 0.05997
9 0.04576 0.11128 0.06553</pre>
{{out}}
Diameters of 1,595 Lunar craters in meters.
<br />
From the Gazetteer of Planetary Nomenclature [https://planetarynames.wr.usgs.gov/],
referred to by the Wikipedia Planetary nomenclature page [https://en.wikipedia.org/wiki/Planetary_nomenclature].
<pre>dig Benford's Craters/D difference
1 0.30103 0.23950 -0.06153
2 0.17609 0.10658 -0.06951
3 0.12494 0.11912 -0.00582
4 0.09691 0.12414 0.02723
5 0.07918 0.11034 0.03116
6 0.06695 0.10596 0.03901
7 0.05799 0.06959 0.01160
8 0.05115 0.06646 0.01531
9 0.04576 0.05831 0.01255</pre>
=={{header|SETL}}==
<syntaxhighlight lang="setl">program benfords_law;
fibos := fibo_list(1000);
expected := [log(1 + 1/d)/log 10 : d in [1..9]];
actual := benford(fibos);
print('d Expected Actual');
loop for d in [1..9] do
print(d, ' ', fixed(expected(d), 7, 5), ' ', fixed(actual(d), 7, 5));
end loop;
proc benford(list);
dist := [];
loop for n in list do
dist(val(str n)(1)) +:= 1;
end loop;
return [d / #list : d in dist];
end proc;
proc fibo_list(n);
a := 1;
b := 1;
fibs := [];
loop while n>0 do
fibs with:= a;
[a, b] := [b, a+b];
n -:= 1;
end loop;
return fibs;
end proc;
end program;</syntaxhighlight>
{{out}}
<pre>d Expected Actual
1 0.30103 0.30100
2 0.17609 0.17700
3 0.12494 0.12500
4 0.09691 0.09600
5 0.07918 0.08000
6 0.06695 0.06700
7 0.05799 0.05600
8 0.05115 0.05300
9 0.04576 0.04500</pre>
=={{header|Sidef}}==
<
var fibonacci = 1000.of {|i| fib(i).digit(
for i in (1..9) {
var num = fibonacci.count_by {|j| j == i }
actuals.append(num / 1000)
Line 3,475 ⟶ 4,204:
"%17s%17s\n".printf("Observed","Expected")
for i in (1..9) {
"%d : %11s %%%15s %%\n".printf(
i, "%.2f".sprintf(100 * actuals[i - 1]),
"%.2f".sprintf(100 * expected[i - 1]),
)
}</
{{out}}
<pre>
Line 3,501 ⟶ 4,230:
The query is the same for any number sequence you care to put in the <tt>benford</tt> table.
<
create table benford (num integer);
Line 3,543 ⟶ 4,272:
-- Tidy up
drop table benford;</
{{out}}
Line 3,560 ⟶ 4,289:
9 rows selected.</pre>
=={{header|Stata}}==
<
set obs 1000
scalar phi=(1+sqrt(5))/2
Line 3,587 ⟶ 4,315:
8 | 53 51.15252245 |
9 | 45 45.75749056 |
+-----------------------------+</
Assuming the data are random, one can also do a goodness of fit [https://en.wikipedia.org/wiki/Pearson%27s_chi-squared_test chi-square test]:
<
chisq=sum((f-p):^2:/p)
chisq
Line 3,598 ⟶ 4,326:
chi2tail(8,chisq)
.9999942179
end</
The p-value is very close to 1, showing that the observed distribution is very close to the Benford law.
Line 3,604 ⟶ 4,332:
The fit is not as good with the sequence (2+sqrt(2))^n:
<
set obs 500
scalar s=2+sqrt(2)
Line 3,635 ⟶ 4,363:
chi2tail(8,chisq)
.0387287805
end</
Now the p-value is less than the usual 5% risk, and one would reject the hypothesis that the data follow the Benford law.
=={{header|Swift}}==
<
/* Reads from a file and returns the content as a String */
Line 3,743 ⟶ 4,470:
let ben = benford(digit: i+1)
print(String(format: "%d\t%.2f\t\t%.2f\t\t%.4f", i+1,ben*100,temp*100,ben-temp))
}</
{{out}}
<pre>$ ./Benford
Line 3,759 ⟶ 4,486:
8 5.12 5.30 -0.0018
9 4.58 4.50 0.0008</pre>
=={{header|Tcl}}==
<
# Count the leading digits (RE matches first digit in each number,
# even if negative)
Line 3,779 ⟶ 4,505:
[expr {log(1+1./$digit)/log(10)*100.0}]]
}
}</
Demonstrating with Fibonacci numbers:
<
for {set a 1;set b [set i 0]} {$i < $n} {incr i} {
lappend result [set b [expr {$a + [set a $b]}]]
Line 3,787 ⟶ 4,513:
return $result
}
benfordTest [fibs 1000]</
{{out}}
<pre>
Line 3,803 ⟶ 4,529:
</pre>
=={{header|True BASIC}}==
{{trans|Chipmunk Basic}}
<syntaxhighlight lang="qbasic">FUNCTION log10(n)
LET log10 = LOG(n)/LOG(10)
END FUNCTION
FUNCTION frequency(n)
LET frequency = (log10(n+1)-log10(n))
END FUNCTION
FUNCTION fibonacci(n)
!https://rosettacode.org/wiki/Fibonacci_sequence#True_BASIC
LET n1 = 0
LET n2 = 1
FOR k = 1 TO ABS(n)
LET sum = n1+n2
LET n1 = n2
LET n2 = sum
NEXT k
IF n < 0 THEN LET fibonacci = n1*((-1)^((-n)+1)) ELSE LET fibonacci = n1
END FUNCTION
CLEAR
LET n = 1000
DIM actual(0)
MAT REDIM actual(n)
FOR nr = 1 TO n
LET num$ = STR$(fibonacci(nr))
LET j = VAL((num$)[1:1])
LET actual(j) = actual(j)+1
NEXT nr
PRINT "First 1000 Fibonacci numbers"
PRINT "Digit Actual freq Expected freq"
FOR i = 1 TO 9
LET freq = frequency(i)*100
PRINT USING "###": i;
PRINT USING " ##.###": actual(i)/10;
PRINT USING " ##.###": freq
NEXT i
END</syntaxhighlight>
=={{header|VBA (Visual Basic for Application)}}==
<syntaxhighlight lang="vb">
Sub BenfordLaw()
Dim BenResult(1 To 9) As Long
BENref = "30,1%|17,6%|12,5%|9,7%|7,9%|6,7%|5,8%|5,1%|4,6%"
For Each c In Selection.Cells
If InStr(1, "-0123456789", Left(c, 1)) > 0 Then
For i = 1 To 9
If CInt(Left(Abs(c), 1)) = i Then BenResult(i) = BenResult(i) + 1: Exit For
Next
End If
Next
Total= Application.Sum(BenResult)
biggest= Len(CStr(BenResult(1)))
txt = "# | Values | Real | Expected " & vbCrLf
For i = 1 To 9
If BenResult(i) > 0 Then
txt = txt & "#" & i & " | " & vbTab
txt = txt & String((biggest - Len(CStr(BenResult(i)))) * 2, " ") & Format(BenResult(i), "0") & " | " & vbTab
txt = txt & String((Len(CStr(Format(BenResult(1) / Total, "##0.0%"))) - Len(CStr(Format(BenResult(i) / Total, "##0.0%")))) * 2, " ") & Format(BenResult(i) / Total, "##0.0%") & " | " & vbTab
txt = txt & Format(Split(BENref, "|")(i - 1), " ##0.0%") & vbCrLf
End If
Next
MsgBox txt, vbOKOnly, "Finish"
End Sub
}
</syntaxhighlight>
=={{header|Visual FoxPro}}==
<
#DEFINE CTAB CHR(9)
#DEFINE COMMA ","
Line 3,845 ⟶ 4,646:
RETURN LOG10(1 + 1/d)
ENDFUNC
</syntaxhighlight>
{{out}}
<pre>
Line 3,860 ⟶ 4,661:
Correlation Coefficient: 0.999908
</pre>
=={{header|V (Vlang)}}==
{{trans|Go}}
{{libheader|Wren-fmt}}
<syntaxhighlight lang="v (vlang)">import math
fn fib1000() []f64 {
mut a, mut b, mut r := 0.0, 1.0, []f64{len:1000}
for i in 0..r.len {
r[i], a, b = b, b, b+a
}
return r
}
fn main() {
show(fib1000(), "First 1000 Fibonacci numbers")
}
fn show(c []f64, title string) {
mut f := [9]int{}
for v in c {
f["$v"[0]-'1'[0]]++
}
println(title)
println("Digit Observed Predicted")
for i, n in f {
println(" ${i+1} ${f64(n)/f64(c.len):9.3f} ${math.log10(1+1/f64(i+1)):8.3f}")
}
}</syntaxhighlight>
{{out}}
<pre>
First 1000 Fibonacci numbers:
Digit Observed Predicted
1 0.301 0.301
2 0.177 0.176
3 0.125 0.125
4 0.096 0.097
5 0.080 0.079
6 0.067 0.067
7 0.056 0.058
8 0.053 0.051
9 0.045 0.046
</pre>
Line 3,865 ⟶ 4,709:
{{trans|Go}}
{{libheader|Wren-fmt}}
<
var fib1000 = Fn.new {
Line 3,902 ⟶ 4,746:
}
show.call(fib1000.call(), "First 1000 Fibonacci numbers:")</
{{out}}
Line 3,918 ⟶ 4,762:
9 0.045 0.046
</pre>
=={{header|XPL0}}==
The program is run like this: benford < fibdig.txt
Luckly, there's no need to show how the first digits of the Fibonacci sequence were obtained.
<syntaxhighlight lang "XPL0">int D, N, Counts(10);
real A, E;
[for D:= 0 to 9 do Counts(D):= 0;
for N:= 1 to 1000 do
[D:= ChIn(1) & $0F; \ASCII to binary
Counts(D):= Counts(D)+1;
];
Text(0, "Digit Actual Expected Difference^m^j");
for D:= 1 to 9 do
[IntOut(0, D); ChOut(0, ^ );
A:= float(Counts(D))/1000.;
RlOut(0, A);
E:= Log(1. + 1./float(D));
RlOut(0, E);
RlOut(0, E-A);
CrLf(0);
];
]</syntaxhighlight>
{{out}}
<pre>
Digit Actual Expected Difference
1 0.30100 0.30103 0.00003
2 0.17700 0.17609 -0.00091
3 0.12500 0.12494 -0.00006
4 0.09600 0.09691 0.00091
5 0.08000 0.07918 -0.00082
6 0.06700 0.06695 -0.00005
7 0.05600 0.05799 0.00199
8 0.05300 0.05115 -0.00185
9 0.04500 0.04576 0.00076
</pre>
=={{header|Yabasic}}==
Using function from https://www.rosettacode.org/wiki/Fibonacci_sequence#Yabasic
<syntaxhighlight lang="vb">n = 1000
dim actual(n)
for nr = 1 to n
num$ = str$(fibonacci(nr))
j = val(left$(num$,1))
actual(j) = actual(j) + 1
next
print "First 1000 Fibonacci numbers"
print "Digit ", "Actual freq ", "Expected freq"
for i = 1 to 9
freq = frequency(i)*100
print i using("###"), " ", (actual(i)/10) using("##.###"), " ", freq using("##.###")
next
end
sub frequency(n)
return (log10(n+1) - log10(n))
end sub
sub log10(n)
return log(n) / log(10)
end sub
sub fibonacci(n)
local n1, n2, k, sum
n1 = 0
n2 = 1
for k = 1 to abs(n)
sum = n1 + n2
n1 = n2
n2 = sum
next k
if n < 0 then
return n1 * ((-1) ^ ((-n) + 1))
else
return n1
end if
end sub</syntaxhighlight>
{{out}}
<pre>First 1000 Fibonacci numbers
Digit Actual freq Expected freq
1 30.100 30.103
2 17.700 17.609
3 12.500 12.494
4 9.600 9.691
5 8.000 7.918
6 6.700 6.695
7 5.600 5.799
8 5.300 5.115
9 4.500 4.576</pre>
=={{header|zkl}}==
{{trans|Go}}
<
(0).pump(1000,List,fcn(ab){ab.append(ab.sum(0.0)).pop(0)}.fp(L(1,1))),
"First 1000 Fibonacci numbers");
Line 3,935 ⟶ 4,871:
(1.0+1.0/i).log10()))
}
}</
{{trans|CoffeeScript}}
<
fcn fibgen(a,b) { return(a,self.fcn.fp(b,a+b)) } //-->L(fib,fcn)
Line 3,957 ⟶ 4,893:
println("Leading digital distribution of the first 1,000 Fibonacci numbers");
println("Digit\tActual\tExpected");
foreach i in ([1..9]){ println("%d\t%.3f\t%.3f".fmt(i,actual[i], expected[i-1])); }</
{{out}}
<pre>
Line 3,972 ⟶ 4,908:
9 0.045 0.046
</pre>
=={{header|ZX Spectrum Basic}}==
{{trans|Liberty BASIC}}
<
20 DIM b(9)
30 LET n=100
Line 3,998 ⟶ 4,933:
1060 NEXT j
1070 RETURN
</syntaxhighlight>
The results obtained are adjusted fairly well, except for the number 8. This occurs with Sinclair BASIC, Sam BASIC and SpecBAS fits.
|