Statistics/Normal distribution: Difference between revisions

 

=={{header|C}}==

* RosettaCode example: Statistics/Normal distribution in C

*

}

return EXIT_FAILURE;

{{out}}

<pre>mean = 0.000477941, stddev = 0.999945

=={{header|C sharp|C#}}==

{{libheader|Math.Net}}

using MathNet.Numerics.Distributions;

using MathNet.Numerics.Statistics;

RunNormal(10000);

}

{{out}}

<pre>Sample size: 100

=={{header|C++}}==

showing features of C++11 here

#include <map>

#include <string>

}

return 0 ;

Output:

<pre>The mean of the distribution is 1 , the standard deviation 1 !

=={{header|D}}==

This uses the Box-Muller method as in the Go entry, and the module from the Statistics/Basic. A ziggurat-based normal generator for the Phobos standard library is in the works.

statistics_basic;

 

writefln("Mean: %8.6f, SD: %8.6f\n", data.meanStdDev[]);

data.map!q{ max(0.0, min(0.9999, a / 3 + 0.5)) }.showHistogram01;

{{out}}

<pre>Mean: 0.000528, SD: 0.502245

0.8: ******

0.9: *</pre>

 

=={{header|Elixir}}==

def normal_distribution(n, w\\5) do

{sum, sum2, hist} = generate(n, w)

