Renamed make-flexp/base' to make-flexpt'

Renamed `dist' struct type to `distribution' ("dist" is too common)
This commit is contained in:
Neil Toronto 2012-12-03 22:43:47 -07:00
parent 6ca52be0ae
commit 055512b4e8
9 changed files with 228 additions and 307 deletions

View File

@ -14,8 +14,8 @@
discrete-dist-probs)
(begin-encourage-inline
(struct: (In Out) discrete-dist-struct dist ([values : (Listof Out)]
[probs : (Listof Positive-Flonum)])
(struct: (In Out) discrete-dist-struct distribution ([values : (Listof Out)]
[probs : (Listof Positive-Flonum)])
#:property prop:custom-print-quotable 'never
#:property prop:custom-write
(λ (v port mode)

View File

@ -1,166 +1,143 @@
#lang racket/base
#lang typed/racket/base
(require racket/performance-hint
racket/promise
"../../flonum.rkt")
(provide
(all-from-out (submod "." typed-defs))
;; Rename transformers
dist-cdf
dist-inv-cdf
dist-min
dist-max)
PDF Sample CDF Inverse-CDF
(struct-out distribution)
(struct-out ordered-dist)
Real-Dist
pdf sample cdf inv-cdf real-dist-prob)
(module typed-defs typed/racket/base
(require racket/performance-hint
racket/promise
"../../flonum.rkt")
(define-type (PDF In)
(case-> (In -> Flonum)
(In Any -> Flonum)))
(define-type (Sample Out)
(case-> (-> Out)
(Integer -> (Listof Out))))
(define-type (CDF In)
(case-> (In -> Flonum)
(In Any -> Flonum)
(In Any Any -> Flonum)))
(define-type (Inverse-CDF Out)
(case-> (Real -> Out)
(Real Any -> Out)
(Real Any Any -> Out)))
(struct: (In Out) distribution ([pdf : (PDF In)]
[sample : (Sample Out)])
#:transparent)
(struct: (In Out) ordered-dist distribution
([cdf : (CDF In)]
[inv-cdf : (Inverse-CDF Out)]
[min : Out]
[max : Out]
[median : (Promise Out)])
#:transparent)
(define-type Real-Dist (ordered-dist Real Flonum))
;; =================================================================================================
(begin-encourage-inline
(provide
PDF Sample CDF Inverse-CDF
(struct-out dist)
(struct-out ordered-dist)
Real-Dist
dist-median
pdf sample cdf inv-cdf real-dist-prob)
(: pdf (All (In Out) (case-> ((distribution In Out) In -> Flonum)
((distribution In Out) In Any -> Flonum))))
(define (pdf d v [log? #f])
((distribution-pdf d) v log?))
(define-type (PDF In)
(case-> (In -> Flonum)
(In Any -> Flonum)))
(: sample (All (In Out) (case-> ((distribution In Out) -> Out)
((distribution In Out) Integer -> (Listof Out)))))
(define sample
(case-lambda
[(d) ((distribution-sample d))]
[(d n) ((distribution-sample d) n)]))
(define-type (Sample Out)
(case-> (-> Out)
(Integer -> (Listof Out))))
(: cdf (All (In Out) (case-> ((ordered-dist In Out) In -> Flonum)
((ordered-dist In Out) In Any -> Flonum)
((ordered-dist In Out) In Any Any -> Flonum))))
(define (cdf d v [log? #f] [1-p? #f])
((ordered-dist-cdf d) v log? 1-p?))
(define-type (CDF In)
(case-> (In -> Flonum)
(In Any -> Flonum)
(In Any Any -> Flonum)))
(define-type (Inverse-CDF Out)
(case-> (Real -> Out)
(Real Any -> Out)
(Real Any Any -> Out)))
(struct: (In Out) dist ([pdf : (PDF In)]
[sample : (Sample Out)])
#:transparent)
(struct: (In Out) ordered-dist dist
([cdf : (CDF In)]
[inv-cdf : (Inverse-CDF Out)]
[min : Out]
[max : Out]
[median : (Promise Out)])
#:transparent)
(define-type Real-Dist (ordered-dist Real Flonum))
;; =================================================================================================
(begin-encourage-inline
(: dist-median (All (In Out) ((ordered-dist In Out) -> Out)))
(define (dist-median d) (force (ordered-dist-median d)))
(: pdf (All (In Out) (case-> ((dist In Out) In -> Flonum)
((dist In Out) In Any -> Flonum))))
(define (pdf d v [log? #f])
((dist-pdf d) v log?))
(: sample (All (In Out) (case-> ((dist In Out) -> Out)
((dist In Out) Integer -> (Listof Out)))))
(define sample
(case-lambda
[(d) ((dist-sample d))]
[(d n) ((dist-sample d) n)]))
(: cdf (All (In Out) (case-> ((ordered-dist In Out) In -> Flonum)
((ordered-dist In Out) In Any -> Flonum)
((ordered-dist In Out) In Any Any -> Flonum))))
(define (cdf d v [log? #f] [1-p? #f])
((ordered-dist-cdf d) v log? 1-p?))
(: inv-cdf (All (In Out) (case-> ((ordered-dist In Out) Real -> Out)
((ordered-dist In Out) Real Any -> Out)
((ordered-dist In Out) Real Any Any -> Out))))
(define (inv-cdf d p [log? #f] [1-p? #f])
((ordered-dist-inv-cdf d) p log? 1-p?))
)
(: real-dist-prob* (Real-Dist Flonum Flonum Any -> Flonum))
;; Assumes a <= b
(define (real-dist-prob* d a b 1-p?)
(define c (dist-median d))
(define cdf (ordered-dist-cdf d))
(define p
(cond [(a . fl= . b) (if 1-p? 1.0 0.0)]
[1-p? (fl+ (cdf a #f #f) (cdf b #f #t))]
[(b . fl<= . c)
;; Both less than the median; use lower tail only
(fl- (cdf b #f #f) (cdf a #f #f))]
[(a . fl>= . c)
;; Both greater than the median; use upper tail only
(fl- (cdf a #f #t) (cdf b #f #t))]
[else
;; Median between a and b; use lower for (a,c] and upper for (c,b]
(define P_x<=a (cdf a #f #f))
(define P_x>b (cdf b #f #t))
(fl+ (fl- 0.5 P_x<=a) (fl- 0.5 P_x>b))]))
(max 0.0 (min 1.0 p)))
(: real-dist-log-prob* (Real-Dist Flonum Flonum Any -> Flonum))
;; Assumes a <= b
(define (real-dist-log-prob* d a b 1-p?)
(define c (dist-median d))
(define cdf (ordered-dist-cdf d))
(define log-p
(cond [(a . fl= . b) (if 1-p? 0.0 -inf.0)]
[1-p? (lg+ (cdf a #t #f) (cdf b #t #t))]
[(b . fl<= . c)
;; Both less than the median; use lower tail only
(define log-P_x<=a (cdf a #t #f))
(define log-P_x<=b (cdf b #t #f))
(cond [(log-P_x<=b . fl< . log-P_x<=a) -inf.0]
[else (lg- log-P_x<=b log-P_x<=a)])]
[(a . fl>= . c)
;; Both greater than the median; use upper tail only
(define log-P_x>a (cdf a #t #t))
(define log-P_x>b (cdf b #t #t))
(cond [(log-P_x>a . fl< . log-P_x>b) -inf.0]
[else (lg- log-P_x>a log-P_x>b)])]
[else
;; Median between a and b; try 1-upper first
(define log-P_x<=a (cdf a #t #f))
(define log-P_x>b (cdf b #t #t))
(define log-p (lg1- (lg+ log-P_x<=a log-P_x>b)))
(cond [(log-p . fl> . (log 0.1)) log-p]
[else
;; Subtracting from 1.0 (in log space) lost bits; split and add instead
(define log-P_a<x<=c
(cond [((fllog 0.5) . fl< . log-P_x<=a) -inf.0]
[else (lg- (fllog 0.5) log-P_x<=a)]))
(define log-P_c<x<=b
(cond [((fllog 0.5) . fl< . log-P_x>b) -inf.0]
[else (lg- (fllog 0.5) log-P_x>b)]))
(lg+ log-P_a<x<=c log-P_c<x<=b)])]))
(min 0.0 log-p))
(: real-dist-prob (case-> (Real-Dist Real Real -> Flonum)
(Real-Dist Real Real Any -> Flonum)
(Real-Dist Real Real Any Any -> Flonum)))
(define (real-dist-prob d a b [log? #f] [1-p? #f])
(let ([a (fl a)] [b (fl b)])
(let ([a (flmin a b)] [b (flmax a b)])
(cond [log? (define p (real-dist-prob* d a b 1-p?))
(cond [(and (p . fl> . +max-subnormal.0) (p . fl< . 0.9)) (fllog p)]
[else (real-dist-log-prob* d a b 1-p?)])]
[else (real-dist-prob* d a b 1-p?)]))))
(: inv-cdf (All (In Out) (case-> ((ordered-dist In Out) Real -> Out)
((ordered-dist In Out) Real Any -> Out)
((ordered-dist In Out) Real Any Any -> Out))))
(define (inv-cdf d p [log? #f] [1-p? #f])
((ordered-dist-inv-cdf d) p log? 1-p?))
)
(require (submod "." typed-defs)
(for-syntax racket/base))
(: real-dist-prob* (Real-Dist Flonum Flonum Any -> Flonum))
;; Assumes a <= b
(define (real-dist-prob* d a b 1-p?)
(define c (force (ordered-dist-median d)))
(define cdf (ordered-dist-cdf d))
(define p
(cond [(a . fl= . b) (if 1-p? 1.0 0.0)]
[1-p? (fl+ (cdf a #f #f) (cdf b #f #t))]
[(b . fl<= . c)
;; Both less than the median; use lower tail only
(fl- (cdf b #f #f) (cdf a #f #f))]
[(a . fl>= . c)
;; Both greater than the median; use upper tail only
(fl- (cdf a #f #t) (cdf b #f #t))]
[else
;; Median between a and b; use lower for (a,c] and upper for (c,b]
(define P_x<=a (cdf a #f #f))
(define P_x>b (cdf b #f #t))
(fl+ (fl- 0.5 P_x<=a) (fl- 0.5 P_x>b))]))
(max 0.0 (min 1.0 p)))
(define-syntax dist-cdf (make-rename-transformer #'ordered-dist-cdf))
(define-syntax dist-inv-cdf (make-rename-transformer #'ordered-dist-inv-cdf))
(define-syntax dist-min (make-rename-transformer #'ordered-dist-min))
(define-syntax dist-max (make-rename-transformer #'ordered-dist-max))
(: real-dist-log-prob* (Real-Dist Flonum Flonum Any -> Flonum))
;; Assumes a <= b
(define (real-dist-log-prob* d a b 1-p?)
(define c (force (ordered-dist-median d)))
(define cdf (ordered-dist-cdf d))
(define log-p
(cond [(a . fl= . b) (if 1-p? 0.0 -inf.0)]
[1-p? (lg+ (cdf a #t #f) (cdf b #t #t))]
[(b . fl<= . c)
;; Both less than the median; use lower tail only
(define log-P_x<=a (cdf a #t #f))
(define log-P_x<=b (cdf b #t #f))
(cond [(log-P_x<=b . fl< . log-P_x<=a) -inf.0]
[else (lg- log-P_x<=b log-P_x<=a)])]
[(a . fl>= . c)
;; Both greater than the median; use upper tail only
(define log-P_x>a (cdf a #t #t))
(define log-P_x>b (cdf b #t #t))
(cond [(log-P_x>a . fl< . log-P_x>b) -inf.0]
[else (lg- log-P_x>a log-P_x>b)])]
[else
;; Median between a and b; try 1-upper first
(define log-P_x<=a (cdf a #t #f))
(define log-P_x>b (cdf b #t #t))
(define log-p (lg1- (lg+ log-P_x<=a log-P_x>b)))
(cond [(log-p . fl> . (log 0.1)) log-p]
[else
;; Subtracting from 1.0 (in log space) lost bits; split and add instead
(define log-P_a<x<=c
(cond [((fllog 0.5) . fl< . log-P_x<=a) -inf.0]
[else (lg- (fllog 0.5) log-P_x<=a)]))
(define log-P_c<x<=b
(cond [((fllog 0.5) . fl< . log-P_x>b) -inf.0]
[else (lg- (fllog 0.5) log-P_x>b)]))
(lg+ log-P_a<x<=c log-P_c<x<=b)])]))
(min 0.0 log-p))
(: real-dist-prob (case-> (Real-Dist Real Real -> Flonum)
(Real-Dist Real Real Any -> Flonum)
(Real-Dist Real Real Any Any -> Flonum)))
(define (real-dist-prob d a b [log? #f] [1-p? #f])
(let ([a (fl a)] [b (fl b)])
(let ([a (flmin a b)] [b (flmax a b)])
(cond [log? (define p (real-dist-prob* d a b 1-p?))
(cond [(and (p . fl> . +max-subnormal.0) (p . fl< . 0.9)) (fllog p)]
[else (real-dist-log-prob* d a b 1-p?)])]
[else (real-dist-prob* d a b 1-p?)]))))

View File

@ -24,16 +24,16 @@
[(d a) (truncated-dist d -inf.0 a)]
[(d a b)
(let*-values ([(a b) (values (fl a) (fl b))]
[(a b) (values (max (dist-min d) (min a b))
(min (dist-max d) (max a b)))])
[(a b) (values (max (ordered-dist-min d) (min a b))
(min (ordered-dist-max d) (max a b)))])
(unsafe-truncated-dist d a b))]))
(: unsafe-truncated-dist (Real-Dist Float Float -> Truncated-Dist))
(define (unsafe-truncated-dist d a b)
(define orig-pdf (dist-pdf d))
(define orig-cdf (dist-cdf d))
(define orig-inv-cdf (dist-inv-cdf d))
(define orig-sample (dist-sample d))
(define orig-pdf (distribution-pdf d))
(define orig-cdf (ordered-dist-cdf d))
(define orig-inv-cdf (ordered-dist-inv-cdf d))
(define orig-sample (distribution-sample d))
(define log-P_a<x<=b (real-dist-prob d a b #t #f))
(define log-P_x<=a (delay (orig-cdf a #t #f)))
(define log-P_x>b (delay (orig-cdf b #t #t)))

View File

@ -53,45 +53,6 @@
(: proc-names (type-name -> arg-types)) ...
(define (proc-names v) (struct-proc-names v)) ...))))]))
(define-syntax (define-distribution-type: stx)
(syntax-case stx (:)
[(_ (type-name T ...) (parent-type-name In Out)
(A B) name ([arg-names : arg-type] ...))
(let ([arg-name-lst (syntax->list #'(arg-names ...))])
(with-syntax* ([struct-type-name (format-id #'type-name "~a-struct" #'type-name)]
[(struct-proc-names ...) (map (λ (arg-name)
(format-id #'type-name
"~a-~a"
#'struct-type-name
arg-name))
arg-name-lst)]
[struct-pred-name (format-id #'type-name "~a?" #'struct-type-name)]
[make-name (format-id #'name "make-~a" #'name)]
[(proc-names ...) (map (λ (arg-name)
(format-id #'name "~a-~a" #'name arg-name))
arg-name-lst)]
[pred-name (format-id #'name "~a?" #'name)]
[format-str
(string-append "(~a "
(string-join (build-list (length arg-name-lst) (λ _ "~v")))
")")])
(syntax/loc stx
(begin
(struct: (A B) struct-type-name parent-type-name ([arg-names : arg-types] ...)
#:property prop:custom-print-quotable 'never
#:property prop:custom-write
(λ (v port write?)
(fprintf port format-str 'name (proc-names v) ...)))
(define-type (type-name T ...) (struct-type-name In Out))
(: proc-names ((type-name T ...) -> arg-types)) ...
(define (proc-names d) (struct-proc-names d)) ...
(define make-name struct-type-name)
(define pred-name struct-pred-name)))))]
[(_ type-name (parent-type-name In Out) name ([arg-names arg-opts ...] ...))
(syntax/loc stx
(define-distribution-type: (type-name) (parent-type-name In Out)
(A B) name ([arg-names arg-opts ...] ...)))]))
;; ===================================================================================================
;; One-sided scale family distributions (e.g. exponential)

View File

@ -11,7 +11,7 @@
(provide flsqrt1pm1
flsinh flcosh fltanh
flasinh flacosh flatanh
make-flexp/base flexpt+ flexpt1p)
make-flexpt flexpt+ flexpt1p)
;; ---------------------------------------------------------------------------------------------------
;; sqrt(1+x)-1
@ -163,8 +163,8 @@
(begin-encourage-inline
(: make-flexp/base (Positive-Exact-Rational -> (Flonum -> Flonum)))
(define (make-flexp/base b)
(: make-flexpt (Positive-Exact-Rational -> (Flonum -> Flonum)))
(define (make-flexpt b)
(define b-hi (fl b))
(define b-lo (fl (- (/ (inexact->exact b-hi) b) 1)))
(cond [(fl= b-lo 0.0) (λ: ([x : Flonum]) (flexpt b-hi x))]

View File

@ -211,7 +211,7 @@
(define pi.128 267257146016241686964920093290467695825/85070591730234615865843651857942052864)
(: flexppi (Flonum -> Flonum))
(define flexppi (make-flexp/base pi.128))
(define flexppi (make-flexpt pi.128))
(: flpsi (Integer Flonum -> Flonum))
(define (flpsi m x)

View File

@ -107,7 +107,7 @@ P. Borwein. An Efficient Algorithm for the Riemann Zeta Function.
(define pi.128 267257146016241686964920093290467695825/85070591730234615865843651857942052864)
(: flexppi (Flonum -> Flonum))
(define flexppi (make-flexp/base pi.128))
(define flexppi (make-flexpt pi.128))
(require math/bigfloat)
(: fleta (Flonum -> Flonum))
@ -124,7 +124,7 @@ P. Borwein. An Efficient Algorithm for the Riemann Zeta Function.
[(fleven? (fltruncate (* 0.5 s))) +inf.0]
[else -inf.0])]
[(s . fl< . -171.0)
(define f (make-flexp/base (/ (* euler.64 pi.64) (inexact->exact (- s)))))
(define f (make-flexpt (/ (* euler.64 pi.64) (inexact->exact (- s)))))
(* (/ 4.0 pi)
(flexp-stirling (- s))
(flsqrt (/ pi (* -0.5 s)))
@ -160,7 +160,7 @@ P. Borwein. An Efficient Algorithm for the Riemann Zeta Function.
[else (- (flexpm1 (fl* (fl- 1.0 s) (fllog 2.0))))]))
(fl/ (fleta s) c)]
[(and (s . fl> . -266.5) (s . fl< . 0.0)) ; s < 0 keeps TR happy
(define f (make-flexp/base (/ (* euler.64 pi.64) (inexact->exact (- s)))))
(define f (make-flexpt (/ (* euler.64 pi.64) (inexact->exact (- s)))))
(* (/ 4.0 pi)
(flexp-stirling (- s))
(flsqrt (/ pi (* -0.5 s)))

View File

@ -12,7 +12,7 @@
@(define untyped-eval (make-untyped-math-eval))
@interaction-eval[#:eval untyped-eval (require racket/list)]
@title[#:tag "dist"]{Probability Distributions}
@title[#:tag "dist" #:style 'toc]{Probability Distributions}
@(author-neil)
@defmodule[math/distributions]
@ -57,7 +57,7 @@ half-open interval (1/2,1], and computes another approximation from random sampl
This plots the pdf and a kernel density estimate of the pdf from random samples:
@interaction[#:eval untyped-eval
(plot (list (function (dist-pdf d) #:color 0 #:style 'dot)
(plot (list (function (distribution-pdf d) #:color 0 #:style 'dot)
(density xs))
#:x-label "x" #:y-label "density of N(2,5)")]
@ -69,7 +69,7 @@ the probabilities within the interval:
(define d-trunc (truncated-dist d -inf.0 5))
(real-dist-prob d-trunc 5 6)
(real-dist-prob d-trunc 0.5 1.0)
(plot (list (function (dist-pdf d-trunc) #:color 0 #:style 'dot)
(plot (list (function (distribution-pdf d-trunc) #:color 0 #:style 'dot)
(density (sample d-trunc 1000)))
#:x-label "x" #:y-label "density of T(N(2,5),-∞,5)")]
@ -132,7 +132,7 @@ When distribution object constructors receive parameters outside their domains,
(cdf (beta-dist 0 0) 1/2)
(inv-cdf (geometric-dist 1.1) 0.2)]
@section{Basic Distribution Types and Operations}
@section{Distribution Types and Operations}
@defform[(PDF In)]{
The type of probability density functions, or @tech{pdfs}, defined as
@ -150,38 +150,6 @@ When given a nonnegative integer @racket[n] as an argument, a sampling procedure
a length-@racket[n] list of independent, random samples.
}
@defstruct*[dist ([pdf (PDF In)] [sample (Sample Out)])]{
The parent type of @tech{distribution objects}. The @racket[In] type parameter is the data type a
distribution accepts as arguments to its @tech{pdf}. The @racket[Out] type parameter is the data
type a distribution returns as random samples.
@examples[#:eval untyped-eval
(dist? (discrete-dist '(a b c)))
(dist? (normal-dist))
((dist-pdf (normal-dist)) 0)
((dist-sample (normal-dist)))]
See @racket[pdf] and @racket[sample] for uncurried forms of @racket[dist-pdf] and
@racket[dist-sample].
}
@defproc[(pdf [d (dist In Out)] [v In] [log? Any #f]) Flonum]{
An uncurried form of @racket[dist-pdf]. When @racket[log?] is not @racket[#f], returns a
log density.
@examples[#:eval untyped-eval
(pdf (discrete-dist '(a b c) '(1 2 3)) 'a)
(pdf (discrete-dist '(a b c) '(1 2 3)) 'a #t)]
}
@defproc*[([(sample [d (dist In Out)]) Out]
[(sample [d (dist In Out)] [n Integer]) (Listof Out)])]{
An uncurried form of @racket[dist-sample].
@examples[#:eval untyped-eval
(sample (exponential-dist))
(sample (exponential-dist) 3)]
}
@section{Ordered Distribution Types and Operations}
@defform[(CDF In)]{
The type of cumulative distribution functions, or @tech{cdfs}, defined as
@racketblock[(case-> (In -> Flonum)
@ -200,17 +168,32 @@ For any function of this type, both optional arguments should default to @racket
interpreted as specified in the description of @racket[inv-cdf].
}
@defstruct*[(ordered-dist dist) ([cdf (CDF In)]
[inv-cdf (Inverse-CDF Out)]
[min Out]
[max Out]
[median (Promise Out)])]{
@defstruct*[distribution ([pdf (PDF In)] [sample (Sample Out)])]{
The parent type of @tech{distribution objects}. The @racket[In] type parameter is the data type a
distribution accepts as arguments to its @tech{pdf}. The @racket[Out] type parameter is the data
type a distribution returns as random samples.
@examples[#:eval untyped-eval
(distribution? (discrete-dist '(a b c)))
(distribution? (normal-dist))
((distribution-pdf (normal-dist)) 0)
((distribution-sample (normal-dist)))]
See @racket[pdf] and @racket[sample] for uncurried forms of @racket[distribution-pdf] and
@racket[distribution-sample].
}
@defstruct*[(ordered-dist distribution) ([cdf (CDF In)]
[inv-cdf (Inverse-CDF Out)]
[min Out]
[max Out]
[median (Promise Out)])]{
The parent type of @italic{ordered} @tech{distribution objects}.
Similarly to @racket[dist], the @racket[In] type parameter is the data type an ordered distribution
accepts as arguments to its @tech{pdf}, and the @racket[Out] type parameter is the data type an
ordered distribution returns as random samples. Additionally, its @tech{cdf} accepts values of type
@racket[In], and its @tech{inverse cdf} returns values of type @racket[Out].
Similarly to @racket[distribution], the @racket[In] type parameter is the data type an ordered
distribution accepts as arguments to its @tech{pdf}, and the @racket[Out] type parameter is the
data type an ordered distribution returns as random samples. Additionally, its @tech{cdf}
accepts values of type @racket[In], and its @tech{inverse cdf} returns values of type
@racket[Out].
@examples[#:eval untyped-eval
(ordered-dist? (discrete-dist '(a b c)))
@ -227,8 +210,24 @@ The parent type of real-valued distributions, such as any distribution returned
@racket[normal-dist]. Equivalent to the type @racket[(ordered-dist Real Flonum)].
}
@defproc[(pdf [d (dist In Out)] [v In] [log? Any #f]) Flonum]{
An uncurried form of @racket[distribution-pdf]. When @racket[log?] is not @racket[#f], returns a
log density.
@examples[#:eval untyped-eval
(pdf (discrete-dist '(a b c) '(1 2 3)) 'a)
(pdf (discrete-dist '(a b c) '(1 2 3)) 'a #t)]
}
@defproc*[([(sample [d (dist In Out)]) Out]
[(sample [d (dist In Out)] [n Integer]) (Listof Out)])]{
An uncurried form of @racket[distribution-sample].
@examples[#:eval untyped-eval
(sample (exponential-dist))
(sample (exponential-dist) 3)]
}
@defproc[(cdf [d (ordered-dist In Out)] [v In] [log? Any #f] [1-p? Any #f]) Flonum]{
An uncurried form of @racket[dist-cdf].
An uncurried form of @racket[ordered-dist-cdf].
When @racket[log?] is @racket[#f], @racket[cdf] returns a probability; otherwise, it returns a log
probability.
@ -238,7 +237,7 @@ When @racket[1-p?] is @racket[#f], @racket[cdf] returns a lower-tail probability
}
@defproc[(inv-cdf [d (ordered-dist In Out)] [p Real] [log? Any #f] [1-p? Any #f]) Out]{
An uncurried form of @racket[dist-inv-cdf].
An uncurried form of @racket[ordered-dist-inv-cdf].
When @racket[log?] is @racket[#f], @racket[inv-cdf] interprets @racket[p] as a probability; otherwise,
it interprets @racket[p] as a log probability.
@ -254,22 +253,6 @@ the two endpoints are swapped first.) The @racket[log?] and @racket[1-p?] argume
meaning of the return value in the same way as the corresponding arguments to @racket[cdf].
}
@deftogether[(@defproc[(dist-cdf [d (ordered-dist In Out)]) (CDF In)]
@defproc[(dist-inv-cdf [d (ordered-dist In Out)]) (Inverse-CDF Out)]
@defproc[(dist-min [d (ordered-dist In Out)]) Out]
@defproc[(dist-max [d (ordered-dist In Out)]) Out]
@defproc[(dist-median [d (ordered-dist In Out)]) Out])]{
The first four are synonyms for @racket[ordered-dist-cdf], @racket[ordered-dist-inv-cdf],
@racket[ordered-dist-min] and @racket[ordered-dist-max]. The last is equivalent to
@racket[(force (ordered-dist-median d))].
@examples[#:eval untyped-eval
(dist-min (normal-dist))
(dist-min (geometric-dist 1/3))
(dist-max (uniform-dist 1 2))
(dist-median (gamma-dist 4 2))]
}
@section{Finite Distribution Families}
@subsection{Unordered Discrete Distributions}
@ -299,7 +282,7 @@ random samples.
(define n 500)
(define h (samples->hash (sample d n)))
(plot (list (discrete-histogram
(map vector xs (map (dist-pdf d) xs))
(map vector xs (map (distribution-pdf d) xs))
#:x-min 0 #:skip 2 #:label "P[x]")
(discrete-histogram
(map vector xs (map (λ (x) (/ (hash-ref h x) n)) xs))
@ -354,11 +337,11 @@ on it are faster.
@examples[#:eval untyped-eval
(define d (bernoulli-dist 0.75))
(map (dist-pdf d) '(0 1))
(map (dist-cdf d) '(0 1))
(map (distribution-pdf d) '(0 1))
(map (ordered-dist-cdf d) '(0 1))
(define d (binomial-dist 1 0.75))
(map (dist-pdf d) '(0 1))
(map (dist-cdf d) '(0 1))]
(map (distribution-pdf d) '(0 1))
(map (ordered-dist-cdf d) '(0 1))]
}
@subsection{Binomial Distributions}
@ -375,9 +358,9 @@ probability of success.
@examples[#:eval untyped-eval
(define d (binomial-dist 15 0.6))
(plot (discrete-histogram
(map vector (build-list 16 values) (build-list 16 (dist-pdf d))))
(map vector (build-list 16 values) (build-list 16 (distribution-pdf d))))
#:x-label "number of successes" #:y-label "probability")
(plot (function-interval (λ (x) 0) (dist-cdf d) -0.5 15.5)
(plot (function-interval (λ (x) 0) (ordered-dist-cdf d) -0.5 15.5)
#:x-label "at-most number of successes" #:y-label "probability")]
}
@ -395,9 +378,9 @@ first success starting from zero.
@examples[#:eval untyped-eval
(define d (geometric-dist 0.25))
(plot (discrete-histogram
(map vector (build-list 16 values) (build-list 16 (dist-pdf d))))
(map vector (build-list 16 values) (build-list 16 (distribution-pdf d))))
#:x-label "first success index" #:y-label "probability")
(plot (function-interval (λ (x) 0) (dist-cdf d) -0.5 15.5)
(plot (function-interval (λ (x) 0) (ordered-dist-cdf d) -0.5 15.5)
#:x-label "at-most first success index" #:y-label "probability"
#:y-max 1)]
}
@ -415,9 +398,9 @@ independent events.
@examples[#:eval untyped-eval
(define d (poisson-dist 6.2))
(plot (discrete-histogram
(map vector (build-list 16 values) (build-list 16 (dist-pdf d))))
(map vector (build-list 16 values) (build-list 16 (distribution-pdf d))))
#:x-label "number of events" #:y-label "probability")
(plot (function-interval (λ (x) 0) (dist-cdf d) -0.5 15.5)
(plot (function-interval (λ (x) 0) (ordered-dist-cdf d) -0.5 15.5)
#:x-label "at-most number of events" #:y-label "probability"
#:y-max 1)]
}
@ -468,14 +451,14 @@ which must both be nonnegative.
(plot (for/list ([α (in-list '(1 2 3 1/2))]
[β (in-list '(1 3 1 1/2))]
[i (in-naturals)])
(function (dist-pdf (beta-dist α β))
(function (distribution-pdf (beta-dist α β))
#:color i #:label (format "Beta(~a,~a)" α β)))
#:x-min 0 #:x-max 1 #:y-max 4 #:y-label "density")
(plot (for/list ([α (in-list '(1 2 3 1/2))]
[β (in-list '(1 3 1 1/2))]
[i (in-naturals)])
(function (dist-cdf (beta-dist α β))
(function (ordered-dist-cdf (beta-dist α β))
#:color i #:label (format "Beta(~a,~a)" α β)))
#:x-min 0 #:x-max 1 #:y-label "probability")]
@ -502,7 +485,7 @@ Represents the Cauchy distribution family parameterized by mode and scale.
(plot (for/list ([m (in-list '(0 -1 0 2))]
[s (in-list '(1 1/2 2.25 0.7))]
[i (in-naturals)])
(function (dist-pdf (cauchy-dist m s))
(function (distribution-pdf (cauchy-dist m s))
#:color i #:label (format "Cauchy(~a,~a)" m s)))
#:x-min -8 #:x-max 8 #:y-label "density"
#:legend-anchor 'top-right)
@ -510,7 +493,7 @@ Represents the Cauchy distribution family parameterized by mode and scale.
(plot (for/list ([m (in-list '(0 -1 0 2))]
[s (in-list '(1 1/2 2.25 0.7))]
[i (in-naturals)])
(function (dist-cdf (cauchy-dist m s))
(function (ordered-dist-cdf (cauchy-dist m s))
#:color i #:label (format "Cauchy(~a,~a)" m s)))
#:x-min -8 #:x-max 8 #:y-label "probability")]
}
@ -529,7 +512,7 @@ Represents the family of distributions whose densities are Dirac delta functions
(pdf (delta-dist) 1)
(plot (for/list ([μ (in-list '(-1 0 1))]
[i (in-naturals)])
(function (dist-cdf (delta-dist μ))
(function (ordered-dist-cdf (delta-dist μ))
#:color i #:style i #:label (format "δ(~a)" μ)))
#:x-min -2 #:x-max 2 #:y-label "probability")]
}
@ -551,14 +534,14 @@ is the reciprocal of mean or scale. Construct exponential distributions from rat
@examples[#:eval untyped-eval
(plot (for/list ([μ (in-list '(2/3 1 2))]
[i (in-naturals)])
(function (dist-pdf (exponential-dist μ))
(function (distribution-pdf (exponential-dist μ))
#:color i #:label (format "Exponential(~a)" μ)))
#:x-min 0 #:x-max 5 #:y-label "density"
#:legend-anchor 'top-right)
(plot (for/list ([μ (in-list '(2/3 1 2))]
[i (in-naturals)])
(function (dist-cdf (exponential-dist μ))
(function (ordered-dist-cdf (exponential-dist μ))
#:color i #:label (format "Exponential(~a)" μ)))
#:x-min 0 #:x-max 5 #:y-label "probability"
#:legend-anchor 'bottom-right)]
@ -583,7 +566,7 @@ which is the reciprocal of scale. Construct gamma distributions from rates using
(plot (for/list ([k (in-list '(1 2 3 9))]
[s (in-list '(2 2 3 1/2))]
[i (in-naturals)])
(function (dist-pdf (gamma-dist k s))
(function (distribution-pdf (gamma-dist k s))
#:color i #:label (format "Gamma(~a,~a)" k s)))
#:x-min 0 #:x-max 15 #:y-label "density"
#:legend-anchor 'top-right)
@ -591,7 +574,7 @@ which is the reciprocal of scale. Construct gamma distributions from rates using
(plot (for/list ([k (in-list '(1 2 3 9))]
[s (in-list '(2 2 3 1/2))]
[i (in-naturals)])
(function (dist-cdf (gamma-dist k s))
(function (ordered-dist-cdf (gamma-dist k s))
#:color i #:label (format "Gamma(~a,~a)" k s)))
#:x-min 0 #:x-max 15 #:y-label "probability"
#:legend-anchor 'bottom-right)]
@ -621,7 +604,7 @@ and scale. In this parameterization, the variance is @racket[(* 1/3 (sqr (* pi s
(plot (for/list ([μ (in-list '(0 -1 0 2))]
[s (in-list '(1 1/2 2.25 0.7))]
[i (in-naturals)])
(function (dist-pdf (logistic-dist μ s))
(function (distribution-pdf (logistic-dist μ s))
#:color i #:label (format "Logistic(~a,~a)" μ s)))
#:x-min -8 #:x-max 8 #:y-label "density"
#:legend-anchor 'top-right)
@ -629,7 +612,7 @@ and scale. In this parameterization, the variance is @racket[(* 1/3 (sqr (* pi s
(plot (for/list ([μ (in-list '(0 -1 0 2))]
[s (in-list '(1 1/2 2.25 0.7))]
[i (in-naturals)])
(function (dist-cdf (logistic-dist μ s))
(function (ordered-dist-cdf (logistic-dist μ s))
#:color i #:label (format "Logistic(~a,~a)" μ s)))
#:x-min -8 #:x-max 8 #:y-label "probability")]
}
@ -652,14 +635,14 @@ which is the square of standard deviation. Construct normal distributions from v
(plot (for/list ([μ (in-list '(0 -1 0 2))]
[σ (in-list '(1 1/2 2.25 0.7))]
[i (in-naturals)])
(function (dist-pdf (normal-dist μ σ))
(function (distribution-pdf (normal-dist μ σ))
#:color i #:label (format "N(~a,~a)" μ σ)))
#:x-min -5 #:x-max 5 #:y-label "density")
(plot (for/list ([μ (in-list '(0 -1 0 2))]
[σ (in-list '(1 1/2 2.25 0.7))]
[i (in-naturals)])
(function (dist-cdf (normal-dist μ σ))
(function (ordered-dist-cdf (normal-dist μ σ))
#:color i #:label (format "N(~a,~a)" μ σ)))
#:x-min -5 #:x-max 5 #:y-label "probability")]
}
@ -685,7 +668,7 @@ before constructing the distribution object.
[b (in-list '(0 1 3))]
[m (in-list '(-2 0 2))]
[i (in-naturals)])
(function (dist-pdf (triangle-dist a b m)) #:color i
(function (distribution-pdf (triangle-dist a b m)) #:color i
#:label (format "Triangle(~a,~a,~a)" a b m)))
#:x-min -3.5 #:x-max 3.5 #:y-label "density")
@ -693,7 +676,7 @@ before constructing the distribution object.
[b (in-list '(0 1 3))]
[m (in-list '(-2 0 2))]
[i (in-naturals)])
(function (dist-cdf (triangle-dist a b m)) #:color i
(function (ordered-dist-cdf (triangle-dist a b m)) #:color i
#:label (format "Triangle(~a,~a,~a)" a b m)))
#:x-min -3.5 #:x-max 3.5 #:y-label "probability")]
@ -724,11 +707,11 @@ Samples are taken by applying the truncated distribution's @tech{inverse cdf} to
(define d (normal-dist))
(define t (truncated-dist d -2 1))
t
(plot (list (function (dist-pdf d) #:label "N(0,1)" #:color 0)
(function (dist-pdf t) #:label "T(N(0,1),-2,1)"))
(plot (list (function (distribution-pdf d) #:label "N(0,1)" #:color 0)
(function (distribution-pdf t) #:label "T(N(0,1),-2,1)"))
#:x-min -3.5 #:x-max 3.5 #:y-label "density")
(plot (list (function (dist-cdf d) #:label "N(0,1)" #:color 0)
(function (dist-cdf t) #:label "T(N(0,1),-2,1)"))
(plot (list (function (ordered-dist-cdf d) #:label "N(0,1)" #:color 0)
(function (ordered-dist-cdf t) #:label "T(N(0,1),-2,1)"))
#:x-min -3.5 #:x-max 3.5 #:y-label "probability")]
}
@ -754,14 +737,14 @@ they are swapped before constructing the distribution object.
(plot (for/list ([a (in-list '(-3 -1 -2))]
[b (in-list '(0 1 3))]
[i (in-naturals)])
(function (dist-pdf (uniform-dist a b)) #:color i
(function (distribution-pdf (uniform-dist a b)) #:color i
#:label (format "Uniform(~a,~a)" a b)))
#:x-min -3.5 #:x-max 3.5 #:y-label "density")
(plot (for/list ([a (in-list '(-3 -1 -2))]
[b (in-list '(0 1 3))]
[i (in-naturals)])
(function (dist-cdf (uniform-dist a b)) #:color i
(function (ordered-dist-cdf (uniform-dist a b)) #:color i
#:label (format "Uniform(~a,~a)" a b)))
#:x-min -3.5 #:x-max 3.5 #:y-label "probability")]

View File

@ -164,7 +164,7 @@ But see @racket[flexpt1p], which is more accurate still.
Like @racket[(flexpt (+ 1.0 x) y)], but accurate for any @racket[x] and @racket[y].
}
@defproc[(make-flexp/base [x Real]) (Flonum -> Flonum)]{
@defproc[(make-flexpt [x Real]) (Flonum -> Flonum)]{
Equivalent to @racket[(λ (y) (flexpt x y))] when @racket[x] is a flonum, but much more
accurate for large @racket[y] when @racket[x] cannot be exactly represented
by a flonum.
@ -181,9 +181,9 @@ the error is low near zero, but grows with distance from the origin:
(eval:result ""))
(eval:alts (flulp-error (flexpt pi y) pi^y)
(eval:result @racketresultfont{43.12619934359266}))]
Using @racket[make-flexp/base], the error is near rounding error everywhere:
Using @racket[make-flexpt], the error is near rounding error everywhere:
@interaction[#:eval untyped-eval
(eval:alts (define flexppi (make-flexp/base (bigfloat->rational pi.bf)))
(eval:alts (define flexppi (make-flexpt (bigfloat->rational pi.bf)))
(eval:result ""))
(eval:alts (flulp-error (flexppi y) pi^y)
(eval:result @racketresultfont{0.8738006564073412}))]