56 lines
1.6 KiB
Text
56 lines
1.6 KiB
Text
import "random" for Random
|
|
import "./fmt" for Fmt
|
|
import "./math" for Nums
|
|
|
|
var rgen = Random.new()
|
|
|
|
// Box-Muller method from Wikipedia
|
|
var normal = Fn.new { |mu, sigma|
|
|
var u1 = rgen.float()
|
|
var u2 = rgen.float()
|
|
var mag = sigma * (-2 * u1.log).sqrt
|
|
var z0 = mag * (2 * Num.pi * u2).cos + mu
|
|
var z1 = mag * (2 * Num.pi * u2).sin + mu
|
|
return [z0, z1]
|
|
}
|
|
|
|
var N = 100000
|
|
var NUM_BINS = 12
|
|
var HIST_CHAR = "■"
|
|
var HIST_CHAR_SIZE = 250
|
|
var bins = List.filled(NUM_BINS, 0)
|
|
var binSize = 0.1
|
|
var samples = List.filled(N, 0)
|
|
var mu = 0.5
|
|
var sigma = 0.25
|
|
for (i in 0...N/2) {
|
|
var rns = normal.call(mu, sigma)
|
|
for (j in 0..1) {
|
|
var rn = rns[j]
|
|
var bn
|
|
if (rn < 0) {
|
|
bn = 0
|
|
} else if (rn >= 1) {
|
|
bn = 11
|
|
} else {
|
|
bn = (rn/binSize).floor + 1
|
|
}
|
|
bins[bn] = bins[bn] + 1
|
|
samples[i*2 + j] = rn
|
|
}
|
|
}
|
|
|
|
Fmt.print("Normal distribution with mean $0.2f and S/D $0.2f for $,d samples:\n", mu, sigma, N)
|
|
System.print(" Range Number of samples within that range")
|
|
for (i in 0...NUM_BINS) {
|
|
var hist = HIST_CHAR * (bins[i] / HIST_CHAR_SIZE).round
|
|
if (i == 0) {
|
|
Fmt.print(" -∞ ..< 0.00 $s $,d", hist, bins[0])
|
|
} else if (i < NUM_BINS - 1) {
|
|
Fmt.print("$4.2f ..< $4.2f $s $,d", binSize * (i-1), binSize * i, hist, bins[i])
|
|
} else {
|
|
Fmt.print("1.00 ... +∞ $s $,d", hist, bins[NUM_BINS - 1])
|
|
}
|
|
}
|
|
Fmt.print("\nActual mean for these samples : $0.5f", Nums.mean(samples))
|
|
Fmt.print("Actual S/D for these samples : $0.5f", Nums.stdDev(samples))
|