diff --git a/collects/math/private/distributions/discrete-dist.rkt b/collects/math/private/distributions/discrete-dist.rkt index cc5a7a5cc7..4b3f6d0c2c 100644 --- a/collects/math/private/distributions/discrete-dist.rkt +++ b/collects/math/private/distributions/discrete-dist.rkt @@ -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) diff --git a/collects/math/private/distributions/dist-struct.rkt b/collects/math/private/distributions/dist-struct.rkt index b2262ef39a..b268534388 100644 --- a/collects/math/private/distributions/dist-struct.rkt +++ b/collects/math/private/distributions/dist-struct.rkt @@ -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_ab) -inf.0] - [else (lg- (fllog 0.5) log-P_x>b)])) - (lg+ log-P_a (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_ab) -inf.0] + [else (lg- (fllog 0.5) log-P_x>b)])) + (lg+ log-P_a (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?)])))) diff --git a/collects/math/private/distributions/truncated-dist.rkt b/collects/math/private/distributions/truncated-dist.rkt index 52fdfa9a13..e7d4ad5d2d 100644 --- a/collects/math/private/distributions/truncated-dist.rkt +++ b/collects/math/private/distributions/truncated-dist.rkt @@ -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_ab (delay (orig-cdf b #t #t))) diff --git a/collects/math/private/distributions/utils.rkt b/collects/math/private/distributions/utils.rkt index c9a3b26a0e..247e7ea428 100644 --- a/collects/math/private/distributions/utils.rkt +++ b/collects/math/private/distributions/utils.rkt @@ -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) diff --git a/collects/math/private/flonum/flonum-more-functions.rkt b/collects/math/private/flonum/flonum-more-functions.rkt index 1a4e133fda..a8d21903e8 100644 --- a/collects/math/private/flonum/flonum-more-functions.rkt +++ b/collects/math/private/flonum/flonum-more-functions.rkt @@ -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))] diff --git a/collects/math/private/functions/psi.rkt b/collects/math/private/functions/psi.rkt index 212052eb74..d57d2c5e81 100644 --- a/collects/math/private/functions/psi.rkt +++ b/collects/math/private/functions/psi.rkt @@ -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) diff --git a/collects/math/private/functions/zeta.rkt b/collects/math/private/functions/zeta.rkt index a4fd030d68..fc91d296bb 100644 --- a/collects/math/private/functions/zeta.rkt +++ b/collects/math/private/functions/zeta.rkt @@ -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))) diff --git a/collects/math/scribblings/math-distributions.scrbl b/collects/math/scribblings/math-distributions.scrbl index 82aa14a8d0..3f62bd1291 100644 --- a/collects/math/scribblings/math-distributions.scrbl +++ b/collects/math/scribblings/math-distributions.scrbl @@ -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")] diff --git a/collects/math/scribblings/math-flonum.scrbl b/collects/math/scribblings/math-flonum.scrbl index b3c89d3cbf..dac1a53154 100644 --- a/collects/math/scribblings/math-flonum.scrbl +++ b/collects/math/scribblings/math-flonum.scrbl @@ -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}))]