[Date Prev][Date Next][Thread Prev][Thread Next][Date Index][Thread Index]

normal distribution



Hi Rich,
try this out, (it took me a while to find the Mathematica stuff).
 
 
;;;-------------------------------------------------------------------------
;;; Gaussian Random Number Generator
;;;
;;; Returns a random number with a "normal" or "Gaussian" Distribution
;;;
;;; Adapted from Mathematica
;;;
;;; Notes:
;;;     1.) Mathematica code has been tested, lisp code has not
;;;     2.) You could speed up the lisp code w/Erran Gat's fpc package
;;;-------------------------------------------------------------------------
 
 
#|
;; Here is Mathematica code to see:
 
;; Returns a random # from ~ -4 to ~ 4 (as default), with mean mu, and sigma =
1.
GaussianRandom[mu_:0, sigma_:1]:= mu + sigma Sqrt[-2 Log[Random[]]] Cos[2Pi
Random[]]//N
 
 
;; Creates a table of 2000 gaussian random #'s
foobar = Sort[Table[GaussianRandom[0,1],{i,2000}]];
 
 
;; Code to display distribution of random #'s
Buckets[aList_] := Module[{
    histCount = Length[aList]-1,
    minItem=Min[aList],
    maxItem=Max[aList],
    },
    Map[Round[((#-minItem)/(maxItem-minItem)) *histCount] &, aList]]
CountBins[aList_] := Module[{
       listLength = Length[aList],
       binCount = Table[0,{i, Length[aList]}]
       },
    Do[binCount[[aList[[i]]+1]]++, {i, 1, listLength}];
       binCount
      ]
Histogram[aList_] := CountBins[Buckets[aList]]
ListPlot[Histogram[foobar],PlotJoined->True]
|#
 
;;Lisp equivalent
 
 
;;Returns a random number between 0 and 1
(defun unit-random (&optional (precision 1000000000))
  (coerce (/ (random precision) precision) 'float))
 
;; Returns a random # with normal distribution from ~ -4 to ~ 4 (as default),
with mean mu, and sigma = 1.
(defun gaussian-random (&optional (mu 0) (sigma 1))
  (+ mu (* sigma (sqrt (* -2.0 (log (unit-random)) (cos (* 2 3.1415926535
(unit-random))))))))