1.1 --- /dev/null	Thu Jan 01 00:00:00 1970 +0000
     1.2 +++ b/src/clojure/contrib/test_contrib/probabilities/examples_monte_carlo.clj	Sat Aug 21 06:25:44 2010 -0400
     1.3 @@ -0,0 +1,73 @@
     1.4 +;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;
     1.5 +;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;
     1.6 +;;
     1.7 +;; Monte-Carlo application examples
     1.8 +;;
     1.9 +;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;
    1.10 +;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;
    1.11 +
    1.12 +(ns
    1.13 +  #^{:author "Konrad Hinsen"
    1.14 +     :skip-wiki true
    1.15 +     :doc "Examples for monte carlo methods"}
    1.16 +  clojure.contrib.probabilities.random.examples-monte-carlo
    1.17 +  (:require [clojure.contrib.generic.collection :as gc])
    1.18 +  (:use [clojure.contrib.probabilities.random-numbers
    1.19 +	 :only (lcg rand-stream)])
    1.20 +  (:use [clojure.contrib.probabilities.finite-distributions
    1.21 +	 :only (uniform)])
    1.22 +  (:use [clojure.contrib.probabilities.monte-carlo
    1.23 +	 :only (random-stream discrete interval normal lognormal exponential
    1.24 +		n-sphere
    1.25 +		sample sample-sum sample-mean sample-mean-variance)]
    1.26 +	:reload)
    1.27 +  (:use [clojure.contrib.monads
    1.28 +	:only (domonad state-m)]))
    1.29 +
    1.30 +; Create a linear congruential generator
    1.31 +(def urng (lcg 259200 7141 54773 1))
    1.32 +
    1.33 +;; Use Clojure's built-in random number generator
    1.34 +;(def urng rand-stream)
    1.35 +
    1.36 +; Sample transformed distributions
    1.37 +(defn sample-distribution
    1.38 +  [n rt]
    1.39 +  (take n (gc/seq (random-stream rt urng))))
    1.40 +
    1.41 +; Interval [-2, 2)
    1.42 +(sample-distribution 10 (interval -2 2))
    1.43 +; Compare with a direct transformation
    1.44 +(= (sample-distribution 10 (interval -2 2))
    1.45 +   (map (fn [x] (- (* 4 x) 2)) (take 10 (gc/seq urng))))
    1.46 +
    1.47 +; Normal distribution
    1.48 +(sample-distribution 10 (normal 0 1))
    1.49 +
    1.50 +; Log-Normal distribution
    1.51 +(sample-distribution 10 (lognormal 0 1))
    1.52 +
    1.53 +; Exponential distribution
    1.54 +(sample-distribution 10 (exponential 1))
    1.55 +
    1.56 +; n-sphere distribution
    1.57 +(sample-distribution 10 (n-sphere 2 1))
    1.58 +
    1.59 +; Discrete distribution
    1.60 +(sample-distribution 10 (discrete (uniform (range 1 7))))
    1.61 +
    1.62 +; Compose distributions in the state monad
    1.63 +(def sum-two-dists
    1.64 +  (domonad state-m
    1.65 +    [r1 (interval -2 2)
    1.66 +     r2 (normal 0 1)]
    1.67 +    (+ r1 r2)))
    1.68 +
    1.69 +(sample-distribution 10 sum-two-dists)
    1.70 +
    1.71 +; Distribution transformations
    1.72 +(sample-distribution  5 (sample 2 (interval -2 2)))
    1.73 +(sample-distribution 10 (sample-sum 10 (interval -2 2)))
    1.74 +(sample-distribution 10 (sample-mean 10 (interval -2 2)))
    1.75 +(sample-distribution 10 (sample-mean-variance 10 (interval -2 2)))
    1.76 +