Enum.each([100,1000,10000], fn n ->

Statistics.normal_distribution(n)

 

{{out}}

=={{header|Factor}}==

{{works with|Factor|0.99 2020-01-23}}

math.statistics random sequences sorting ;

 

[ [ CHAR: * ] "" replicate-as ] dip ! make the bar

"% 5.2f: %s %d\n" printf ! print a line of the histogram

{{out}}

<pre style="font-size:80%; height: 120ex; overflow: scroll">

{{works with|Fortran|95 and later}}

Using the Marsaglia polar method

implicit none

 

end do

{{out}}

<pre>

 

=={{header|FreeBASIC}}==

 

Const pi As Double = 3.141592653589793

Print

Print "Press any key to quit"

Sample output:

{{out}}

=={{header|Go}}==

Box-Muller method shown here. Go has a normally distributed random function in the standard library, as shown in the Go [[Random numbers]] solution. It uses the ziggurat method.

 

import (

fmt.Println(strings.Repeat("*", p/scale))

}

Output:

<pre>

 

=={{header|Haskell}}==

import Data.Maybe (fromMaybe)

import Text.Printf (printf)

main = do

runTest 1000

 

{{out}}

=={{header|J}}==

'''Solution'''

rnorm01=. (2 o. 2p1 * runif01) * [: %: _2 * ^.@runif01 NB. random normal number generator (Box-Muller)

 

mean=: +/ % # NB. mean

stddev=: (<:@# %~ +/)&.:*:@(- mean) NB. standard deviation

'''Example Usage'''

(mean , stddev) DataSet

0.000781667 1.00154

require 'plot'

 

=={{header|Java}}==

{{trans|D}}

{{works with|Java|8}}

import static java.util.Arrays.stream;

import java.util.Locale;

.toArray());

}

<pre>Mean: -0.001870, SD: 0.500539

0.0: **

'''Preliminaries'''

# Input: a number without an exponent

# Output: a string holding the reformatted number so that there are at least `left` characters

def sample_mean_and_variance:

.mean = (.sum/.n)

'''The Task'''

def parameters: {

N: 100000,

histogram ;

 

{{out}}

<pre>

=={{header|Julia}}==

Julia has the builtin package "Distributions" to generate random numbers from a standard distribution (Normal, Chisq etc.).

 

data = rand(Normal(0, 1), 1000)

@printf("range = (%2.2f, %2.2f\n)", minimum(data), maximum(data))

h = plot(x=data, Geom.histogram)

 

{{out}}

=={{header|Kotlin}}==

{{trans|FreeBASIC}}

 

val rand = java.util.Random()

val sampleSizes = intArrayOf(100, 1_000, 10_000, 100_000)

for (sampleSize in sampleSizes) normalStats(sampleSize)

 

{{out}}

 

=={{header|Lasso}}==

if(#a->size) => {

local(mean = (with n in #a sum #n) / #a->size)

histogram(#n)

'\r\r'

 

{{out}}

=={{header|Liberty BASIC}}==

Uses LB Statistics/Basic

 

end

v =rnd( 1)

normalDist =( -2 *log( u))^0.5 *cos( 2 *3.14159265 *v)

100000 data terms used.

Largest term was 4.12950792 & smallest was -4.37934139

=={{header|Lua}}==

Lua provides math.random() to generate uniformly distributed random numbers. The function gaussian() shown here uses math.random() to generate normally distributed random numbers with given mean and variance.

return math.sqrt(-2 * variance * math.log(math.random())) *

math.cos(2 * math.pi * math.random()) + mean

print("Mean:", mean(t) .. ", expected " .. average)

print("StdDev:", std(t) .. ", expected " .. math.sqrt(variance) .. "\n")

{{out}}

<pre>Mean: 50.008328894275, expected 50

=={{header|Maple}}==

Maple has a built-in for sampling directly from [http://www.maplesoft.com/support/help/Maple/view.aspx?path=Statistics/Distributions/Normal Normal] distributions:

n := 100000:

X := Sample( Normal(0,1), n );

Mean( X );

StandardDeviation( X );

 

=={{header|Mathematica}}/{{header|Wolfram Language}}==

SampleNormal[n_] := (Print[#//Length, " numbers, Mean : ", #//Mean, ", StandardDeviation : ", #//StandardDeviation];

Histogram[#, BarOrigin -> Left,Axes -> False])& [(Table[(-2*Log[x])^0.5*Cos[2*Pi*x], {n} ]]

SampleNormal[ 10000 ]

->10000 numbers, Mean : -0.0122308, StandardDeviation : 1.00646

[[File:Mma_NormalDistribution.png]]

 

=={{header|MATLAB}} / {{header|Octave}}==

x = randn(N,1);

mean(x)

std(x)

[nn,xx] = hist(x,100);

 

=={{header|Nim}}==

Here is a way to generate random values following normal distribution from random values following uniform distribution. It uses the Basic form of the Box-Muller transform.

 

 

proc drawHistogram(ns: seq[float]) =

 

echo &"μ = {z.mean:.12f} σ = {z.standardDeviation:.12f}"

 

{{out}}

=={{header|PARI/GP}}==

{{works with|PARI/GP|2.4.3 and above}}

\\ Could easily be extended with a second normal at very little cost.

};

for(i=1,#h,for(j=1,h[i]\sz,print1("#"));print());

};

 

For versions before 2.4.3, define

my(pr=32*ceil(default(realprecision)*log(10)/log(4294967296))); \\ Current precision

random(2^pr)*1.>>pr

and use <code>rreal()</code> in place of <code>random(1.)</code>.

 

A PARI implementation:

rnormal(long prec)

{

ret = gerepileupto(ltop, ret);

return ret;

Use <code>mpsincos</code> and caching to generate two values at nearly the same cost.

 

 

From [http://www.freepascal.org/docs-html/rtl/math/randg.html Free Pascal Docs unit math]

{$IFDEF FPC}

{$MOde objFPC}

mySol := getSol(@rnorm,cMean,cStdDiv,1000000);

Histo(mySol,cHistCnt,cColLen);

{{out}}<pre>

function randg

=={{header|Perl}}==

{{trans|Raku}}

use List::Util qw(sum reduce min max);

 

print "$i\t$t1$t2\n";

}

{{out}}

<pre style="height:35ex">32 ⎸

=={{header|Phix}}==

{{trans|Liberty_BASIC}}

<span style="color: #008080;">with</span> <span style="color: #008080;">javascript_semantics</span>

<span style="color: #008080;">procedure</span> <span style="color: #000000;">sample</span><span style="color: #0000FF;">(</span><span style="color: #004080;">integer</span> <span style="color: #000000;">n</span><span style="color: #0000FF;">)</span>

<span style="color: #000000;">sample</span><span style="color: #0000FF;">(</span><span style="color: #000000;">100000</span><span style="color: #0000FF;">)</span>

{{Out}}

<pre>

</pre>

{{trans|Lua}}

<span style="color: #008080;">with</span> <span style="color: #008080;">javascript_semantics</span>

<span style="color: #008080;">function</span> <span style="color: #000000;">gaussian</span><span style="color: #0000FF;">(</span><span style="color: #004080;">atom</span> <span style="color: #000000;">mean</span><span style="color: #0000FF;">,</span> <span style="color: #004080;">atom</span> <span style="color: #000000;">variance</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;">"StdDev: %g, expected %g\n"</span><span style="color: #0000FF;">,{</span><span style="color: #000000;">std</span><span style="color: #0000FF;">(</span><span style="color: #000000;">t</span><span style="color: #0000FF;">),</span><span style="color: #7060A8;">sqrt</span><span style="color: #0000FF;">(</span><span style="color: #000000;">variance</span><span style="color: #0000FF;">)})</span>

<span style="color: #000000;">showHistogram</span><span style="color: #0000FF;">(</span><span style="color: #000000;">t</span><span style="color: #0000FF;">)</span>

{{Out}}

<pre>

 

=={{header|PureBasic}}==

ProcedureReturn Random(resolution) / resolution

EndProcedure

Print(#CRLF$ + #CRLF$ + "Press ENTER to exit"): Input()

CloseConsole()

Sample output:

<pre>100000 data terms used.

=={{header|Python}}==

This uses the external [http://matplotlib.org/ matplotlib] package as well as the built-in standardlib function [http://docs.python.org/2/library/random.html?highlight=gauss#random.gauss random.gauss].

import matplotlib.pyplot as plt

import random

% (mn, sd, max(data), min(data), size))

 

 

{{out}}

 

Generate normal random numbers from uniform random numbers, using the [https://en.wikipedia.org/wiki/Box%E2%80%93Muller_transform Box-Muller transform]. Both x and y are normally distributed, and independent.

u <- sqrt(-2*log(runif(n)))

v <- 2*pi*runif(n)

x <- u*cos(v)

y <- v*sin(v)

 

 

R can generate random normal distributed numbers using the [https://stat.ethz.ch/R-manual/R-devel/library/stats/html/Normal.html rnorm] function:

x <- rnorm(n, mean=0, sd=1)

mean(x)

sd(x)

 

=={{header|Racket}}==

compute statistics and plot a histogram.

[[File:histogram-racket.png|thumb|right]]

#lang racket

(require math (planet williams/science/histogram-with-graphics))

(for ([x data]) (histogram-increment! h x))

(histogram-plot h "Normal distribution μ=50 σ=4")

 

The other part of the task was to produce normal distributed numbers from a unit distribution.

version of [http://planet.plt-scheme.org/package-source/schematics/random.plt/1/0/random.ss code]

originally written by Sebastian Egner.

#lang racket

(require math)

(set! next (* scale v2))

(+ μ (* σ scale v1))])))))))

 

=={{header|Raku}}==

(formerly Perl 6)

{{works with|Rakudo|2018.03}}

my $r = sqrt -2 * log rand;

my $Θ = τ * rand;

say $i, "\t", '█' x $full, @subbar[$part];

}

{{out}}

<pre>"m" => 50.006107405837142e0

The REXX language doesn't have any "higher math" BIF functions like &nbsp; SIN, COS, LN, LOG, SQRT, EXP, POW, etc,

<br>so we hoi polloi programmers have to roll our own.

numeric digits 20 /*use twenty decimal digs for accuracy.*/

parse arg n seed . /*obtain optional arguments from the CL*/

sqrt: procedure; parse arg x; if x=0 then return 0; d= digits(); m.= 9; numeric digits; numeric form; h= d+6

parse value format(x,2,1,,0) 'E0' with g 'E' _ .; g=g*.5'e'_%2; do j=0 while h>9; m.j=h; h=h%2+1; end /*j*/

This REXX program makes use of &nbsp; '''scrsize''' &nbsp; REXX program (or BIF) which is used to determine the screen size of the terminal (console); &nbsp; this is to aid in maximizing the width of the horizontal histogram.

 

{{works with|Ruby|2.7}}

'''The Implementation'''

# Uses the Marsaglia polar method.

 

end

end

'''The Task'''

{{libheader|enumerable-statistics}}

 

def show_stats_and_histogram(data, bins)

show_stats_and_histogram(100.times.map { gen_normal.rand }, 40)

show_stats_and_histogram(10000.times.map { gen_normal.rand }, 60)

{{out}}

<pre style="height: 200ex; overflow: scroll">

 

=={{header|Run BASIC}}==

s = 100000

h$ = "============================================================="

print using("##",b);" ";using("#####",bins(b));" ";left$(h$,(bins(b) / mb) * 90)

next b

{{out}}

<pre>

{{libheader|rand}}

{{libheader|rand_distr}}

use math::{histogram::Histogram, traits::ToIterator};

use rand;

print_histogram(&data, 80, 40, '-');

}

{{out}}

<pre style="height: 120ex; overflow: scroll">

 

=={{header|SAS}}==

n=100000;

twopi=2*constant('pi');

y 0.000023995 1.0019996

z 0.0012857 1.0056536

 

=={{header|Sidef}}==

{{trans|Raku}}

 

func normdist (m, σ) {

var part = (8 * (x - full))

say (i, "\t", '█' * full, subbar[part])

{{out}}

<pre>

Pairs of normal numbers are generated from pairs of uniform numbers using the [https://en.wikipedia.org/wiki/Box%E2%80%93Muller_transform Box-Muller method]. A normal density is added to the histogram for comparison. See '''[http://www.stata.com/help.cgi?histogram histogram]''' in Stata help. A [https://en.wikipedia.org/wiki/Q%E2%80%93Q_plot Q-Q plot] is also drawn.

 

set obs 100000

gen u=runiform()

hist y, normal

hist z, normal

 

=={{header|Tcl}}==

# Uses the Box-Muller transform to compute a pair of normal random numbers

proc tcl::mathfunc::nrand {mean stddev} {

stats 10000

puts ""

Sample output:

<pre>

 

=={{header|VBA}}==

Dim s() As Variant, bins(71) As Single

ReDim s(20000)

Debug.Print

Next i

<pre>sample size 20000 mean-5,26306310478751E-03 standard deviation 1,00355037427319

-3,60--3,50

{{libheader|Wren-fmt}}

{{libheader|Wren-math}}

 

var rgen = Random.new()

}

Fmt.print("\nActual mean for these samples : $0.5f", Nums.mean(samples))

 

{{out}}

=={{header|zkl}}==

{{trans|Go}}

const PI2=(0.0).pi*2;;

rnd:=(0.0).random.fp(1); // random number in [0,1), using partial application

b:=(v + SIG)*BINS/SIG/2;

if(0<=b<BINS) h[b]+=1;

Partial application: rnd() --> (0.0).random(1). Basically, the fp method fixes the call parameters, which are then used when the partial thing is run.

println("Samples: %,d".fmt(N));

println("Mean: ", m:=sum/N);

println("Stddev: ", (sumSq/N - m*m).sqrt());

{{out}}

<pre>