The point of writing a generate model is typically not to simply run
it to generate data, but rather to do inference on it—
A sampler is an object that contains a probabilistic program and produces samples from the posterior distribution of its result expression. Samplers can be either unweighted (where each sample is equally representative of the posterior distribution) or weighted (where each sample comes with a factor that corrects for the ratio of its sampling frequency and its probability in the posterior distribution).
exhaustive enumeration (up to optional threshold), exponential in number of random variables, cannot handle sampling from continuous random variables
The examples in the following sections use the following function definitions:
> (define (count-heads n) (if (zero? n) 0 (+ (if (flip) 1 0) (count-heads (sub1 n)))))
> (define (geom) (if (flip) 0 (add1 (geom))))
A sampler is also an applicable object. That is, if s is a sampler, then evaluating (s) generates a sample.
(generate-samples s n [ f #:burn burn #:thin thin]) → (vectorof any/c) s : weighted-sampler? n : exact-nonnegative-integer? f : (-> any/c any/c) = values burn : exact-nonnegative-integer? = 0 thin : exact-nonnegative-integer? = 0
The sampler is first called burn times and the results are discarded. In addition, the sampler is called thin times before every sample to be retained.
(generate-weighted-samples s n [ f #:burn burn #:thin thin]) → (vectorof (cons/c any/c (>=/c 0))) s : weighted-sampler? n : exact-nonnegative-integer? f : (-> any/c any/c) = values burn : exact-nonnegative-integer? = 0 thin : exact-nonnegative-integer? = 0
(sampler->discrete-dist sampler n [ f #:burn burn #:thin thin]) → discrete-dist? sampler : weighted-sampler? n : exact-positive-integer? f : (-> any/c any/c) = (lambda (x) x) burn : exact-nonnegative-integer? = 0 thin : exact-nonnegative-integer? = 0
The samplers supported by this language consist of simple samplers and a more complicated and flexible Metropolis-Hastings sampler framework.
(rejection-sampler def/expr ... result-expr)
The sampler is implemented using rejection sampling—
(importance-sampler def/expr ... result-expr)
Metropolis-Hastings (MH) is an algorithm framework for producing a correlated sequence of samples where each sample is based on the previous. The algorithm is parameterized by the mechanism for proposing a new state given the previous state; given a proposal, the MH algorithm accepts or rejects it based on how the proposal was generated and the relative likelihood of the proposed state.
The mh-sampler form implements the Metropolis-Hastings framework, and the proposal mechanisms are implemented by a variety of MH transition types.
(mh-sampler maybe-transition def/expr ... result-expr)
| #:transition transition-expr
The transition-expr determines the mechanism used to propose new states. If absent, (single-site) is used.
A zone pattern matches a zone if the two values are equal? or if the zone pattern is #f.
(slice [ #:scale scale-factor #:zone zone-pattern]) → mh-transition? scale-factor : (>/c 0) = 1 zone-pattern : any/c = #f
The scale-factor argument controls the width parameter used to find the slice bounds.
Note: unlike traditional Gibbs sampling, this transition picks a choice at random rather than perturbing all choices round-robin.
txs : (listof mh-transition?) weights : (listof (>/c 0)) = ( (lambda (tx) 1) txs)
The initial value is proposal:drift.
Random choices made during the execution of a model can be grouped into zones using the with-zone form. Specialized transitions can then be assigned to choices based on their membership in zones.
(with-zone zone-expr body ...+)
A particle set consists of a collection of particles, each of which contains a current state estimate and a likelihood weight. Particle sets are updated via a (stochastic) state transformer that produces a new state for each particle. Observations performed by the state transformer adjust the particles’ weights.
n : exact-nonnegative-integer? init-state : any/c
ps : particles? update-state : (-> any/c any/c)
ps : particles? n : exact-nonnegative-integer? = (particles-count ps) algorithm : (or/c 'multinomial 'residual #f) = 'multinomial
ps : particles?
ps : particles?
In general, it is only sensible to call this function when the weights are known to be equal, such as after calling particles-resample.
The enumerate form works by exploring all possibilities using
the technique described in [EPP]. If limit-expr is
absent or evaluates to #f, then exploration ceases only when
all paths have been explored; if any path is infinite, then
enumerate fails to terminate. If limit-expr is
given and evaluates to a probability limit, then exploration
ceases when the error of the normalized result distribution would be
less than limit—
Only discrete and integer-valued distributions can be sampled with enumerate, and infinite-range distributions (such as geometric-dist) and infinitely-deep recursive functions (such as geometric) require the use of #:limit for termination.
Use #:limit to control when exploration should be cut off.
> (enumerate #:limit 1e-06 (geom))
Use #:normalize? #f to get the result probabilities without normalizing by the acceptance rate:
> (enumerate #:normalize? #f (define (drop-coin?) (flip 0.9)) (define (drunk-flips n) (cond [(zero? n) #t] [(drop-coin?) 'failed] [else (and (flip) (drunk-flips (sub1 n)))])) (define A (drunk-flips 10)) (observe/fail (not (eq? A 'failed))) (eq? A #t))
(discrete-measure [#f 1.0] [#t 1.855468750000192e-12])