(aux hansei): reification, reflection and exploration
Introduction
Quoting from Embedded domain-specific languages for probabilistic programming:
❝
Broadly speaking, probabilistic programming languages are to express computations with degrees of uncertainty, which comes from the imprecision in input data, lack of the complete knowledge or is inherent in the domain. More precisely, the goal of probabilistic programming languages is to represent and automate reasoning about probabilistic models, which describe uncertain quantities -- random variables -- and relationships among them.
Here we focus on the dsl Hansei Embedded probabilistic domain-specific language HANSEI and the corresponding published paper Embedded Probabilistic Programming, In proceedings of the IFIP working conference on domain-specific languages, ed. Walid Taha. LNCS 5658, Springer, 2009, pp. 360-384., by Oleg Kiselyov and Chung-chieh Shan; moreover, an application to logic can be found in HANSEI as a Declarative Logic Programming Language.
Implementation
The first implementation of our language uses the probability
monad that represents a stochastic computation as a lazy search
tree. That is, our implementation uses the type constructor pV defined below.
type 'a vc = V of 'a | C of (unit -> 'a pV) and 'a pV = (prob * 'a vc) list[(p, V v)] denotes a deterministic successful outcomevwith the
probability mass p. A branch of the form V v is a leaf node
that describes a possible successful outcome, whereas a branch
of the form C thunk is not yet explored.
The intended meaning of a search tree of type 'a pV is a discrete probability distribution over values of type 'a.(module
(aux hansei)
(op/times
op/plus
op/subtract
op/divide
op/greater
probcc-map
probcc-verbose
probcc-τ
probcc-value
probcc-unit
probcc-mem-empty
probcc-thread-local
probcc-mem-new-loc
probcc-mem-ref
probcc-mem-set
probcc-mem-prompt
probcc-mem-resume
probcc-reify/0
probcc-distribution
probcc-impossible
probcc-reflect
probcc-bernoulli
probcc-coin
probcc-uniform
probcc-uniform/range
probcc-uniform/either
probcc-geometric
probcc-when
probcc-observe
probcc-geometric/count
probcc-geometric/bounded
probcc-letlazy
probcc-letlazy/nesting
probcc-memo
probcc-wref-new
probcc-wref-get
probcc-wref-set!
probcc-explore
probcc-normalize/total
probcc-normalize
probcc-reify
probcc-reify/exact/a
probcc-reify/exact
probcc-variable-elimination
λ-probcc-bucket
probcc-prob-of
probcc-at-least
probcc-map-dist
probcc-random-selector
probcc-dist-selector
probcc-nearly-one
probcc-explore/shallow
probcc-first-success
probcc-explore/bounded
probcc-sample-dist/rejection
probcc-sample-runner/times
probcc-sample-dist
probcc-sample-rejection/a
probcc-sample-rejection
probcc-sample-importance/a
probcc-sample-importance
probcc-statistics
probcc-timeit
probcc-reify/part
probcc-soft-cut
probcc-bind
probcc-leaves
probcc-dfs)
(import
scheme
(only (scheme base) make-parameter)
(chicken base)
(only (chicken sort) sort)
(only (chicken time) cpu-time)
(only (chicken format) printf)
(only (chicken flonum) fpsqrt)
(only (chicken fixnum) fx= fx< fx+)
(only (chicken condition)
make-property-condition
make-composite-condition
abort)
(only (chicken memory representation) record-instance? record->vector)
(only srfi-1 fold)
(aux base)
(aux continuation delimited)
(aux hansei pmap)
(only (aux hansei random)
ocaml-random-global
ocaml-random-init!
make-ocaml-random-source)
(only (aux hakaru10 random) rng-source? rng-double!))
(define op/times (make-parameter *))
(define op/plus (make-parameter +))
(define op/subtract (make-parameter -))
(define op/divide
(make-parameter
(λ (m n)
(if (zero? n)
(/ (exact->inexact m) (exact->inexact n))
(exact->inexact (/ m n))))))
(define op/greater (make-parameter >))
(define (probcc-ratio m n)
(if (zero? n) (/ (exact->inexact m) (exact->inexact n)) (/ m n)))
(define (probcc-signal kind loc msg . args)
(abort (make-composite-condition
(make-property-condition
'exn
'location
loc
'message
msg
'arguments
args)
(make-property-condition 'hansei 'kind kind))))
(define (probcc-map f lst)
(if (null? lst)
'()
(let* ((h (f (car lst))) (t (probcc-map f (cdr lst)))) (cons h t))))
(define probcc-verbose (make-parameter #f))
(define-syntax-rule (probcc-τ p body ...) `((C ,(τ body ...)) ,p))
(define-syntax-rule (probcc-value p body ...) `((V ,(begin body ...)) ,p))
(define (probcc-unit v) (list (probcc-value 1 v)))
(define probcc-mem-empty (pmap-create fx-compare))
(define probcc-thread-local probcc-mem-empty)
(define probcc-mem-counter 0)
(define (probcc-mem-new-loc)
(let1 (loc probcc-mem-counter) (set! probcc-mem-counter (add1 loc)) loc))
(define (probcc-mem-ref loc mem #!optional default)
(pmap-find
loc
mem
(or default
(τ (probcc-signal 'memory 'probcc-mem-ref "unbound location" loc)))))
(define (probcc-mem-set loc v mem) (pmap-add loc v mem))
(define (probcc-mem-prompt thunk)
(let* ((mem probcc-thread-local) (v (thunk)))
(set! probcc-thread-local mem)
v))
(define (probcc-mem-resume curr k v)
(let1 (mem probcc-thread-local)
(set! probcc-thread-local curr)
(let1 (r (k v)) (set! probcc-thread-local mem) r)))
(define (probcc-reify/0 model)
(let* ((mem probcc-thread-local)
(tree (resetcc
(set! probcc-thread-local probcc-mem-empty)
(probcc-unit (model)))))
(set! probcc-thread-local mem)
tree))
(define (probcc-distribution pairs)
(let1 (curr probcc-thread-local)
(letcc/shift
k
(map (λ1-match/first
((,v ,p) (probcc-τ p (probcc-mem-resume curr k v))))
pairs))))
(define (probcc-impossible) (probcc-distribution '()))
(define (probcc-reflect choices)
(let1 (curr probcc-thread-local)
(letcc/shift
k
(letrec ((make-choices (λ (pv) (map F pv)))
(F (λ1-match/first
(((V ,v) ,p)
(probcc-τ p (probcc-mem-resume curr k v)))
(((C ,t) ,p) (probcc-τ p (make-choices (t)))))))
(make-choices choices)))))
(define (probcc-bernoulli t f p)
(probcc-distribution `((,t ,p) (,f ,((op/subtract) 1 p)))))
(define (probcc-coin p) (probcc-bernoulli #t #f p))
(define (probcc-uniform n)
(cond ((not (integer? n))
(probcc-signal 'malformed 'probcc-uniform "not an integer count" n))
((> n 0)
(letrec ((p (/ 1 n))
(plus (op/plus))
(subtract (op/subtract))
(loop (λ (pacc acc i)
(if (zero? i)
(probcc-distribution
(cons `(,i ,(subtract 1 pacc)) acc))
(loop (plus pacc p) (cons `(,i ,p) acc) (sub1 i))))))
(loop 0 '() (sub1 n))))
(else (probcc-impossible))))
(define (probcc-uniform/range low high)
(+ low (probcc-uniform (- high low))))
(define (probcc-uniform/either lst)
(list-ref lst (probcc-uniform (length lst))))
(define (probcc-geometric p s f)
(letrec ((subtract (op/subtract))
(loop (λ (n)
(list (probcc-τ p (probcc-unit (cons s n)))
(probcc-τ (subtract 1 p) (loop (cons f n)))))))
(probcc-reflect (loop '()))))
(define-syntax-rule
(probcc-when test body ...)
(cond (test body ...) (else (probcc-impossible))))
(define (probcc-observe test thunk)
(let1 (v (thunk)) (if (test v) v (probcc-impossible))))
(define (probcc-geometric/count p)
(letrec ((subtract (op/subtract))
(loop (λ (n)
(list (probcc-τ p (probcc-unit n))
(probcc-τ (subtract 1 p) (loop (add1 n)))))))
(probcc-reflect (loop 0))))
(define (probcc-geometric/bounded n p)
(letrec ((times (op/times))
(subtract (op/subtract))
(expi (λ (x n)
(cond ((= n 0) 1)
((= n 1) x)
(else (times x (expi x (sub1 n)))))))
(loop (λ (pp i)
(if (> i n)
'()
(cons (list i pp)
(loop (times (subtract 1 p) pp) (add1 i)))))))
(let1 (norm (subtract 1 (expi (subtract 1 p) (add1 n))))
(probcc-distribution (loop (probcc-ratio p norm) 0)))))
(define (probcc-letlazy e)
(let1 (loc (probcc-mem-new-loc))
(τ (probcc-mem-ref
loc
probcc-thread-local
(τ (let1 (v (e))
(set! probcc-thread-local
(probcc-mem-set loc v probcc-thread-local))
v))))))
(define (probcc-letlazy/nesting e)
(let1 (loc (probcc-mem-new-loc))
(set! probcc-thread-local
(probcc-mem-set loc `(MT ,e) probcc-thread-local))
(letrec ((deref (τ (match/first
(probcc-mem-ref loc probcc-thread-local (τ #f))
((MV ,v) v)
((MT ,e)
(let1 (v (e))
(set! probcc-thread-local
(probcc-mem-set loc `(MV ,v) probcc-thread-local))
v))
(#f
(let1 (mem probcc-thread-local)
(letcc/shift
k
(list (probcc-τ
1
(let1 (v (deref))
(probcc-mem-prompt
(τ (set! probcc-thread-local
(probcc-mem-set loc `(MV ,v) mem))
(k v)))))))))))))
deref)))
(define (probcc-memo f)
(let1 (loc (probcc-mem-new-loc))
(set! probcc-thread-local
(probcc-mem-set loc pmap-empty probcc-thread-local))
(letrec ((table (τ (probcc-mem-ref
loc
probcc-thread-local
(τ (probcc-signal
'memory
'probcc-memo
"memo table not in this world's memory"
loc))))))
(λ (x)
(pmap-find
x
(table)
(τ (let* ((v (f x)) (table* (pmap-add x v (table))))
(set! probcc-thread-local
(probcc-mem-set loc table* probcc-thread-local))
v)))))))
(define (probcc-wref-new v)
(let1 (loc (probcc-mem-new-loc))
(set! probcc-thread-local (probcc-mem-set loc v probcc-thread-local))
loc))
(define (probcc-wref-get r) (probcc-mem-ref r probcc-thread-local))
(define (probcc-wref-set! r v)
(set! probcc-thread-local (probcc-mem-set r v probcc-thread-local)))
(define (probcc-pmap->pv ans mass tail)
(pmap-foldi (λ (v p a) (cons (probcc-value (mass p) v) a)) ans tail))
(define (probcc-sort/mass pv)
(let1 (greater (op/greater))
(sort pv (λ (a b) (greater (cadr a) (cadr b))))))
(define (probcc-total plus ch) (fold (λ (e pa) (plus pa (cadr e))) 0 ch))
(define (probcc-scale ch ptotal)
(probcc-map
(λ1-match/first ((,x ,p) (list x (probcc-ratio p ptotal))))
ch))
(define (probcc-explore maxdepth choices)
(letrec ((times (op/times))
(plus (op/plus))
(loop (λ (p depth down choices ans susp)
(match/first
choices
(() (cons ans susp))
((((V ,v) ,pt) unquote rest)
(loop p
depth
down
rest
(pmap-insert-with plus v (times pt p) ans)
susp))
((((C ,t) ,pt) unquote rest)
(cond (down
(let1 (down* (or (not maxdepth) (< depth maxdepth)))
(match/first
(loop (times pt p)
(add1 depth)
down*
(t)
ans
susp)
((,ans* unquote susp*)
(loop p depth down rest ans* susp*)))))
(else
(let1 (susp* (cons `((C ,t) ,(times pt p)) susp))
(loop p depth down rest ans susp*)))))))))
(match/first
(loop 1 0 #t choices pmap-empty '())
((,ans unquote susp)
(probcc-sort/mass (probcc-pmap->pv ans identity susp))))))
(define (probcc-normalize/total choices)
(let* ((plus (op/plus))
(divide (op/divide))
(total (fold (λ (e acc) (plus (cadr e) acc)) 0 choices)))
(list total
(probcc-map
(λ1-match/first ((,slot ,p) (list slot (divide p total))))
choices))))
(define (probcc-normalize choices) (cadr (probcc-normalize/total choices)))
(define ((probcc-reify depth) model)
(probcc-explore depth (probcc-reify/0 model)))
(define probcc-reify/exact/a (probcc-reify +inf.0))
(define-syntax-rule
(probcc-reify/exact body ...)
(probcc-reify/exact/a (τ body ...)))
(define (probcc-variable-elimination f)
(λ args (probcc-reflect (probcc-reify/exact (apply f args)))))
(define-syntax-rule
(λ-probcc-bucket args body ...)
(letrec ((F (λ args body ...))
(B (λ-memo bargs (probcc-reify/exact (apply F bargs)))))
(o probcc-reflect B)))
(define (probcc-nan-inside? x)
(cond ((number? x) (or (nan? (real-part x)) (nan? (imag-part x))))
((pair? x)
(or (probcc-nan-inside? (car x)) (probcc-nan-inside? (cdr x))))
((vector? x)
(let loop ((i 0))
(and (fx< i (vector-length x))
(or (probcc-nan-inside? (vector-ref x i))
(loop (fx+ i 1))))))
((record-instance? x) (probcc-nan-inside? (record->vector x)))
(else #f)))
(define (probcc-structural=? a b)
(and (fx= 0 (pmap-compare a b)) (not (probcc-nan-inside? a))))
(define (probcc-prob-of v pv)
(let loop ((pv* pv))
(match/first
pv*
(() 0)
((((V ,w) ,p) unquote rest)
(if (probcc-structural=? w v) p (loop rest)))
((,_ unquote rest) (loop rest)))))
(define (probcc-at-least prob v pv) (>= (probcc-prob-of v pv) prob))
(define (probcc-map-dist f pv)
(let* ((plus (op/plus))
(pm (fold (λ (e pm)
(match/first
e
(((V ,v) ,p) (pmap-insert-with plus (f v) p pm))
(else (probcc-signal
'malformed
'probcc-map-dist
"a branch not yet explored"
e))))
pmap-empty
pv)))
(pmap-foldi (λ (v p a) (cons (probcc-value p v) a)) pm '())))
(define probcc-random-global-source
(make-ocaml-random-source ocaml-random-global))
(define (probcc-random-selector seed-or-source)
(let1 (src (cond ((rng-source? seed-or-source) seed-or-source)
(else
(ocaml-random-init! ocaml-random-global seed-or-source)
probcc-random-global-source)))
(λ (choices)
(let* ((plus (op/plus))
(ptotal (probcc-total plus choices))
(r (rng-double! src 0.0 ptotal)))
(let loop ((pcum 0) (rest choices))
(match/first
rest
(()
(probcc-signal
'selector
'probcc-random-selector
"Choice selection: can't happen"
choices))
(((,x ,p) unquote rest*)
(let1 (pcum* (plus pcum p))
(if (< r pcum*)
(list x ptotal)
(loop pcum* rest*))))))))))
(define (probcc-dist-selector choices)
(let* ((ptotal (probcc-total (op/plus) choices))
(x (probcc-distribution (probcc-scale choices ptotal))))
(list x ptotal)))
(define probcc-nearly-one (- 1.0 9.9999999999999995e-08))
(define (probcc-explore/shallow maxdepth choices)
(letrec ((times (op/times))
(plus (op/plus))
(loop (λ (pc depth ans acc lst)
(match/first
lst
(() (cons ans acc))
((((V ,v) ,p) unquote rest)
(loop pc
depth
(pmap-insert-with plus v (times p pc) ans)
acc
rest))
(((,c unquote rest) ⊣ (>= depth maxdepth))
(loop pc depth ans (cons c acc) rest))
((((C ,t) ,p) unquote rest)
(match/first
(loop (times pc p) (add1 depth) ans '() (t))
((,ans* unquote ch)
(let* ((ptotal (probcc-total plus ch))
(acc* (cond ((zero? ptotal) acc)
((< ptotal probcc-nearly-one)
(let1 (ch* (probcc-scale ch ptotal))
(cons `((C ,(τ ch*)) ,(times p ptotal)) acc)))
(else (cons `((C ,(τ ch)) ,p) acc)))))
(loop pc depth ans* acc* rest)))))))))
(match/first
(loop 1 0 pmap-empty '() choices)
((,ans unquote acc) (probcc-pmap->pv ans identity acc)))))
(define (probcc-first-success choices)
(let1 (times (op/times))
(let loop ((choices* choices))
(match/first
choices*
(() '())
((((V ,_) ,_) unquote _) choices*)
((((C ,t) ,pt) unquote rest)
(let1 (ch (probcc-map
(λ1-match/first
((,slot ,p) (list slot (times pt p))))
(t)))
(loop (append rest ch))))))))
(define (probcc-explore/bounded maxsize choices)
(letrec ((times (op/times))
(plus (op/plus))
(loop (λ (explore? pc low high q qsize lst)
(match/first
lst
(() (next low high q qsize))
((((V ,_) ,p) unquote rest)
(let1 (pe (times pc p))
(loop explore?
pc
(plus low pe)
(plus high pe)
q
qsize
rest)))
((((C ,t) ,p) unquote rest)
(if explore?
(let1 (q* (pmap-insert-with append (times pc p) (list t) q))
(loop explore? pc low high q* (add1 qsize) rest))
(loop explore?
pc
low
(plus high (times pc p))
q
qsize
rest))))))
(next (λ (low high q qsize)
(if (zero? qsize)
(list low high)
(let1 (explore? (< qsize maxsize))
(match/first
((if explore?
pmap-delete-find-max
pmap-delete-find-min)
q)
(((,p ,t unquote ts) ,q*)
(loop explore?
p
low
high
(if (null? ts) q* (pmap-add p ts q*))
(sub1 qsize)
(t)))))))))
(loop #t 1 0 0 (pmap-create number-compare) 0 choices)))
(define (probcc-sample-dist/rejection selector n choices)
(letrec ((times (op/times))
(plus (op/plus))
(loop (λ (pcontrib ans ch)
(match/first
ch
((((V ,v) ,p))
(pmap-insert-with plus v (times p pcontrib) ans))
(() ans)
((((C ,t) ,p)) (loop (times p pcontrib) ans (t)))
(else (match/first
(selector ch)
((,th ,ptotal)
(loop (times pcontrib ptotal)
ans
(list (list th 1)))))))))
(driver
(λ (ans i)
(if (zero? i) ans (driver (loop 1 ans choices) (sub1 i))))))
(let1 (ans (driver pmap-empty n))
(when (probcc-verbose)
(printf "rejection_sample: done ~a worlds\n" n))
(probcc-sort/mass
(probcc-pmap->pv ans (λ (p) (probcc-ratio p n)) '())))))
(define ((probcc-sample-runner/times n) z step)
(let loop ((z* z) (i n))
(if (zero? i) (list z* n) (loop (step z*) (sub1 i)))))
(define (probcc-sample-dist selector runner choices)
(letrec ((times (op/times))
(plus (op/plus))
(add (λ (v p ans) (pmap-insert-with plus v p ans)))
(look-ahead
(λ (pcontrib)
(λ (e answers)
(let ((ans (car answers)) (acc (cdr answers)))
(match/first
e
(((V ,v) ,p)
(cons (add v (times p pcontrib) ans) acc))
(((C ,t) ,p)
(match/first
(t)
(() answers)
((((V ,v) ,p1))
(cons (add v (times (times p p1) pcontrib) ans)
acc))
(,ch
(let1 (ptotal (probcc-total plus ch))
(cons ans
(cons (if (< ptotal probcc-nearly-one)
(list (probcc-scale ch ptotal) (times p ptotal))
(list ch p))
acc)))))))))))
(loop (λ (pcontrib ans ch)
(match/first
ch
((((V ,v) ,p)) (add v (times p pcontrib) ans))
(() ans)
((((C ,t) ,p)) (loop (times p pcontrib) ans (t)))
(else (match/first
(fold (look-ahead pcontrib) (list ans) ch)
((,ans*) ans*)
((,ans* unquote cch)
(match/first
(selector cch)
((,th ,ptotal)
(loop (times pcontrib ptotal) ans* th)))))))))
(toploop
(λ (pcontrib ans cch)
(match/first
(selector cch)
((,th ,ptotal) (loop (times pcontrib ptotal) ans th)))))
(driver
(λ (pcontrib vals cch)
(match/first
(runner pmap-empty (λ (ans) (toploop pcontrib ans cch)))
((,ans ,ns)
(let1 (ans* (pmap-foldi
(λ (v p ans) (add v (times ns p) ans))
vals
ans))
(when (probcc-verbose)
(printf "sample_importance: done ~a worlds\n" ns))
(probcc-sort/mass
(probcc-pmap->pv
ans*
(λ (p) (probcc-ratio p ns))
'())))))))
(make-threads
(λ (pcontrib ans ch)
(match/first
(fold (look-ahead pcontrib) (list ans) ch)
((,ans*)
(probcc-sort/mass (probcc-pmap->pv ans* identity '())))
((,ans* (,ch1 ,p))
(make-threads (times pcontrib p) ans* ch1))
((,ans* unquote cch)
(driver pcontrib ans* (reverse cch)))))))
(make-threads 1 pmap-empty choices)))
(define (probcc-sample-rejection/a selector n thunk)
(probcc-sample-dist/rejection selector n (probcc-reify/0 thunk)))
(define-syntax-rule
(probcc-sample-rejection selector n body ...)
(probcc-sample-rejection/a selector n (τ body ...)))
(define (probcc-sample-importance/a selector n thunk)
(probcc-sample-dist
selector
(probcc-sample-runner/times n)
(probcc-explore/shallow 3 (probcc-reify/0 thunk))))
(define-syntax-rule
(probcc-sample-importance selector n body ...)
(probcc-sample-importance/a selector n (τ body ...)))
(define (probcc-statistics seed1 seed2 sampler)
(let* ((times (op/times))
(plus (op/plus))
(subtract (op/subtract))
(add (λ (new old)
(list (plus (car old) (car new))
(plus (cadr old) (cadr new)))))
(answers
(let loop ((seed seed1) (answers pmap-empty))
(if (> seed seed2)
answers
(loop (add1 seed)
(fold (λ (e answers)
(match/first
e
(((V ,v) ,p)
(pmap-insert-with
add
v
(list p (times p p))
answers))
(else (probcc-signal
'malformed
'probcc-statistics
"a branch not yet explored"
e))))
answers
(sampler seed))))))
(n (add1 (- seed2 seed1)))
(stats (pmap-foldi
(λ (v s a)
(match/first
s
((,p ,p2)
(let* ((variance (/ (subtract p2 (/ (times p p) n)) n))
(std (fpsqrt (exact->inexact variance))))
(cons (list v (/ p n) std) a)))))
answers
'())))
(sort stats
(let1 (greater (op/greater))
(λ (a b) (greater (cadr a) (cadr b)))))))
(define (probcc-timeit thunk)
(let* ((seconds
(λ ()
(call-with-values
cpu-time
(λ (user system) (/ (+ user system) 1000.0)))))
(start (seconds))
(r (thunk)))
(printf "\nTime spent: ~a sec\n" (- (seconds) start))
r))
(define (probcc-bind m f)
(probcc-map
(λ1-match/first
(((V ,x) ,p) (probcc-τ p (f x)))
(((C ,t) ,p) (probcc-τ p (probcc-bind (t) f))))
m))
(define (probcc-reify/part bound model)
(foldr (λ (e acc)
(match/first e (((V ,v) ,p) (cons (list v p) acc)) (else acc)))
'()
(probcc-explore bound (probcc-reify/0 model))))
(define (probcc-soft-cut p q r)
(let1 (t (probcc-first-success (probcc-reify/0 p)))
(if (null? t) (r) (q (probcc-reflect t)))))
(define (probcc-leaves choices)
(let L ((choices* choices) (count 0))
(let1 (F (λ (probpair acc)
(match/first
probpair
(((V ,_) ,_) (add1 acc))
(((C ,t) ,_) (L (t) acc)))))
(foldr F count choices*))))
(define (probcc-dfs choices)
(letrec ((times (op/times))
(D (λ (choices) (apply append (probcc-map M choices))))
(M (λ1-match/first
(((C ,t) ,p)
(D (probcc-map
(λ1-match/first ((,slot ,pi) `(,slot ,(times p pi))))
(t))))
(,probpair (list probpair)))))
(D choices))))
Reification and reflection
The model is ordinary Scheme code; a choice is a call to probcc-distribution, which captures its continuation up to the
nearest resetcc and answers a list of branches, one thunk per
value, each resuming the continuation with its value (DSL p.14 l.627-628). probcc-reify/0 runs a model under resetcc and so turns it into a lazy search tree (DSL p.15 l.637):
let reify0 m = reset (fun () -> pv_unit (m ()))(aux continuation delimited), where probM.ml goes through a request type and its interpreter (deviation 2):❝
We pay the price of embedding stochastic computations only when calling the two stochastic operators dist and reify0. The deterministic parts of a stochastic program run at full speed, as if no randomness were present and no inference were in progress.
The lawn model of DSL p.12 l.539-546, for example, is a function of no
arguments: rain and a sprinkler may wet the grass, which is observed wet, and (probcc-reify/exact/a lawn) answers rain at 0.2838 and no rain
at 0.3219999999999999, before normalization (test/probcc/grass-model reifies the same body):
(define (lawn)
(let* ((rain (probcc-coin 0.29999999999999999))
(sprinkler (probcc-coin 0.5))
(grass-is-wet
(or (and (probcc-coin 0.90000000000000002) rain)
(and (probcc-coin 0.80000000000000004) sprinkler)
(probcc-coin 0.10000000000000001))))
(probcc-when grass-is-wet rain)))
probcc-reflect is the inverse of probcc-reify/0 (DSL p.16 l.694-700): it chooses following a tree, so that reifying a reflected
tree gives an equivalent tree, and reflecting a reified model expresses the same
distribution (DSL p.17 l.728-730; test/reflect/inverse). Reifying
a sub-model, flattening its tree with probcc-explore and reflecting it is variable elimination (probcc-variable-elimination, DSL p.17 l.748-771): flips_xor' 10 explores 9 trees of depth 2
instead of the 1024 leaves of flips_xor 10; memoizing the eliminated function on its arguments is bucket elimination (λ-probcc-bucket, DSL p.18 l.776-783; test/variable-elimination).
Exploration and its order
probcc-explore (I:11-28, DSL p.9 l.393-408) flattens a tree to a
given depth, or to the end when the depth is #f or +inf.0 (OCaml's None): the values reached are
collated, their masses summed, and the branches below the depth are kept as
suspended entries. Upstream threads the answer, a persistent map, and the list of
suspensions through its loop, and so does this port: a thunk forced by the explorer
may choose in an enclosing model, whose every world then resumes the loop, so a
mutable table would be shared by all of them (test/explore/persistent).
The answer is sorted by descending mass, with a stable sort of upstream's order
(deviation 5): equal masses keep the values in descending pmap-compare order, the stand-in for OCaml's compare (deviation 6), then the suspended branches in reverse
discovery order. So the order of the tables is deterministic (test/explore/order, test/uniform/edges).
World-local memory
Every world of a model has its own memory, which the lazy variables, the memo tables and the world-local references live in. It is probM's first-class store (P:40-75):
❝
The first-class storage is immutable and pure functional. Only the pointer to it is mutable. This is critical as different threads may have different views of the same `memory' and one thread's updating its memory should not affect other threads that share it.
❝
We implement first-class stores as immutable maps, and associate one with each delimited continuation.
probcc-thread-local is that pointer, a persistent map from
fixnum locations to values, changed only with set!: probcc-reify/0 starts the model in an empty memory and
restores the caller's when the first segment returns (P:104-107); probcc-distribution and probcc-reflect read the memory before capturing the continuation, and every branch resumes in
it (P:118-148). A parameter or a dynamic-wind would not do:
aux's shift enters the dynamic-wind frames of the captured continuation again on every resume, so a value
cached in one world would show in another. A parameterize inside a model is independent of the memory (test/parameterize/memory).
Laziness
Events take place in another order than they are observed, and a value drawn
long before it is looked at duplicates the search tree in between (DSL p.19-20
l.858-865). probcc-letlazy (P:321-330) answers a thunk that
evaluates its body at most once per world, when first called, keeping the value in
the world-local memory: "lazy let insertion at the point of reification"
(P:292-298). The built-in laziness cannot do it:
❝So
It is incorrect to implement letlazy using OCaml's built-in laziness, because the latter uses the mutable state of the underlying hardware, which does not `snap back' when the inference procedure backtracks to an earlier point in the execution of the probabilistic program.
delay, force and δ are never used for model values. probcc-letlazy is not nesting-safe, as upstream: an inner reification starts with an empty
memory, so a lazy variable of the outer world is evaluated again inside, and the
effects are reordered (N:65-81). probcc-letlazy/nesting (P:342-370) stores the unevaluated body where it is created and lets an inner
world ask its parent for it (test/letlazy, test/nested).Memoization
probcc-memo (P:449-479) generalizes the lazy variables to
stochastic functions:
❝Its table lives in the world-local memory, so a memo must be made inside the reification that uses it (deviation 16;
The function (memo f) :: a -> Prob b is then a proper function: applied to an argument, it gives a sample from its output distribution. When applied to the argument again, it gives the same sample. The memo facility is the generalization of letlazy: indeed, letlazy is the particular case of memo when the function to memoize has the argument unit.
test/memo). The table of λ-probcc-bucket is global instead, shared by
every inference, which is right only for a function that reads no world-local state
(deviation 11).Nested inference
An inference procedure is an ordinary function, so a model may call one, and the choices it makes are choices of the enclosing model:
❝An inner
Because a random choice made by an inference algorithm is expressed with dist just like any other random choice, any inference procedure, such as the importance-sampling algorithm in §4, can reason about itself or any other inference procedure.
probcc-reify/0 pushes its own delimiter and an empty
memory; the inner explorer forces its thunks in the extent of the outer model, so a
choice directly inside them (probcc-dist-selector, a lazy variable asking its parent) captures the inner explorer's loop, which
therefore threads persistent state, as every inference loop of the core does (plan
§5.7; test/nested, test/nested/thresholds). EPP p.4 l.374-379 memoizes the outer biased "as usual" and uses it at the inner level; with probM's memo that raises Not_found, and only letlazy_nesting works, here as upstream.Inference about inference
PE:69-77 asks how likely an estimate made by sampling is to be at least 0.3:
❝With
Suppose we choose a coin as before, then estimate p by flipping the coin twice and dividing the count of true by 2 (that is, by using rejection sampling with count 2). What is the probability that our estimate is at least 0.3? It is 7/8, because the only way for us to estimate below 0.3 is to choose a fair coin and get false from it twice.
probcc-dist-selector, the sampler's random choices are
choices of the outer model, which exact inference enumerates: the outer table is
exactly 7/8 and 1/8 (test/nested/inference-about-inference).Samplers and selectors
A selector chooses one entry of a list of (x p) entries and answers (x ptotal) (deviation 18). probcc-random-selector draws from OCaml's generator, which (aux hansei random) transcribes, and its seed re-initialises the
one global state, as Random.init does (P:221-222), so seeded
runs replay upstream's records bit for bit (see the suite of (aux hansei random)); probcc-dist-selector chooses through the enclosing model instead (P:238-241).
probcc-sample-dist/rejection (I:121-135) follows a path from the
root to a leaf per sample. probcc-sample-dist (I:137-198) is importance sampling with look-ahead:
❝On the drunk coin of DSL p.18-19, the and of ten tosses of a coin lost with probability 0.9, rejection sampling never sees the ten heads, whereas importance sampling estimates their probability, 9.77e-14 exactly (
The key is to explore the tree shallowly before committing to any random choice among its branches. Whenever this look-ahead comes across a branch that ends in failure, that branch is eliminated from the available choices. Similarly, whenever the look-ahead comes across a branch that ends in success, the successful result is registered as a sample, and that branch is again eliminated from the available choices.
test/drunk-coin), at 7e-14 from 5000 samples with seed 1 (the replays of the random suite). EPP
p.6 l.547-590 gives the algorithm step by step.(define (drunk-coin)
(let* ((x (probcc-coin 0.5)) (lost (probcc-coin 0.90000000000000002)))
(if lost (probcc-impossible) x)))
(define (dcoin-and n)
(if (= n 1) (drunk-coin) (and (drunk-coin) (dcoin-and (sub1 n)))))
Two CHICKEN hazards
Two properties of CHICKEN break models that resume a continuation more than once, which every model of these suites avoids (plan §5.6):
mapbuilds its result destructively (srfi-1's too), so a continuation captured inside its procedure must not be resumed:(probcc-reify/exact (map (λ (x) (probcc-coin 1/2)) '(1 2)))answers the single entry(#t #t #f #f #t #f)at mass 1, the four worlds sharing one list built in place, instead of four lists at 1/4. Models useprobcc-map, left to right as OCaml'sList.map(deviation 20);list-tabulateandfilter-mapare re-entrant but call their procedure from the last element to the first (test/map/reentrant);- under
csi, a continuation captured while the operands of a call with five or more operands are evaluated accumulates them when resumed, even with plaincall/cc:(probcc-reify/exact (list 'a 'b 'c 'd (probcc-coin 1/2)))answers(a b c d #t)and(a b c d #t #f); quasiquote templates and constructors are calls too. Models bind every choice withlet*first. Compiled code is not affected.
OCaml evaluates the components of a tuple and the arguments of a call from right
to left, let x = e1 and y = e2 in and the scrutinee of match (e1, e2) with from left to right; the ported models
sequence their choices with let* in that order, which the
depth-bounded and the replayed oracles depend on (plan D19).
Deviations
The departures from upstream, each tagged in the source with its number as deviation N (plan §9):
- the representation: entries
((V v) p)and((C thunk) p)and pairs(v p)put the probability second; weights may be exact, and the operations on them are the parametersop/times..op/greater(T:4-15); - the direct reifier: no
pReq,interporread_answer, the body of the reset answers the tree already interpreted (P:24, 97-108); - a single untagged delimiter: a user's
resetccinside a model captures HANSEI's choices (P:26); shift, notshift0, which coincide since no shift body shifts again (P:120, 142, 361);- the output order of exploration: descending mass, stable over upstream's
order;
+inf.0or#fforNone;probcc-reify/partkeeps the explorer's order (I:28; LV:68-70); pmap-comparefor OCaml'scompare: symbols are ordered by name, not by declaration, and unorderable values raise (PM:101);probcc-uniform 1makes a choice node of one branch,n <= 0fails instead of raising, a count that is not an integer raises(exn hansei malformed), and the masses are the exact 1/n (P:156-165);probcc-uniform/rangeis half open, low .. high - 1 (P:173-174);probcc-geometricanswers the trace of the trials;probcc-geometric/countis upstream's (P:181-185);- normalization turns exact weights into flonums with one rounding (I:299-301);
- the bucket of
λ-probcc-bucketis the globalλ-memo, wrong when the function reads world-local state (DSL p.18 l.776-783); probcc-variable-eliminationis variadic (P:267-268);probcc-letlazy/nestingprints no debugging output (P:347-363);- the samplers print their
donelines only underprobcc-verbose(I:131, 186); - the memory stores values directly in a map of fixnum locations, with one global pointer, single-threaded (P:40-64);
- a missing memory location or memo table raises
(exn hansei memory)instead ofNot_found, and memo tables are ordered bypmap-compare(P:471); - the world-local references are reconstructed, since
WorldLocalMemis not shipped, andprobcc-wref-set!takes the reference first (TI:121-126); - a selector answers the list
(x ptotal); a seed re-initialises the global OCaml-compatible state, and an rng-source may be given instead (P:221-233; T:20); - the ratios
p /. ptotalkeep exact weights exact, and a zero total gives NaN or an infinity as upstream (P:205, 240; I:47, 132, 160, 187); probcc-mapmust replacemaparound stochastic code (new);probcc-statisticsis sorted, not in hash order (I:292);- the parser combinators are prefixed
parsec-;parsec-oris binary,parsec-altandparsec-stream-overtake lists, and the unit value is(void)(PS:59-66, 106-113); probcc-prob-ofskips branches not yet explored (N:292-293);- lazy lists are
'()and(head-thunk . tail), promoted to a library, andprobcc-lobservecompares heads withequal?(SL:325-326; M2:116); - OCaml's tuples are lists (D18);
- the exact
1and0stand for upstream's1.0and0.0in deterministic branches, folds and missing values (P:100, 182, 198, 361; I:108, 128, 300; N:293); - the standard deviation of
probcc-statisticsisfpsqrtof a flonum, so a negative rounded variance gives+nan.0as OCaml (I:289-292); probcc-timeitprints the seconds with~a, not%g(I:307);probcc-map-distis promoted from example code, and a branch not yet explored raises(exn hansei malformed)(BL:456-460);- the errors of
(aux hansei pmap)are(exn pmap)conditions raised withabort, andpmap-findtakes a default thunk (PM:66, 76, 85, 141); - the errors of
(aux hansei random)are(exn hansei)conditions of kindrandom(OR:91-94); ocaml-random-float!draws three words, as the OCaml that made HANSEI's records did, where OCaml 4.14 draws two (OR:167-172);- a count of
probcc-lazy-permutation/insthat is not a non-negative exact integer raises(exn hansei malformed)(LP:281-286).
Kept as upstream, with a comment: one global RNG state, re-seeded by every random_selector seed; letlazy is not nesting-safe; a memo must be created inside the reification that uses
it; bounded_explore counts every leaf, whatever its value; random_selector raises "can't happen" past the last entry; first_success is breadth-first through append.
Upstream discrepancies: EPP p.4 says biased can be memoized as usual and used at the inner level, which raises with P's memo (only letlazy_nesting works); M1:203-205 re-sums the exact r, not the sample; PD:122 has a syntax typo; TI depends on the unreleased WorldLocalMem; DSL p.5's normalized lawn differs in the last
digits from what today's code computes; and the sampling records of sorted.ml, incite.ml, hmm.ml, colored_balls.ml, bloodtype.ml, blip_naive.ml and blip_ailog.ml do not replay under the three-draw generator, which replays those of samplingM.ml, nested.ml, slazy.ml, parsec.ml and music1a.ml: they come from earlier code or generators, and the
suites check those models exactly, against today's upstream under OCaml, or
statistically in the long suite.
test/probcc/coin-model: pass
Joint distribution of tossing two biased coins, where head has probability 0.59999999999999998 to appear.
(define (test/probcc/coin-model _)
(⊦= '(((V ((x #t) (y #t))) 0.35999999999999999)
((V ((x #t) (y #f))) 0.23999999999999999)
((V ((x #f) (y #t))) 0.23999999999999999)
((V ((x #f) (y #f))) 0.16))
(probcc-normalize
(probcc-reify/exact
(let* ((p 0.59999999999999998)
(x (probcc-coin p))
(y (probcc-coin p)))
`((x ,x) (y ,y)))))))
((eta 0.0) (memory #(12582912 4125944 1048576)) (stdout "") (stderr ""))
test/probcc/coin-model/when: pass
Slightly variation of the previous test, here it has been observed that at least one head appeared.
(define (test/probcc/coin-model/when _)
(⊦= '(((V (#t #t)) 0.42857142857142899)
((V (#t #f)) 0.28571428571428598)
((V (#f #t)) 0.28571428571428598))
(probcc-normalize
(probcc-reify/exact
(let* ((p 0.59999999999999998)
(x (probcc-coin p))
(y (probcc-coin p)))
(probcc-when (or x y) (list x y)))))))
((eta 0.0) (memory #(12582912 4127744 1048576)) (stdout "") (stderr ""))
test/probcc/grass-model: pass
❝
The canonical example is the grass model, with three random variables representing the events of rain, of a switched-on sprinkler and wet grass. The (a priori) probabilities of the first two events are judged to be 30% and 50% correspondingly. Probabilities are non-negative real numbers that may be regarded as weights on non-deterministic choices. Rain almost certainly (90%) wets the grass. The sprinkler also makes the grass wet, in 80% of the cases. The grass may also be wet for some other reason. The modeler gives such an unaccounted event 10% of a chance. This model is often depicted as a directed acyclic graph (DAG) -- so-called Bayesian, or belief network -- with nodes representing random variables and edges conditional dependencies. Associated with each node is a distribution (such as Bernoulli distribution: the flip of a biased coin), or a function that computes a distribution from the node's inputs (such as the noisy disjunction nor). The sort of reasoning we wish to perform on the model is finding out the probability distribution of some of its random variables. For example, we can work out from the model that the probability of the grass being wet is 60.6%. Such reasoning is called probabilistic inference. Often we are interested in the distribution conditioned on the fact that some random variables have been observed to hold a particular value. In our example, having observed that the grass is wet, we want to find out the chance it was raining on that day.
The solution to this problem shows the probability distribution of raining, provided that has been observed a wet grass:
(((V (rain #f)) 0.53152855727963011) ((V (rain #t)) 0.46847144272036984))
(define (test/probcc/grass-model _)
(define result
(probcc-reify/exact
(let* ((rain (probcc-coin 0.29999999999999999))
(sprinkler (probcc-coin 0.5))
(grass-is-wet
(or (and (probcc-coin 0.90000000000000002) rain)
(and (probcc-coin 0.80000000000000004) sprinkler)
(probcc-coin 0.10000000000000001))))
(probcc-when grass-is-wet `(rain ,rain)))))
(⊦= (list (probcc-value 0.32200000000000001 '(rain #f))
(probcc-value 0.2838 '(rain #t)))
result))
((eta 0.0) (memory #(12582912 4129952 1048576)) (stdout "") (stderr ""))
test/probcc/grass-model/complete: pass
If we remove the assumption that has been observed a wet grass, then we have the joint probability distribution of all variables:
(define (test/probcc/grass-model/complete _)
(define result
(probcc-reify/exact
(let* ((rain (probcc-coin 0.29999999999999999))
(sprinkler (probcc-coin 0.5))
(grass-is-wet
(or (and (probcc-coin 0.90000000000000002) rain)
(and (probcc-coin 0.80000000000000004) sprinkler)
(probcc-coin 0.10000000000000001))))
`((rain ,rain) (sprinkler ,sprinkler) (grass-is-wet ,grass-is-wet)))))
(⊦= '(((V ((rain #f) (sprinkler #f) (grass-is-wet #f))) 0.315)
((V ((rain #f) (sprinkler #t) (grass-is-wet #t))) 0.28699999999999998)
((V ((rain #t) (sprinkler #t) (grass-is-wet #t))) 0.14729999999999999)
((V ((rain #t) (sprinkler #f) (grass-is-wet #t))) 0.13650000000000001)
((V ((rain #f) (sprinkler #t) (grass-is-wet #f))) 0.063)
((V ((rain #f) (sprinkler #f) (grass-is-wet #t))) 0.035000000000000003)
((V ((rain #t) (sprinkler #f) (grass-is-wet #f))) 0.0135)
((V ((rain #t) (sprinkler #t) (grass-is-wet #f)))
0.0027000000000000001))
result)
(⊦= (probcc-normalize result) result))
((eta 0.0) (memory #(12582912 4134120 1048576)) (stdout "") (stderr ""))
test/uniform/range: pass
(define (test/uniform/range _)
(⊦= '(((V 1) 1/7)
((V 2) 1/7)
((V 3) 1/7)
((V 4) 1/7)
((V 5) 1/7)
((V 6) 1/7)
((V 7) 1/7))
(sort (probcc-reify/exact (probcc-uniform/range 1 8))
(λ (a b) (< (cadr (car a)) (cadr (car b)))))))
((eta 0.0) (memory #(12582912 4135920 1048576)) (stdout "") (stderr ""))
test/dfs: pass
(define (test/dfs _)
(⊦= '(((V #t) 1/2) ((V #f) 1/2))
(probcc-dfs (probcc-reify/0 (τ (probcc-coin 1/2)))))
(⊦= '(((V (#t #t)) 1/6)
((V (#t #f)) 1/3)
((V (#f #t)) 1/6)
((V (#f #f)) 1/3))
(probcc-dfs
(probcc-reify/0
(τ (let* ((a (probcc-coin 1/2)) (b (probcc-coin 1/3)))
(list a b)))))))
((eta 0.0) (memory #(12582912 4137120 1048576)) (stdout "") (stderr ""))
test/parameterize: pass
(define (test/parameterize _)
(define q (make-parameter 'outer))
(⊦= '(((V (#f inner)) 2/3) ((V (#t inner)) 1/3))
(probcc-reify/exact
(parameterize ((q 'inner)) (let1 (x (probcc-coin 1/3)) (list x (q))))))
(⊦= 'outer (q)))
((eta 0.0) (memory #(12582912 4139400 1048576)) (stdout "") (stderr ""))
test/geometric: pass
(define (test/geometric _)
(define result
((probcc-reify 5) (τ (probcc-geometric 0.84999999999999998 's 'f))))
(define t6 (cadr (car (sixth result))))
(define t7 (cadr (car (seventh result))))
(define t8 (cadr (car (eighth result))))
(⊦= `(((V (s)) 0.84999999999999998)
((V (s f)) 0.1275)
((V (s f f)) 0.019125)
((V (s f f f)) 0.0028687500000000002)
((V (s f f f f)) 0.00043031249999999998)
((C ,t6) 6.4546875000000002e-05)
((C ,t7) 9.6820312500000105e-06)
((C ,t8) 1.70859375e-06))
result))
((eta 0.0) (memory #(12582912 4144376 1048576)) (stdout "") (stderr ""))
test/pmap: pass
The persistent AVL maps of (aux hansei pmap), ExtLib's PMap as augmented by HANSEI with insert_with (cf. PM:106-199): 200 insertions in a scrambled order list their keys in
ascending order, adding to or removing from a map leaves it unchanged, pmap-insert-with stores (f new old), and every node keeps the AVL invariant with upstream's tolerance 2.
(define (test/pmap _)
(define (avl? m)
(letrec ((H (λ (t)
(cond ((null? t) 0)
(else
(let ((hl (H (vector-ref t 0)))
(hr (H (vector-ref t 3))))
(and hl
hr
(<= (abs (- hl hr)) 2)
(= (vector-ref t 4) (add1 (max hl hr)))
(vector-ref t 4))))))))
(and (pmap? m) (H (vector-ref (record->vector m) 2)) #t)))
(define squares (map (λ (k) (cons k (* k k))) (iota 200)))
(define m200
(let loop ((i 0) (m pmap-empty))
(if (= i 200)
m
(let1 (k (remainder (* 37 i) 200))
(loop (add1 i) (pmap-add k (* k k) m))))))
(define m201 (pmap-add 200 'x m200))
(define m199 (pmap-remove 100 m201))
(define m3 (alist->pmap '((5 . b) (1 . a) (9 . c))))
(⊦ equal? squares (pmap->alist m200))
(⊦ equal? #t (avl? m200))
(⊦ equal? #t (avl? m201))
(⊦ equal? #t (avl? m199))
(⊦ equal? 200 (length (pmap->alist m200)))
(⊦ equal? 201 (length (pmap->alist m201)))
(⊦ equal? 200 (length (pmap->alist m199)))
(⊦ equal? #f (pmap-member? 100 m199))
(⊦ equal? #t (pmap-member? 100 m201))
(⊦ equal? #f (pmap-member? 200 m200))
(⊦ equal? 'x (pmap-find 200 m199))
(⊦ equal?
'(n . o)
(pmap-find 'k (pmap-insert-with cons 'k 'n (pmap-add 'k 'o pmap-empty))))
(⊦ equal? 'n (pmap-find 'k (pmap-insert-with cons 'k 'n pmap-empty)))
(⊦ equal?
'((1 . a) ((5 . b) (9 . c)))
(match/first
(pmap-delete-find-min m3)
((,b ,m) (list b (pmap->alist m)))))
(⊦ equal?
'((9 . c) ((1 . a) (5 . b)))
(match/first
(pmap-delete-find-max m3)
((,b ,m) (list b (pmap->alist m)))))
(⊦ equal? '((1 . a) (5 . b) (9 . c)) (pmap->alist m3)))
((eta 0.0) (memory #(12582912 4156320 1048576)) (stdout "") (stderr ""))
test/pmap/api: pass
The rest of the interface of (aux hansei pmap) (cf. PMI): pmap-find with and without a default thunk, the traversals in
ascending key order, the round trip through association lists, a map ordered by a
custom comparison, and the (exn pmap) conditions raised on a
missing key or an empty map.
(define (test/pmap/api _)
(define m (alist->pmap '((b . 2) (a . 1) (c . 3))))
(define (collect-for-each m)
(let1 (acc '())
(pmap-for-each (λ (k v) (set! acc (cons (list k v) acc))) m)
(reverse acc)))
(⊦ equal? 2 (pmap-find 'b m))
(⊦ equal? 'none (pmap-find 'z m (τ 'none)))
(⊦raises (exn pmap) (pmap-find 'z m))
(⊦raises (exn) (with-exception-handler (λ (e) 'x) (τ (pmap-find 'z m))))
(⊦ equal? '(a . 1) (pmap-find-min m))
(⊦raises (exn pmap) (pmap-find-min pmap-empty))
(⊦ equal? #t (pmap-member? 'c m))
(⊦ equal? #f (pmap-member? 'd m))
(⊦ equal? '((a 1) (b 2) (c 3)) (collect-for-each m))
(⊦ equal?
'((a . 10) (b . 20) (c . 30))
(pmap->alist (pmap-map (λ (v) (* 10 v)) m)))
(⊦ equal? '((a a 1) (b b 2) (c c 3)) (pmap->alist (pmap-map/key list m)))
(⊦ equal? '(3 2 1) (pmap-fold cons m '()))
(⊦ equal?
'(c 3 b 2 a 1)
(pmap-foldi (λ (k v acc) (cons k (cons v acc))) m '()))
(⊦ equal?
'((a . 1) (b . 2) (c . 3))
(pmap->alist (alist->pmap (pmap->alist m))))
(⊦ equal? '((a . 9)) (pmap->alist (alist->pmap '((a . 1) (a . 9)))))
(⊦ equal? #t (pmap-empty? pmap-empty))
(⊦ equal? #f (pmap-empty? (pmap-add 1 1 pmap-empty)))
(⊦ equal? #t (pmap? m))
(⊦ equal? #f (pmap? '()))
(⊦ equal?
'(-1 0 1)
(list (fx-compare 1 2) (fx-compare 2 2) (fx-compare 3 2)))
(⊦ equal?
'(-1 0 1)
(list (number-compare 0.5 1)
(number-compare 1/2 0.5)
(number-compare 2 1.5)))
(⊦ equal?
'((3 . c) (2 . b) (1 . a))
(pmap->alist
(alist->pmap '((1 . a) (3 . c) (2 . b)) (λ (a b) (fx-compare b a)))))
(⊦raises (exn pmap) (pmap-delete-find-min pmap-empty))
(⊦raises (exn pmap) (pmap-delete-find-max pmap-empty)))
((eta 0.0) (memory #(12582912 4161288 1048576)) (stdout "") (stderr ""))
test/pmap/compare: pass
pmap-compare stands in for OCaml's polymorphic compare, the default order of PMap.empty (cf. PM:101): types are ranked, numbers compare by value with an exact number
before an equal inexact one, -0. equals 0., symbols compare by name, and procedures cannot be compared, as OCaml's compare: functional value.
(define (test/pmap/compare _)
(⊦ equal?
'(#t 1 1.0 2 #\c "b" a (1) #(0))
(sort (list "b" 'a 2 #t #\c '(1) #(0) 1.0 1)
(λ (a b) (negative? (pmap-compare a b)))))
(⊦ equal? 0 (pmap-compare -0.0 0.0))
(⊦raises (exn pmap) (pmap-compare car cdr)))
((eta 0.0) (memory #(12582912 4163000 1048576)) (stdout "") (stderr ""))
test/explore/persistent: pass
The explorer threads its answer through its loop in a persistent map, as explore threads (ans, susp) (cf. I:16-28):
a thunk it forces may choose in the enclosing model, here the coin of the second
branch, so every world of that model resumes the same loop. With a mutable table
shared by the worlds, the #t of the first world leaked into the
second one, whose inner table read (((V a) 1/2) ((V #f) 1/2) ((V #t) 1/2)).
(define (test/explore/persistent _)
(⊦ equal?
'(((V (((V a) 1/2) ((V #t) 1/2))) 1/2)
((V (((V a) 1/2) ((V #f) 1/2))) 1/2))
(probcc-reify/exact
(probcc-explore
#f
(list (probcc-value 1/2 'a)
(probcc-τ 1/2 (probcc-unit (probcc-coin 1/2))))))))
((eta 0.0) (memory #(12582912 4161976 1048576)) (stdout "") (stderr ""))
test/explore/order: pass
The order of the explorer's output (deviation 5): descending mass, and equal masses
in the order of upstream's PMap.foldi (I:28), descending pmap-compare order of the values, by a stable sort.
(define (test/explore/order _)
(⊦ equal?
'(((V (#t #t)) 1/4) ((V (#t #f)) 1/4) ((V (#f #t)) 1/4) ((V (#f #f)) 1/4))
(probcc-reify/exact
(let* ((a (probcc-coin 1/2)) (b (probcc-coin 1/2))) (list a b)))))
((eta 0.001) (memory #(12582912 4163912 1048576)) (stdout "") (stderr ""))
test/explore/depth: pass
Exploring the lawn model to a depth (cf. DSL p.12 l.539-546, p.16 l.712-719,
I:16-28): its tree starts with the choice of rain; explored to depth 1, the
choices of rain and of the sprinkler are made, and the eight branches of the next
choice, the coin of 0.9, are left unexplored with their masses; from depth 4 on, as
with no bound (#f, OCaml's None), the tree is explored completely. A branch left unexplored keeps its thunk (eq?, I:26), and is not wrapped in a new one.
(define (test/explore/depth _)
(define tree (probcc-reify/0 lawn))
(define depth1 (probcc-explore 1 tree))
(define t0 (τ '()))
(define inner `(((C ,t0) 1/2)))
(⊦= '(0.29999999999999999 0.69999999999999996) (probcc-map cadr tree))
(⊦ equal? 8 (length depth1))
(⊦ equal? #t (every (λ1-match/first (((C ,_) ,_) #t) (else #f)) depth1))
(⊦= '(0.014999999999999999
0.014999999999999999
0.035000000000000003
0.035000000000000003
0.13500000000000001
0.13500000000000001
0.315
0.315)
(sort (probcc-map cadr depth1) <))
(⊦= '(((V #f) 0.32200000000000001) ((V #t) 0.2838)) (probcc-explore 4 tree))
(⊦= '(((V #f) 0.32200000000000001) ((V #t) 0.2838)) (probcc-explore #f tree))
(⊦ eq? t0 (car (cdr (car (car (probcc-explore 0 `(((C ,(τ inner)) 1)))))))))
((eta 0.0) (memory #(12582912 4149696 1048576)) (stdout "") (stderr ""))
test/explore/trees: pass
The five trees of DSL p.8-9 l.364-374 represent the same distribution: the first
is a flat probability table, the others have unexplored branches, the last one a
failure node '(); exploring each of them completely gives the
table.
(define (test/explore/trees _)
(define trees
(list '(((V #t) 0.59999999999999998) ((V #f) 0.40000000000000002))
'(((V #t) 0.40000000000000002)
((V #f) 0.40000000000000002)
((V #t) 0.20000000000000001))
`(((V #t) 0.40000000000000002)
((C ,(τ '(((V #f) 0.40000000000000002)
((V #t) 0.20000000000000001))))
1.0))
`(((V #t) 0.40000000000000002)
((C ,(τ '(((V #f) 0.5) ((V #t) 0.25)))) 0.80000000000000004))
`(((V #t) 0.40000000000000002)
((C ,(τ `(((V #f) 0.5) ((V #t) 0.25) ((C ,(τ '())) 0.25))))
0.80000000000000004))))
(⊦= (make-list
5
'(((V #t) 0.59999999999999998) ((V #f) 0.40000000000000002)))
(map (λ (tree) (probcc-explore #f tree)) trees)))
((eta 0.0) (memory #(12582912 4152344 1048576)) (stdout "") (stderr ""))
test/reflect/inverse: pass
probcc-reflect is the inverse of probcc-reify/0 (cf. P:136-148, DSL p.16 l.694-700, p.17 l.728-733, NIPS l.65-67): reifying the
reflection of a table, or of the lazy tree of the lawn model, gives the
distribution they denote.
(define (test/reflect/inverse _)
(define table
'(((V a) 0.20000000000000001) ((V b) 0.5) ((V c) 0.29999999999999999)))
(define tree (probcc-reify/0 lawn))
(⊦= (probcc-explore #f table) (probcc-reify/exact (probcc-reflect table)))
(⊦= (probcc-explore #f tree) (probcc-reify/exact (probcc-reflect tree)))
(⊦= '(((V b) 0.5) ((V c) 0.29999999999999999) ((V a) 0.20000000000000001))
(probcc-reify/exact (probcc-reflect table))))
((eta 0.0) (memory #(12582912 4151576 1048576)) (stdout "") (stderr ""))
test/normalize: pass
probcc-normalize/total is normalize (cf. I:296-301): the total mass, summed from the left, and the distribution
divided by it; probcc-normalize is the list alone, as the normalize of PE:22-24, here on the extended grass model of
PE:8-28. Exact masses are divided exactly and rounded once, and a total of zero
gives NaN masses, as upstream, instead of raising.
(define (test/normalize _)
(define (grass-cloudy)
(let* ((cloudy (probcc-coin 0.5))
(rain (probcc-coin
(if cloudy 0.80000000000000004 0.20000000000000001)))
(sprinkler (probcc-coin (if cloudy 0.10000000000000001 0.5)))
(wet-roof (and (probcc-coin 0.69999999999999996) rain))
(wet-grass
(or (and (probcc-coin 0.90000000000000002) rain)
(and (probcc-coin 0.90000000000000002) sprinkler))))
(probcc-when wet-grass rain)))
(define posterior (probcc-reify/exact/a grass-cloudy))
(⊦= '(((V #t) 0.45810000000000001) ((V #f) 0.189)) posterior)
(⊦= '(0.64710000000000001
(((V #t) 0.70792767732962447) ((V #f) 0.29207232267037547)))
(probcc-normalize/total posterior))
(⊦= (cadr (probcc-normalize/total posterior)) (probcc-normalize posterior))
(⊦ equal?
`(3/10 (((V a) ,(exact->inexact 1/3)) ((V b) ,(exact->inexact 2/3))))
(probcc-normalize/total '(((V a) 1/10) ((V b) 1/5))))
(⊦ equal? #t (nan? (cadr (car (probcc-normalize '(((V a) 0))))))))
((eta 0.0) (memory #(12582912 4151600 1048576)) (stdout "") (stderr ""))
test/uniform/edges: pass
probcc-uniform (cf. P:156-165) gives the value 0 the mass 1 - (n - 1)/n, the exact 1/n here, so the masses tie and the values come in descending order. Deviation 7: (probcc-uniform 1) is a choice of one branch, where upstream
answers 0 without choosing, a count of zero fails, where upstream raises, and a
count that is not an integer raises (exn hansei) of kind malformed, where OCaml's type rules it out.
(define (test/uniform/edges _)
(⊦ equal?
'(((V 2) 1/3) ((V 1) 1/3) ((V 0) 1/3))
(probcc-reify/exact (probcc-uniform 3)))
(⊦ equal?
#t
(match/first
(probcc-reify/0 (τ (probcc-uniform 1)))
((((C ,_) 1)) #t)
(else #f)))
(⊦ equal? '(((V 0) 1)) (probcc-reify/exact (probcc-uniform 1)))
(⊦ equal? '() (probcc-reify/exact (probcc-uniform 0)))
(⊦ equal? '() (probcc-reify/exact (probcc-uniform/either '())))
(⊦raises (exn hansei) (probcc-reify/exact (probcc-uniform 5/2)))
(⊦raises (exn hansei) (probcc-reify/exact (probcc-uniform 2.5)))
(⊦raises (exn hansei) (probcc-reify/exact (probcc-uniform/range 0 2.5)))
(⊦ equal?
'(((V 1.0) 0.5) ((V 0.0) 0.5))
(probcc-reify/exact (probcc-uniform 2.0))))
((eta 0.001) (memory #(12582912 4156624 1048576)) (stdout "") (stderr ""))
test/map/reentrant: pass
probcc-map is OCaml's List.map, left to right, and a continuation captured inside its procedure may be resumed:
CHICKEN's map builds its result destructively, and (probcc-reify/exact (map (λ (x) (probcc-coin 1/2)) '(1 2))) gives the single entry (#t #t #f #f #t #f) at mass 1, the four worlds sharing one list built in place (plan §5.6).
(define (test/map/reentrant _)
(⊦ equal?
'(((V (#t #t)) 1/4) ((V (#t #f)) 1/4) ((V (#f #t)) 1/4) ((V (#f #f)) 1/4))
(probcc-reify/exact (probcc-map (λ (_) (probcc-coin 1/2)) '(1 2))))
(⊦ equal?
'(((V (#t #t #f #f #t #f)) 1))
(probcc-reify/exact (map (λ (_) (probcc-coin 1/2)) '(1 2)))))
((eta 0.0) (memory #(12582912 4155960 1048576)) (stdout "") (stderr ""))
test/geometric/count: pass
probcc-geometric/count is upstream's geometric (cf. P:181-185), the number of failures before the first success; explored to
depth 4 (cf. X:63-66), the values 0 to 3 come with upstream's masses, bit for bit,
and three branches are left unexplored.
(define (test/geometric/count _)
(define result
(probcc-explore
4
(probcc-reify/0 (τ (probcc-geometric/count 0.84999999999999998)))))
(⊦ equal?
'(((V 0) 0.84999999999999998)
((V 1) 0.1275)
((V 2) 0.019125000000000003)
((V 3) 0.0028687500000000011))
(take result 4))
(⊦ equal?
#t
(every (λ1-match/first (((C ,_) ,_) #t) (else #f)) (drop result 4)))
(⊦ equal?
'(1.1390625000000009e-05 6.4546875000000043e-05 0.00043031250000000025)
(sort (probcc-map cadr (drop result 4)) <)))
((eta 0.0) (memory #(12582912 4161864 1048576)) (stdout "") (stderr ""))
test/geometric/bounded: pass
probcc-geometric/bounded (cf. P:187-206) restricts the geometric
distribution to 0 .. n and normalizes it; with upstream's
operations in upstream's order, the masses are bit for bit those of X:72-77,
BN:189-190 and BL:332-334 (listed here by descending mass), and exact masses stay
exact.
(define (test/geometric/bounded _)
(⊦ equal?
'(((V 0) 0.85000021784575897)
((V 1) 0.12750003267686386)
((V 2) 0.019125004901529581)
((V 3) 0.0028687507352294374)
((V 4) 0.00043031261028441569)
((V 5) 6.454689154266237e-05)
((V 6) 9.6820337313993568e-06)
((V 7) 1.4523050597099037e-06))
(probcc-reify/exact (probcc-geometric/bounded 7 0.84999999999999998)))
(⊦ equal?
'(((V 0) 0.86956521739130432) ((V 1) 0.13043478260869568))
(probcc-reify/exact (probcc-geometric/bounded 1 0.84999999999999998)))
(⊦ equal?
'(((V 0) 0.85043053045604333)
((V 1) 0.12756457956840653)
((V 2) 0.019134686935260981)
((V 3) 0.0028702030402891474))
(probcc-reify/exact (probcc-geometric/bounded 3 0.84999999999999998)))
(⊦ equal?
'(((V 0) 8/15) ((V 1) 4/15) ((V 2) 2/15) ((V 3) 1/15))
(probcc-reify/exact (probcc-geometric/bounded 3 1/2))))
((eta 0.0) (memory #(12582912 4158968 1048576)) (stdout "") (stderr ""))
test/observe: pass
probcc-observe (cf. P:213-214) asserts the evidence that the value
of a thunk passes a test, here the wet grass of the lawn model of X:48-56, whose
masses it gives bit for bit.
(define (test/observe _)
(define (grass-model)
(let* ((rain (probcc-coin 0.29999999999999999))
(sprinkler (probcc-coin 0.5)))
(probcc-observe
(λ (wet) wet)
(τ (or (and (probcc-coin 0.90000000000000002) rain)
(and (probcc-coin 0.80000000000000004) sprinkler)
(probcc-coin 0.10000000000000001))))
rain))
(⊦ equal?
'(((V #f) 0.3219999999999999) ((V #t) 0.2838))
(probcc-reify/exact/a grass-model))
(⊦ equal?
'(((V 2) 1/4) ((V 0) 1/4))
(probcc-reify/exact (probcc-observe even? (τ (probcc-uniform 4))))))
((eta 0.0) (memory #(12582912 4160096 1048576)) (stdout "") (stderr ""))
test/letlazy: pass
probcc-letlazy (cf. P:321-330) evaluates its thunk at most once
per world, when it is first needed: a lazy coin gives the same value twice in a
world; in tests4 (cf. SL:73-93) the expensive function is never
called in the world u = 0, and four times in the world u = 1, split by the coin before x is forced.
(define (test/letlazy _)
(define calls '())
(define (expensive x) (set! calls (cons x calls)) x)
(define (tests4)
(let* ((u (probcc-distribution '((0 1/2) (1 1/2))))
(x (probcc-letlazy
(τ (expensive
(probcc-distribution
`((,(+ u 10) 1/2) (,(+ u 20) 1/2))))))))
(cond ((= u 0) #f)
(else
(let* ((c (probcc-coin 1/2)) (b (x)) (a (x))) (list a b c))))))
(⊦ equal?
'(((V (#f #t #t)) 1/3)
((V (#f #f #f)) 1/3)
((V (#t #t #t)) 1/6)
((V (#t #f #f)) 1/6))
(probcc-reify/exact
(let* ((x (probcc-letlazy (τ (probcc-coin 1/2)))) (c (probcc-coin 1/3)))
(list c (x) (x)))))
(⊦ equal?
'(((V #f) 1/2)
((V (21 21 #t)) 1/8)
((V (21 21 #f)) 1/8)
((V (11 11 #t)) 1/8)
((V (11 11 #f)) 1/8))
(probcc-reify/exact/a tests4))
(⊦ equal? 4 (length calls))
(⊦ equal? #t (every (λ (x) (and (memv x '(11 21)) #t)) calls)))
((eta 0.0) (memory #(12582912 4162520 1048576)) (stdout "") (stderr ""))
test/memo: pass
probcc-memo (cf. P:449-479) makes a stochastic function a proper
function in every world: test2m and test3m (cf. SL:694-729) give the same value for the same argument in a world, also
through another memoized function, and the grass model of PE:33-45, all of whose
variables are memoized, gives upstream's masses bit for bit. The table lives in the
world-local memory, which a reification starts empty (P:106): a memo created outside
the reification it is used in raises (deviation 16).
(define (test/memo _)
(define (f* x) (probcc-distribution `((,x 1/2) (,(add1 x) 1/2))))
(define (test2m)
(let* ((f (probcc-memo f*)) (c (probcc-coin 1/2)))
(if c
(let* ((e (f 2)) (d (f 2)) (b (f 1)) (a (f 1))) (list c a b d e))
(let* ((e (f 1)) (d (f 1)) (b (f 2)) (a (f 2))) (list c a b d e)))))
(define (test3m)
(let* ((f (probcc-memo f*))
(g (probcc-memo
(λ (x)
(let1 (u (f x))
((probcc-distribution
`((,(τ u) 1/2) (,(τ (f x)) 1/2))))))))
(c (probcc-coin 1/2)))
(if c
(let* ((e (g 1)) (d (f 1)) (b (g 1)) (a (f 1))) (list c a b d e))
(let* ((e (f 1)) (d (g 1)) (b (f 1)) (a (g 1))) (list c a b d e)))))
(define (grass-model-memo)
(let* ((cloudy (probcc-memo (λ (_) (probcc-coin 0.5))))
(rain (probcc-memo
(λ (_)
(probcc-coin
(if (cloudy '())
0.80000000000000004
0.20000000000000001)))))
(sprinkler
(probcc-memo
(λ (_)
(probcc-coin (if (cloudy '()) 0.10000000000000001 0.5)))))
(wet-roof
(probcc-memo
(λ (_) (and (probcc-coin 0.69999999999999996) (rain '())))))
(wet-grass
(probcc-memo
(λ (_)
(or (and (probcc-coin 0.90000000000000002) (rain '()))
(and (probcc-coin 0.90000000000000002)
(sprinkler '())))))))
(if (wet-grass '()) (rain '()) (probcc-impossible))))
(⊦ equal?
'(((V (#t 2 2 3 3)) 1/8)
((V (#t 2 2 2 2)) 1/8)
((V (#t 1 1 3 3)) 1/8)
((V (#t 1 1 2 2)) 1/8)
((V (#f 3 3 2 2)) 1/8)
((V (#f 3 3 1 1)) 1/8)
((V (#f 2 2 2 2)) 1/8)
((V (#f 2 2 1 1)) 1/8))
(probcc-reify/exact/a test2m))
(⊦ equal?
'(((V (#t 2 2 2 2)) 1/4)
((V (#t 1 1 1 1)) 1/4)
((V (#f 2 2 2 2)) 1/4)
((V (#f 1 1 1 1)) 1/4))
(probcc-reify/exact/a test3m))
(⊦ equal?
'(((V #t) 0.45810000000000006) ((V #f) 0.18900000000000003))
(probcc-reify/exact/a grass-model-memo))
(⊦raises
(exn hansei)
(probcc-reify/exact
(let1 (f (probcc-memo f*)) (probcc-reify/exact (f 1))))))
((eta 0.0) (memory #(12582912 4168248 1048576)) (stdout "") (stderr ""))
test/nested: pass
Nested inference (cf. N:15-130): the distribution of the distribution of a coin
that is fair or always true, testn11, testn21 and testn21'; with a lazy variable of the outer world, probcc-letlazy is evaluated again inside, where the inner
reification starts with an empty memory, and the answer is wrong, as upstream's testnl10 (N:65-81); probcc-letlazy/nesting asks the parent world instead, and testnl11, testnl12 and testnl13 are right, also through two levels of nesting.
(define (test/nested _)
(⊦ equal?
'(((V (((V #t) 1))) 1/2) ((V (((V #t) 1/2) ((V #f) 1/2))) 1/2))
(probcc-reify/exact
(let1 (coin (probcc-distribution
`((,(τ (probcc-coin 1/2)) 1/2) (,(τ #t) 1/2))))
(probcc-reify/exact/a coin))))
(⊦ equal?
'(((V (((V #t) 1))) 1/2) ((V (((V #t) 1/2) ((V #f) 1/2))) 1/2))
(probcc-reify/exact
(let1 (coin-choice (probcc-coin 1/2))
(probcc-reify/exact (or (probcc-coin 1/2) coin-choice)))))
(⊦ equal?
'(((V (((V #t) 1))) 1/2) ((V (((V #t) 1/2) ((V #f) 1/2))) 1/2))
(probcc-reify/exact
(let1 (coin-choice (probcc-coin 1/2))
(probcc-reify/exact (or coin-choice (probcc-coin 1/2))))))
(⊦ equal?
'(((V (((V #t) 3/4) ((V #f) 1/4))) 1))
(probcc-reify/exact
(let1 (coin-choice (probcc-letlazy (τ (probcc-coin 1/2))))
(probcc-reify/exact (or (coin-choice) (probcc-coin 1/2))))))
(⊦ equal?
'(((V (((V #t) 1))) 1/2) ((V (((V #t) 1/2) ((V #f) 1/2))) 1/2))
(probcc-reify/exact
(let1 (coin-choice (probcc-letlazy/nesting (τ (probcc-coin 1/2))))
(probcc-reify/exact (or (coin-choice) (probcc-coin 1/2))))))
(⊦ equal?
'(((V (((V #t) 1))) 1/2) ((V (((V #t) 1/2) ((V #f) 1/2))) 1/2))
(probcc-reify/exact
(let1 (coin-choice (probcc-letlazy/nesting (τ (probcc-coin 1/2))))
(probcc-reify/exact (or (probcc-coin 1/2) (coin-choice))))))
(⊦ equal?
'(((V (#t #t (((V (#t #t #t)) 1/2) ((V (#f #t #t)) 1/2)))) 1/2)
((V (#f #f (((V (#t #f #f)) 1/2) ((V (#f #f #f)) 1/2)))) 1/2))
(probcc-reify/exact
(let* ((coin-choice (probcc-letlazy/nesting (τ (probcc-coin 1/2))))
(v1 (probcc-reify/exact
(let* ((v (probcc-coin 1/2))
(c (coin-choice))
(b (coin-choice)))
(list v b c))))
(b (coin-choice))
(a (coin-choice)))
(list a b v1))))
(⊦ equal?
'(((V (#t
(((V (#t (((V #t) 1)) #t)) 1/2) ((V (#f (((V #t) 1)) #t)) 1/2))))
1/2)
((V (#f
(((V (#t (((V #f) 3/4) ((V #t) 1/4)) #f)) 1/2)
((V (#f (((V #f) 3/4) ((V #t) 1/4)) #f)) 1/2))))
1/2))
(probcc-reify/exact
(let* ((coin-choice (probcc-letlazy/nesting (τ (probcc-coin 1/2))))
(v (probcc-reify/exact
(let* ((w (probcc-coin 1/2))
(inner (probcc-reify/exact
(or (coin-choice) (probcc-coin 1/4))))
(c (coin-choice)))
(list w inner c))))
(a (coin-choice)))
(list a v)))))
((eta 0.0) (memory #(12582912 4169776 1048576)) (stdout "") (stderr ""))
test/wref: pass
World-local references (cf. TI:121-126, whose WorldLocalMem is not shipped, deviation 17): a reference set in a world holds its value in that
world only.
(define (test/wref _)
(⊦ equal?
'(((V (#t 1)) 1/2) ((V (#f 2)) 1/2))
(probcc-reify/exact
(let* ((r (probcc-wref-new 0)) (c (probcc-coin 1/2)))
(probcc-wref-set! r (if c 1 2))
(list c (probcc-wref-get r)))))
(⊦raises
(exn hansei)
(probcc-reify/exact
(let1 (r (probcc-wref-new 0))
(probcc-reify/exact (probcc-wref-get r))))))
((eta 0.0) (memory #(12582912 4174184 1048576)) (stdout "") (stderr ""))
test/parameterize/memory: pass
The world-local memory is a variable changed with set!, not a parameter, so it is not affected by the dynamic-wind frames that a resume enters again: a lazy variable forced inside a parameterize keeps its value outside, a nested reification inside
a parameterize sees the binding, and a memo called inside a parameterize answers the same value outside.
(define (test/parameterize/memory _)
(define q (make-parameter 'outer))
(⊦ equal?
'(((V (#f #t #t outer)) 1/3)
((V (#f #f #f outer)) 1/3)
((V (#t #t #t inner)) 1/6)
((V (#t #f #f inner)) 1/6))
(probcc-reify/exact
(let* ((x (probcc-letlazy (τ (probcc-coin 1/2))))
(c (probcc-coin 1/3))
(a (if c
(parameterize ((q 'inner)) (list (x) (q)))
(list (x) (q))))
(b (x)))
(list c (car a) b (cadr a)))))
(⊦ equal?
'(((V (#t (((V (#t #t)) 1/2) ((V (#t #f)) 1/2)))) 1/2)
((V (#f (((V (#f #t)) 1/2) ((V (#f #f)) 1/2)))) 1/2))
(probcc-reify/exact
(let1 (c (probcc-coin 1/2))
(parameterize
((q c))
(list c
(probcc-reify/exact
(let1 (d (probcc-coin 1/2)) (list (q) d))))))))
(⊦ equal?
'(((V (#t #t outer)) 1/2) ((V (#f #f outer)) 1/2))
(probcc-reify/exact
(let* ((f (probcc-memo (λ (_) (probcc-coin 1/2))))
(a (parameterize ((q 'inner)) (f 'k)))
(b (f 'k)))
(list a b (q)))))
(⊦ equal? 'outer (q)))
((eta 0.0) (memory #(12582912 4170280 1048576)) (stdout "") (stderr ""))
test/map-dist: pass
probcc-map-dist (cf. BL:456-460) maps the values of a flat table
and collates the equal results, in upstream's order, the results descending.
(define (test/map-dist _)
(⊦ equal?
'(((V 2) 0.25) ((V 1) 0.75))
(probcc-map-dist length '(((V (a)) 0.25) ((V (b)) 0.5) ((V (a b)) 0.25))))
(⊦raises (exn hansei) (probcc-map-dist length `(((C ,(τ '())) 1/2)))))
((eta 0.0) (memory #(12582912 4174184 1048576)) (stdout "") (stderr ""))
test/prob-of: pass
probcc-prob-of and probcc-at-least (cf. N:292-304, PE:49-61): the mass of a value in a table, 0 when it is missing,
values compared as OCaml's = (-0. is 0., NaN is never found), and a branch not yet
explored is skipped (deviation 23); a coin that is fair or
always true gives true with probability at least 0.3 in every world (PE:62-66).
(define (test/prob-of _)
(define table `(((C ,(τ '())) 1/2) ((V #f) 1/4) ((V #t) 1/4)))
(⊦ equal? 1/4 (probcc-prob-of #t table))
(⊦ equal? 0 (probcc-prob-of 'x table))
(⊦ equal? 1/2 (probcc-prob-of -0.0 '(((V 0.0) 1/2))))
(⊦ equal?
1
(probcc-prob-of
-0.0
(probcc-reify/exact (if (probcc-coin 1/2) -0.0 0.0))))
(⊦ equal? 0 (probcc-prob-of +nan.0 '(((V +nan.0) 1/2))))
(⊦ equal? 0 (probcc-prob-of '(1 +nan.0) '(((V (1 +nan.0)) 1/2))))
(⊦ equal? 1/2 (probcc-prob-of '#(a (0.0 b)) '(((V #(a (-0.0 b))) 1/2))))
(⊦ equal? #t (probcc-at-least 1/4 #f table))
(⊦ equal? #f (probcc-at-least 1/3 #f table))
(⊦ equal?
'(((V #t) 1))
(probcc-reify/exact
(let* ((biased (probcc-coin 1/2))
(coin (τ (or (probcc-coin 1/2) biased))))
(probcc-at-least
0.29999999999999999
#t
(probcc-reify/exact/a coin))))))
((eta 0.0) (memory #(12582912 4173848 1048576)) (stdout "") (stderr ""))
test/nested/inference-about-inference: pass
Inference about inference (cf. PE:69-89; EPP p.4 l.345-379): a coin is fair or
always true; its probability of true, estimated by two rejection samples drawn with probcc-dist-selector, that is through the non-determinism of the
outer reification, is at least 0.3 with probability 7/8, whether the bias is chosen
before or lazily, with probcc-letlazy/nesting; the exact
estimate is in test/prob-of.
(define (test/nested/inference-about-inference _)
(define (coin-at-least biased)
(let1 (coin (τ (or (probcc-coin 1/2) (biased))))
(probcc-at-least
0.29999999999999999
#t
(probcc-sample-rejection/a probcc-dist-selector 2 coin))))
(⊦ equal?
'(((V #t) 7/8) ((V #f) 1/8))
(probcc-reify/exact
(let1 (biased (probcc-coin 1/2)) (coin-at-least (τ biased)))))
(⊦ equal?
'(((V #t) 7/8) ((V #f) 1/8))
(probcc-reify/exact
(coin-at-least (probcc-letlazy/nesting (τ (probcc-coin 1/2)))))))
((eta 0.0) (memory #(12582912 4175272 1048576)) (stdout "") (stderr ""))
test/nested/thresholds: pass
The thresholds of N:307-345, with exact tenths from -0.5 to 1.5: the probability that the estimate of a coin, fair or always true, gives true a mass of at least p, by exact inference (1, 1/2, 0) and by two rejection samples (1, 7/8, 5/8, 0), the bias chosen before or lazily.
(define (test/nested/thresholds _)
(define ((test-nest p inferencer))
(let1 (x (probcc-coin 1/2))
(probcc-at-least p #t (inferencer (τ (or (probcc-coin 1/2) x))))))
(define ((test-nest-lazy p inferencer))
(let1 (x (probcc-letlazy/nesting (τ (probcc-coin 1/2))))
(probcc-at-least p #t (inferencer (τ (or (probcc-coin 1/2) (x)))))))
(define (sampler m) (probcc-sample-rejection/a probcc-dist-selector 2 m))
(define (true-mass model) (probcc-prob-of #t (probcc-reify/exact/a model)))
(define thresholds (probcc-map (λ (n) (/ n 10)) (iota 21 -5)))
(⊦ equal?
(probcc-map
(λ (p)
(let ((expected-exact (cond ((<= p 1/2) 1) ((<= p 1) 1/2) (else 0)))
(expected-sample
(cond ((<= p 0) 1)
((<= p 1/2) 7/8)
((<= p 1) 5/8)
(else 0))))
(list p
expected-exact
expected-exact
expected-sample
expected-sample)))
thresholds)
(probcc-map
(λ (p)
(let* ((a (true-mass (test-nest p probcc-reify/exact/a)))
(b (true-mass (test-nest-lazy p probcc-reify/exact/a)))
(c (true-mass (test-nest p sampler)))
(d (true-mass (test-nest-lazy p sampler))))
(list p a b c d)))
thresholds)))
((eta 0.001) (memory #(12582912 4185680 1048576)) (stdout "") (stderr ""))
test/dist-selector: pass
probcc-dist-selector (cf. P:238-241) chooses through the
non-determinism of the enclosing reification, so a nested rejection sampler of n
samples answers the distribution of its estimates (N:145-160): one sample is the
identity, two give three tables; it answers the chosen entry and the total mass of
the entries, here as a distribution.
(define (test/dist-selector _)
(define (tflip2-noshared) (and (probcc-coin 1/2) (probcc-coin 1/2)))
(⊦ equal?
'(((V (((V #f) 1))) 3/4) ((V (((V #t) 1))) 1/4))
(probcc-reify/exact
(probcc-sample-rejection/a probcc-dist-selector 1 tflip2-noshared)))
(⊦ equal?
'(((V (((V #f) 1))) 9/16)
((V (((V #t) 1/2) ((V #f) 1/2))) 3/8)
((V (((V #t) 1))) 1/16))
(probcc-reify/exact
(probcc-sample-rejection/a probcc-dist-selector 2 tflip2-noshared)))
(⊦ equal?
'(((V (#t 4)) 3/4) ((V (#f 4)) 1/4))
(probcc-reify/exact (probcc-dist-selector '((#t 3) (#f 1))))))
((eta 0.0) (memory #(12582912 4179256 1048576)) (stdout "") (stderr ""))
test/first-success: pass
probcc-first-success (cf. I:55-64) opens the branches of a tree
breadth first, the masses multiplied along the way, until its first entry is a leaf,
and answers that tree; a tree without leaves gives the empty tree.
(define (test/first-success _)
(⊦ equal?
'(((V a) 1/2) ((V b) 1/2))
(probcc-first-success `(((C ,(τ '(((V b) 1)))) 1/2) ((V a) 1/2))))
(⊦ equal? '() (probcc-first-success '()))
(⊦ equal?
'()
(probcc-first-success
(probcc-reify/0 (τ (probcc-coin 1/2) (probcc-impossible)))))
(⊦ equal?
'(((V #t) 1/6) ((V #f) 1/2))
(probcc-first-success
(probcc-reify/0
(τ (if (probcc-coin 1/3) (probcc-impossible) (probcc-coin 1/4)))))))
((eta 0.0) (memory #(12582912 4180864 1048576)) (stdout "") (stderr ""))
test/soft-cut: pass
probcc-soft-cut p q r (cf. LV:386-391, PS:156-161) is Prolog's
soft cut p *-> q ; r: q gets every value of p when p has one, r runs otherwise.
(define (test/soft-cut _)
(⊦ equal?
'(((V (q #t)) 1/2) ((V (q #f)) 1/2))
(probcc-reify/exact
(probcc-soft-cut (τ (probcc-coin 1/2)) (λ (x) (list 'q x)) (τ 'r))))
(⊦ equal?
'(((V r) 1))
(probcc-reify/exact
(probcc-soft-cut (τ (probcc-impossible)) (λ (x) (list 'q x)) (τ 'r)))))
((eta 0.0) (memory #(12582912 4182912 1048576)) (stdout "") (stderr ""))
test/explore/shallow: pass
probcc-explore/shallow (cf. I:35-53) explores a tree to a depth
without flattening it: exploring the result fully gives the exact distribution of
the lawn model, of geometric_bounded 7 0.85 and of dcoin_and 5, whatever the depth, up to rounding; the entries are
upstream's, in its order, which the importance sampler's selections follow.
(define (test/explore/shallow _)
(define (dcoin-and-5) (dcoin-and 5))
(define (geometric-7) (probcc-geometric/bounded 7 0.84999999999999998))
(define (shallow-then-explore k model)
(probcc-explore #f (probcc-explore/shallow k (probcc-reify/0 model))))
(⊦ equal?
'((#t #t #t) (#t #t #t) (#t #t #t))
(probcc-map
(λ (model)
(let1 (exact (probcc-explore #f (probcc-reify/0 model)))
(probcc-map
(λ (k)
(let1 (got (shallow-then-explore k model))
(and (equal? (map car exact) (map car got))
(every (λ (a b)
(< (abs (- (cadr a) (cadr b)))
1.0000000000000001e-15))
exact
got))))
'(0 1 3))))
(list lawn geometric-7 dcoin-and-5)))
(define t (τ '(((V b) 1/4) ((V c) 1/4))))
(⊦ equal?
`(((V b) 1/2) ((V a) 1/4) ((C ,t) 1/2))
(probcc-explore/shallow
0
`(((V a) 1/8) ((V b) 1/2) ((C ,t) 1/2) ((V a) 1/8))))
(define u (τ `(((V c) 1/4) ((V b) 1/4) ((C ,t) 1/2))))
(⊦ equal?
`(((C ,t) 1))
(match/first
(probcc-explore/shallow 1 `(((C ,u) 1)))
((((V c) 1/4) ((V b) 1/4) ((C ,t*) 1/2)) (t*)))))
((eta 0.001) (memory #(12582912 4185504 1048576)) (stdout "") (stderr ""))
test/explore/bounded: pass
probcc-explore/bounded (cf. I:71-109) bounds the probability of
evidence, bit for bit as SM:80-82 and SM:469-480: geometric_bounded is solved with a queue of one branch; for dcoin_and 10 the bounds close as the queue grows, and the lower bound is the exact 0 while no
leaf is reached (deviation 26).
(define (test/explore/bounded _)
(define (u) (probcc-reify/0 (τ (unless (dcoin-and 10) (probcc-impossible)))))
(⊦ equal?
'(0.99949400615938178 0.99949400615938178)
(probcc-explore/bounded
1
(probcc-reify/0
(τ (when (> (probcc-geometric/bounded 7 0.84999999999999998) 3)
(probcc-impossible))))))
(⊦ equal? '(0 1.0) (probcc-explore/bounded 1 (u)))
(⊦ equal? '(0 0.049999999999999989) (probcc-explore/bounded 3 (u)))
(⊦ equal?
'(9.7656249999999776e-14 9.7656249999999776e-14)
(probcc-explore/bounded 5 (u))))
((eta 0.001) (memory #(12582912 4186224 1048576)) (stdout "") (stderr ""))
test/sample/look-ahead-exact: pass
With lazy flips (DSL p.20 l.866-909), each flip is observed right after it is made, so the look-ahead of importance sampling finds the probability 2^-20 of 20 true flips exactly, with one sample; eager flips are produced before any is observed, and exact inference explores the whole tree.
(define (test/sample/look-ahead-exact _)
(define (flips/lazy p n)
(if (zero? n)
'()
(let* ((x (probcc-letlazy (τ (probcc-coin p))))
(xs (probcc-letlazy (τ (flips/lazy p (sub1 n))))))
(cons x xs))))
(define (trues/lazy n xs)
(match/first
xs
(() (zero? n))
((,x unquote xs*) (and (x) (trues/lazy (sub1 n) (xs*))))))
(define (flips p n)
(if (zero? n)
'()
(let* ((x (probcc-coin p)) (xs (flips p (sub1 n)))) (cons x xs))))
(⊦ equal?
'(((V ok) 9.5367431640625e-07))
(probcc-sample-importance
(probcc-random-selector 17)
1
(probcc-when (trues/lazy 20 (flips/lazy 0.5 20)) 'ok)))
(⊦ equal?
'(((V ok) 1/4096))
(probcc-reify/exact (probcc-when (every identity (flips 1/2 12)) 'ok))))
((eta 0.016) (memory #(12582912 4189344 1048576)) (stdout "") (stderr ""))
test/sample/pre-explored: pass
The importance sampler first explores the tree to depth 3 and looks one step ahead (cf. I:188-197): these models of SM:13-206 are solved before any sample is drawn, so their results are exact, bit for bit as upstream's, and the selector, which here raises, is never called.
(define (test/sample/pre-explored _)
(define (no-selector choices) (error "a pre-explored tree draws" choices))
(define (tflip2-shared) (let1 (v (probcc-coin 0.5)) (and v v)))
(define (tflip2-noshared) (and (probcc-coin 0.5) (probcc-coin 0.5)))
(define (talarm)
(let* ((earthquake (probcc-coin 0.01))
(burglary (probcc-coin 0.10000000000000001)))
(if earthquake
(if burglary
(probcc-coin 0.98999999999999999)
(probcc-coin 0.20000000000000001))
(if burglary (probcc-coin 0.97999999999999998) (probcc-coin 0.01)))))
(define (tiblt2)
(let* ((x (probcc-observe
(λ (r) (eq? (car r) 'a))
(τ (probcc-distribution
'(((a b) 0.01)
((a c) 0.02)
((d e) 0.96999999999999997))))))
(y (case (cadr x)
((b)
(probcc-distribution
'((#t 0.90000000000000002) (#f 0.10000000000000001))))
((c)
(probcc-distribution
'((#t 0.59999999999999998) (#f 0.40000000000000002))))
(else
(probcc-distribution
'((#t 0.20000000000000001) (#f 0.80000000000000004)))))))
y))
(define (tfooibl3)
(probcc-observe
(λ (x) (eq? x #t))
(τ (eq? 'a
(if (probcc-distribution
'((#t 0.10000000000000001) (#f 0.90000000000000002)))
(probcc-distribution
'((a 0.69999999999999996) (b 0.29999999999999999)))
(probcc-distribution
'((a 0.20000000000000001) (b 0.80000000000000004))))))))
(define (tmusictr42)
(cadr (probcc-observe
(λ (r) (eq? (car r) 'a))
(τ (if (probcc-coin 0.01)
(let1 (q (probcc-coin 0.29999999999999999)) (list 'a q))
(list 'b #t))))))
(⊦ equal?
'(((V #t) 0.5) ((V #f) 0.5))
(probcc-sample-importance/a no-selector 100 tflip2-shared))
(⊦ equal?
'(((V #f) 0.75) ((V #t) 0.25))
(probcc-sample-importance/a no-selector 100 tflip2-noshared))
(⊦ equal?
'(((V #f) 0.89128000000000007) ((V #t) 0.10872000000000001))
(probcc-sample-importance/a no-selector 100 talarm))
(⊦ equal?
'(((V #t) 0.021000000000000001) ((V #f) 0.0090000000000000011))
(probcc-sample-importance/a no-selector 100 tiblt2))
(⊦ equal?
'(((V #t) 0.25))
(probcc-sample-importance/a no-selector 100 tfooibl3))
(⊦ equal?
'(((V #f) 0.0069999999999999993) ((V #t) 0.0030000000000000001))
(probcc-sample-importance/a no-selector 1 tmusictr42)))
((eta 0.001) (memory #(12582912 4184912 1048576)) (stdout "") (stderr ""))
test/sample/dist-selector-importance: pass
Bounded-rational theory of mind (cf. IN:344-426): agent S chooses what to say,
knowing that agent A estimates the state of the world by importance sampling of two
samples, through probcc-dist-selector; saying only not all gives the maximal expected utility, -0.
(define ((agent-A-model form))
(let* ((s0 (probcc-letlazy (τ (probcc-coin 1/2))))
(s1 (probcc-letlazy (τ (probcc-coin 1/2))))
(sc (probcc-letlazy (τ (probcc-coin 1/2))))
(state (λ (s) (case s ((S_b0) (s0)) ((S_b1) (s1)) ((S_conj) (sc))))))
(probcc-when
(language form state)
(let* ((action (probcc-distribution '((0 1) (1 1) (2 1) (3 1))))
(u (utility state action)))
(list u action)))))
(define ((agent-S-model inferencer))
(let* ((formvs
(probcc-distribution (probcc-map (λ (v) (list v 1)) all-forms)))
(form (λ (fv) (and (memq fv formvs) #t)))
(action-A
(match/first
(argmax (inferencer (agent-A-model form)))
(((,_ ,x)) x)
(() (probcc-uniform 4))
(,lst
(let1 (p (/ 1 (length lst)))
(probcc-distribution
(probcc-map
(λ1-match/first ((,_ ,x) (list x p)))
lst))))))
(u (utility (λ (s) (eq? s 'S_b1)) action-A)))
(list u formvs)))
(define (test/sample/dist-selector-importance _)
(define (inferencer m) (probcc-sample-importance/a probcc-dist-selector 2 m))
(⊦ equal?
'((#t (F_notall)))
(probcc-map
(λ1-match/first ((,u ,a) (list (= u 0) a)))
(argmax (probcc-reify/exact/a (agent-S-model inferencer))))))
((eta 0.002) (memory #(12582912 4197240 1048576)) (stdout "") (stderr ""))
test/statistics: pass
probcc-statistics (cf. I:279-294) gives the mean and the
standard deviation of the mass of each value over the samplers of a range of seeds,
by descending mean (deviation 21); the deviation is a flonum, NaN when rounding makes
the variance negative, as upstream (deviation 27).
(define (test/statistics _)
(⊦ equal?
'((a 1/5 0.081649658092772609))
(probcc-statistics 1 3 (λ (s) (list (probcc-value (/ s 10) 'a)))))
(⊦ equal?
#t
(flonum?
(caddr (car (probcc-statistics
1
3
(λ (s) (list (probcc-value 0.10000000000000001 'a))))))))
(⊦ equal?
'((b 3/4 0.0) (a 1/4 0.0))
(probcc-statistics 1 2 (λ (s) '(((V a) 1/4) ((V b) 3/4)))))
(⊦raises (exn hansei) (probcc-statistics 1 1 (λ (s) `(((C ,(τ '())) 1))))))
((eta 0.0) (memory #(12582912 4201320 1048576)) (stdout "") (stderr ""))
test/verbose: pass
The samplers print their done lines (I:131, 186) only when probcc-verbose is true (deviation 14).
(define (test/verbose _)
(⊦ equal?
""
(with-output-to-string
(τ (probcc-sample-rejection
(probcc-random-selector 1)
10
(probcc-coin 1/2)))))
(⊦ equal?
"rejection_sample: done 10 worlds\n"
(with-output-to-string
(τ (parameterize
((probcc-verbose #t))
(probcc-sample-rejection
(probcc-random-selector 1)
10
(probcc-coin 1/2))))))
(⊦ equal?
"sample_importance: done 10 worlds\n"
(with-output-to-string
(τ (parameterize
((probcc-verbose #t))
(probcc-sample-importance
(probcc-random-selector 1)
10
(dcoin-and 10)))))))
((eta 0.001) (memory #(12582912 4193128 1048576)) (stdout "") (stderr ""))
test/variable-elimination: pass
Variable elimination (DSL p.17 l.751-771; X:305-329): the xor of 10 flips explores the 1024 leaves of a full binary tree, while reifying and reflecting the recursive call leaves 4 at the top; the bucket elimination of DSL p.18 l.776-783 memoizes the reified calls.
(define (test/variable-elimination _)
(define (flips-xor n)
(if (= n 1)
(probcc-coin 1/2)
(let* ((r (flips-xor (sub1 n))) (f (probcc-coin 1/2))) (not (eq? f r)))))
(define (flips-xor* n)
(if (= n 1)
(probcc-coin 1/2)
(let* ((r (probcc-reflect (probcc-reify/exact (flips-xor* (sub1 n)))))
(f (probcc-coin 1/2)))
(not (eq? f r)))))
(define xor-bucket
(λ-probcc-bucket
(n)
(if (= n 1)
(probcc-coin 1/2)
(let* ((r (xor-bucket (sub1 n))) (f (probcc-coin 1/2)))
(not (eq? f r))))))
(⊦ equal? '(((V #t) 1/2) ((V #f) 1/2)) (probcc-reify/exact (flips-xor 10)))
(⊦ equal? 1024 (probcc-leaves (probcc-reify/0 (τ (flips-xor 10)))))
(⊦ equal? '(((V #t) 1/2) ((V #f) 1/2)) (probcc-reify/exact (flips-xor* 10)))
(⊦ equal? 4 (probcc-leaves (probcc-reify/0 (τ (flips-xor* 10)))))
(⊦ equal? '(((V #t) 1/2) ((V #f) 1/2)) (probcc-reify/exact (xor-bucket 10))))
((eta 0.0050000000000000001)
(memory #(12582912 4201000 1048576))
(stdout "")
(stderr ""))
test/drunk-coin: pass
The drunk coin of SM:224-238 and DSL p.19 l.814-817: the exact distribution of the conjunction of 10 tosses, bit for bit as upstream.
(define (drunk-coin)
(let* ((x (probcc-coin 0.5)) (lost (probcc-coin 0.90000000000000002)))
(if lost (probcc-impossible) x)))
(define (dcoin-and n)
(if (= n 1) (drunk-coin) (and (drunk-coin) (dcoin-and (sub1 n)))))
(define (test/drunk-coin _)
(⊦ equal?
'(((V #f) 0.052631578947363269) ((V #t) 9.7656249999999776e-14))
(probcc-reify/exact (dcoin-and 10))))
((eta 0.0) (memory #(12582912 4198176 1048576)) (stdout "") (stderr ""))
test/reify/part: pass
probcc-reify/part (cf. LV:68-70) runs a model to a depth and
keeps its values, as (v p) pairs: the boolean lists of RT p.5 l.340-365, of which
depth 3 reaches the empty list and the lists of one element.
(define (test/reify/part _)
(define (a-list)
(if (probcc-coin 1/2)
'()
(let* ((t (a-list)) (h (probcc-coin 1/2))) (cons h t))))
(⊦ equal? '((() 1/2) ((#t) 1/8) ((#f) 1/8)) (probcc-reify/part 3 a-list))
(⊦ equal? '((#t 1/2) (#f 1/2)) (probcc-reify/part #f (τ (probcc-coin 1/2)))))
((eta 0.0) (memory #(12582912 4199800 1048576)) (stdout "") (stderr ""))
test/bind: pass
probcc-bind (cf. DSL p.10 l.415-425 pv_bind) with probcc-unit is the search-tree monad, whose laws hold up
to the equivalence of trees explored fully.
(define (test/bind _)
(define (f x) (probcc-reify/0 (τ (let1 (c (probcc-coin 1/3)) (list x c)))))
(define (m1) (probcc-coin 1/2))
(define (m2) (let* ((a (probcc-coin 1/2)) (b (probcc-coin 1/4))) (and a b)))
(⊦ equal?
(probcc-reify/exact (probcc-reflect (f (m1))))
(probcc-explore #f (probcc-bind (probcc-reify/0 m1) f)))
(⊦ equal?
(probcc-reify/exact (probcc-reflect (f (m2))))
(probcc-explore #f (probcc-bind (probcc-reify/0 m2) f)))
(⊦ equal?
(probcc-explore #f (f 'x))
(probcc-explore #f (probcc-bind (probcc-unit 'x) f))))
((eta 0.0) (memory #(12582912 4201448 1048576)) (stdout "") (stderr ""))
test/timeit-runner-memory: pass
probcc-timeit (cf. I:303-308) prints the processor time a thunk
takes and answers its value; probcc-sample-runner/times (cf. I:142-146) applies a step n times; probcc-mem-ref (cf. P:52-54) calls its default thunk on a missing location, and raises without
one.
(define (test/timeit-runner-memory _)
(define v #f)
(define out (with-output-to-string (τ (set! v (probcc-timeit (τ 42))))))
(⊦ equal? 42 v)
(⊦ equal? "\nTime spent: " (substring out 0 13))
(⊦ equal? " sec\n" (substring out (- (string-length out) 5)))
(⊦ equal? '(3 3) ((probcc-sample-runner/times 3) 0 add1))
(⊦ equal? 'd (probcc-mem-ref (probcc-mem-new-loc) probcc-mem-empty (τ 'd)))
(⊦raises
(exn hansei)
(probcc-mem-ref (probcc-mem-new-loc) probcc-mem-empty)))
((eta 0.0) (memory #(12582912 4197824 1048576)) (stdout "") (stderr ""))
(aux hansei random): replaying the samples of HANSEI bit for bit
Replaying upstream's records
HANSEI's samplers draw from OCaml's Random module, through random_selector seed, which re-initialises its global state (P:221-222). (aux hansei random) transcribes that generator, the lagged Fibonacci
generator of 55 words seeded with MD5 RFC 1321, with the xor bits and the three-draw Random.float of the OCaml releases HANSEI's records were made with (OCaml 4.14 draws two words), the
only variant that replays the records of samplingM.ml and nested.ml bit for bit; see the plan, D6 and D20, and Embedded probabilistic domain-specific language HANSEI.
Why three draws: OCaml 4.14's Random.float (OR:167-172) draws two
30-bit words, (r1 /. 2^30 +. r2) /. 2^30, and with it the rejection
sampling of dcoin_and 10 over 10000 samples, seed 17, gives 0.054 where SM records 0.0514; three words, ((r0 /. 2^30 +. r1) /. 2^30 +. r2) /. 2^30, replay every assertion
of samplingM.ml, nested.ml, slazy.ml, parsec.ml and music1a.ml under an OCaml build of HANSEI patched so (deviation 32).
The global state is upstream's too: a selector made with a seed re-initialises the one
state every selector draws from, which N:30-31 needs (D6).
The expected values come from three sources: upstream's own assertions; runs of
upstream HANSEI under OCaml 4.14.2 with delimcc, two-draw and three-draw; and a Python
transcription of the samplers over OCaml's generator, written while planning the port
and kept in test/hansei-python/ (its README says how to run it under python3 -I; it is not part of make test).
test/md5: pass
md5-digest is OCaml's Digest.string, which Random.full_init uses to seed the state: the test suite of RFC 1321 (A.5), and messages of
55 to 120 bytes around the padding boundaries, whose digests come from OCaml 4.14.2.
(define (test/md5 _)
(⊦ equal? "d41d8cd98f00b204e9800998ecf8427e" (md5-hex ""))
(⊦ equal? "0cc175b9c0f1b6a831c399e269772661" (md5-hex "a"))
(⊦ equal? "900150983cd24fb0d6963f7d28e17f72" (md5-hex "abc"))
(⊦ equal? "f96b697d7cb7938d525a2f31aaf161d0" (md5-hex "message digest"))
(⊦ equal?
"c3fcd3d76192e4007dfb496cca67e13b"
(md5-hex "abcdefghijklmnopqrstuvwxyz"))
(⊦ equal?
"d174ab98d277d9f5a5611c2c9f419d9f"
(md5-hex
"ABCDEFGHIJKLMNOPQRSTUVWXYZabcdefghijklmnopqrstuvwxyz0123456789"))
(⊦ equal?
"57edf4a22be3c955ac49da2e2107b67a"
(md5-hex (apply string-append (make-list 8 "1234567890"))))
(⊦ equal? "04364420e25c512fd958a70738aa8f72" (md5-hex (make-string 55 #\x)))
(⊦ equal? "668a72d5ba17f08e62dabcafad6db14b" (md5-hex (make-string 56 #\x)))
(⊦ equal? "7dc2ca208106a2f703567bdff99d8981" (md5-hex (make-string 63 #\x)))
(⊦ equal? "c1bb4f81d892b2d57947682aeb252456" (md5-hex (make-string 64 #\x)))
(⊦ equal? "b2ea3db8221fa9ae01ee0adf42970ba0" (md5-hex (make-string 119 #\y)))
(⊦ equal? "723628a6315542e744e4bc6e2f9c84d2" (md5-hex (make-string 120 #\y)))
(⊦ equal?
#u8(144 1 80 152 60 210 79 176 214 150 63 125 40 225 127 114)
(md5-digest #u8(97 98 99))))
((eta 0.0) (memory #(12582912 2650224 1048576)) (stdout "") (stderr ""))
test/ocaml/init: pass
make-ocaml-random-state seed is OCaml's Random.State.make [| seed |] (OR:40-64): the first and last words
after seeding, as OCaml 4.14.2 prints them through Obj.magic (and as the Python transcription of the plan, §10.3, computes them). Without a seed
it is OCaml's default state, the one Random.bits () draws from in a
program that never seeds; a seed that is not an exact integer, or a vector of them,
raises (exn hansei) of kind random.
(define (test/ocaml/init _)
(define (words s) (vector-ref (record->vector s) 1))
(define (idx s) (vector-ref (record->vector s) 2))
(define (first-words s)
(list (vector-ref (words s) 0)
(vector-ref (words s) 1)
(vector-ref (words s) 2)
(vector-ref (words s) 3)
(vector-ref (words s) 4)))
(define s1 (make-ocaml-random-state 1))
(⊦ equal?
'(613788242 197116386 1014951421 468579290 805446940)
(first-words s1))
(⊦ equal?
'(145365742 323272762 161958640 284739705 672753100)
(list (vector-ref (words s1) 50)
(vector-ref (words s1) 51)
(vector-ref (words s1) 52)
(vector-ref (words s1) 53)
(vector-ref (words s1) 54)))
(⊦ equal? 0 (idx s1))
(⊦ equal?
'(43684900 826807282 853359337 631525392 682719058)
(first-words (make-ocaml-random-state 17)))
(⊦ equal?
'(76722835 608589901 287828993 653199777 415478086)
(first-words (make-ocaml-random-state 0)))
(⊦ equal?
'(640528891 982617058 790696909 229908377 1042402058)
(first-words (make-ocaml-random-state -5)))
(⊦ equal?
'(972831248 467146799 206477288 695949822 606390225)
(first-words (make-ocaml-random-state 123456789)))
(⊦ equal?
'(987910699 495797812 364182224 414272206 318284740)
(first-words (make-ocaml-random-state)))
(⊦ equal? #t (ocaml-random-state? ocaml-random-global))
(⊦ equal? #f (ocaml-random-state? (words s1)))
(⊦raises (exn hansei) (make-ocaml-random-state 'seed))
(⊦raises (exn hansei) (make-ocaml-random-state 1.0))
(⊦raises (exn hansei) (make-ocaml-random-state #(1 a))))
((eta 0.002) (memory #(12582912 2654472 1048576)) (stdout "") (stderr ""))
test/ocaml/bits: pass
ocaml-random-bits! is OCaml's Random.State.bits (OR:76-83), the xor variant of the lagged Fibonacci step: the first five draws after Random.init for several seeds, after Random.full_init, and from the default state, as OCaml 4.14.2 draws
them.
(define (test/ocaml/bits _)
(define (draws n f)
(if (zero? n) '() (let1 (x (f)) (cons x (draws (sub1 n) f)))))
(define (bits s) (draws 5 (τ (ocaml-random-bits! s))))
(⊦ equal?
'(111472438 155495486 1029890784 56915826 660037571)
(bits (make-ocaml-random-state 1)))
(⊦ equal?
'(230024973 197978033 495359752 1036278348 906260834)
(bits (make-ocaml-random-state 17)))
(⊦ equal?
'(218040752 178628190 911923154 33712936 1041571162)
(bits (make-ocaml-random-state 0)))
(⊦ equal?
'(638983422 1072006411 37028206 270829641 521483706)
(bits (make-ocaml-random-state -5)))
(⊦ equal?
'(1047213900 619322777 232370739 831219731 841638859)
(bits (make-ocaml-random-state 123456789)))
(⊦ equal?
'(89809344 994326685 290180182 366831641 332180439)
(bits (make-ocaml-random-state)))
(⊦ equal?
'(218040752 178628190 911923154 33712936 1041571162)
(bits (make-ocaml-random-state #())))
(⊦ equal?
'(311998834 678888426 938209723 870508417 248299610)
(bits (make-ocaml-random-state #(1 2 3))))
(⊦ equal?
'(586659661 422753529 1044869830)
(let1 (s (make-ocaml-random-state 17))
(draws 997 (τ (ocaml-random-bits! s)))
(draws 3 (τ (ocaml-random-bits! s))))))
((eta 0.002) (memory #(12582912 2689088 1048576)) (stdout "") (stderr ""))
test/ocaml/float: pass
The doubles of a fresh source are OCaml's three-draw rawfloat, ((r0 / 2^30 + r1) / 2^30 + r2) / 2^30, bit for bit (OCaml
4.14.2 running the three-draw formula over Random.bits): the one
variant that replays HANSEI's records (plan D6, D20).
(define (test/ocaml/float _)
(define (doubles seed)
(let* ((src (make-ocaml-random-source seed))
(a (rng-double! src 0.0 1.0))
(b (rng-double! src 0.0 1.0))
(c (rng-double! src 0.0 1.0)))
(list a b c)))
(⊦ equal?
'(0.95916053666250445 0.047719400948572799 0.64524492170467784)
(doubles 1))
(⊦ equal?
'(0.46133971976524352 0.18637992054598515 0.699662755317079)
(doubles 17))
(⊦ equal?
'(0.84929461978968279 0.62951315564107035 0.4331830841471257)
(doubles 0))
(⊦ equal?
'(0.03448520507512965 0.97297789481065189 0.34460494547462833)
(doubles -5))
(⊦ equal?
'(0.00046133971976524351 0.055913976163795544 1.7491568882926976)
(let* ((s (make-ocaml-random-state 17))
(a (ocaml-random-float! s 0.001))
(b (ocaml-random-float! s 0.29999999999999999))
(c (ocaml-random-float! s 2.5)))
(list a b c)))
(⊦ equal?
0.95916053666250445
(rng-double! (make-ocaml-random-source 1) 0 1)))
((eta 0.001) (memory #(12582912 2658416 1048576)) (stdout "") (stderr ""))
test/ocaml/source: pass
make-ocaml-random-source makes an OCaml state an rng-source of (aux hakaru10 random): rng-source-copy replays it, rng-integer! is lo + Random.int (hi -
lo + 1) (OR:86-94, values from OCaml 4.14.2), and a source over a given state, such
as ocaml-random-global, draws from that state, not from a copy.
(define (test/ocaml/source _)
(define (draws n f)
(if (zero? n) '() (let1 (x (f)) (cons x (draws (sub1 n) f)))))
(define (ints seed n)
(let1 (s (make-ocaml-random-state seed))
(draws 5 (τ (ocaml-random-int! s n)))))
(define src (make-ocaml-random-source 17))
(define copy (rng-source-copy src))
(define used (make-ocaml-random-state 17))
(define global-copy (ocaml-random-state-copy ocaml-random-global))
(define global-source (make-ocaml-random-source ocaml-random-global))
(⊦ equal?
(draws 1000 (τ (rng-double! src 0.0 1.0)))
(draws 1000 (τ (rng-double! copy 0.0 1.0))))
(⊦ equal?
(draws 1000 (τ (rng-integer! src 0 99)))
(draws 1000 (τ (rng-integer! copy 0 99))))
(⊦ equal?
'(-3 -2 -1 0 1 2 3 4)
(sort (delete-duplicates (draws 1000 (τ (rng-integer! src 4 -3)))) <))
(⊦ equal?
'(9 7 5 7 2)
(let1 (s (make-ocaml-random-source 1))
(draws 5 (τ (rng-integer! s 1 10)))))
(⊦ equal? '(8 6 4 6 1) (ints 1 10))
(⊦ equal? '(438 486 784 826 571) (ints 1 1000))
(⊦ equal? '(3 3 2 8 4) (ints 17 10))
(⊦ equal? '(973 33 752 348 834) (ints 17 1000))
(⊦ equal? '(0 6 6 2 4) (ints 1 7))
(⊦ equal?
'(230024973 197978033 495359752 1036278348 906260834)
(ints 17 1073741823))
(⊦raises (exn hansei) (ocaml-random-int! used 0))
(⊦raises (exn hansei) (ocaml-random-int! used 1073741824))
(⊦raises (exn hansei) (rng-integer! src 0 1.0))
(⊦ equal?
'(230024973 197978033 495359752 1036278348 906260834)
(begin
(draws 100 (τ (ocaml-random-bits! used)))
(let1 (s (ocaml-random-init! used 17))
(draws 5 (τ (ocaml-random-bits! s))))))
(⊦ equal?
(ocaml-random-float! global-copy 1.0)
(rng-double! global-source 0.0 1.0))
(⊦ equal?
(ocaml-random-bits! global-copy)
(ocaml-random-bits! ocaml-random-global)))
((eta 0.0050000000000000001)
(memory #(12582912 2661168 1048576))
(stdout "")
(stderr ""))
test/replay/sampling: pass
The sampling records of samplingM.ml (SM:16-269, 398-401),
replayed bit for bit with probcc-random-selector: rejection and
importance sampling of shared and independent flips, the burglar alarm, the geometric
distribution, the observation of tiblt2 and of the Music TR, and
the drunk coin, whose dcoin_and_true 10 rejection sampling
accepts no world in 10000. The tables are in the order of probcc-explore (deviation 5), and the models draw in OCaml's order (D19).
(define (test/replay/sampling _)
(⊦ equal?
'(((V #f) 0.52000000000000002) ((V #t) 0.47999999999999998))
(probcc-sample-rejection/a (probcc-random-selector 1) 100 tflip2-shared))
(⊦ equal?
'(((V #f) 0.76000000000000001) ((V #t) 0.23999999999999999))
(probcc-sample-rejection/a
(probcc-random-selector 1)
100
tflip2-noshared))
(⊦ equal?
'(((V #f) 0.85999999999999999) ((V #t) 0.14000000000000001))
(probcc-sample-rejection/a (probcc-random-selector 1) 100 talarm))
(⊦ equal?
'(((V #f) 0.89128000000000007) ((V #t) 0.10872000000000001))
(probcc-sample-importance/a (probcc-random-selector 1) 100 talarm))
(⊦ equal?
'(((V 0) 0.79000000000000004)
((V 1) 0.17000000000000001)
((V 2) 0.040000000000000001))
(probcc-sample-rejection
(probcc-random-selector 17)
100
(probcc-geometric/count 0.84999999999999998)))
(⊦ equal?
'(((V 0) 0.84999999999999998)
((V 1) 0.1275)
((V 2) 0.019574999999999981)
((V 3) 0.0024750000000000011)
((V 4) 0.0004500000000000001))
(probcc-sample-importance
(probcc-random-selector 17)
100
(probcc-geometric/count 0.84999999999999998)))
(⊦ equal?
'(((V 0) 0.84999999999999998)
((V 1) 0.1275)
((V 2) 0.019080000000000406)
((V 3) 0.0028349999999999968)
((V 4) 0.00047250000000000032)
((V 5) 9.0000000000000019e-05)
((V 6) 2.2500000000000005e-05))
(probcc-sample-importance
(probcc-random-selector 17)
1000
(probcc-geometric/count 0.84999999999999998)))
(⊦ equal?
'(((V #f) 0.02))
(probcc-sample-rejection/a (probcc-random-selector 1) 100 tiblt2))
(⊦ equal?
'(((V #t) 0.14000000000000001))
(probcc-sample-importance/a (probcc-random-selector 1) 100 tmusictr41))
(⊦ equal?
'(((V #f) 0.050000000000000003))
(probcc-sample-rejection (probcc-random-selector 17) 100 (dcoin-and 10)))
(⊦ equal?
'(((V #f) 0.051400000000000001))
(probcc-sample-rejection
(probcc-random-selector 17)
10000
(dcoin-and 10)))
(⊦ equal?
'(((V #f) 0.052301965599999993))
(probcc-sample-importance (probcc-random-selector 17) 100 (dcoin-and 10)))
(⊦ equal?
'(((V #f) 0.052635478842930135) ((V #t) 6.9999999999999841e-14))
(probcc-sample-importance (probcc-random-selector 1) 5000 (dcoin-and 10)))
(⊦ equal?
'(((V #f) 0.052641559436170146) ((V #t) 5.9999999999999871e-14))
(probcc-sample-importance
(probcc-random-selector 17)
5000
(dcoin-and 10)))
(⊦ equal?
'()
(probcc-sample-rejection
(probcc-random-selector 1)
10000
(or (dcoin-and 10) (probcc-impossible)))))
((eta 0.17000000000000001)
(memory #(12582912 2662720 1048576))
(stdout "")
(stderr ""))
test/selector: pass
probcc-random-selector (cf. P:221-233) chooses an entry with the
probability of its mass and answers it with the total mass (deviation 18). Its seed
re-initialises OCaml's global Random state, as upstream's Random.init, so selectors made with the same seed share one
stream; an rng-source of (aux hakaru10 random) (OCaml's, the StdGen of Haskell, CHICKEN's) may replace the seed.
(define (test/selector _)
(define (draws n f)
(if (zero? n) '() (let1 (x (f)) (cons x (draws (sub1 n) f)))))
(define choices (probcc-map (λ (i) (list i 1/10)) (iota 10)))
(define (picks selector) (draws 5 (τ (car (selector choices)))))
(define s1 (probcc-random-selector 1))
(define first-picks (picks s1))
(define s2 (probcc-random-selector 1))
(⊦ equal? first-picks (picks s2))
(⊦ equal? (begin (probcc-random-selector 1) first-picks) (picks s1))
(define src (make-ocaml-random-source 1))
(define src-copy (rng-source-copy src))
(define global (ocaml-random-state-copy ocaml-random-global))
(⊦ equal? first-picks (picks (probcc-random-selector src)))
(⊦ equal?
(ocaml-random-bits! global)
(ocaml-random-bits! ocaml-random-global))
(⊦ equal? first-picks (picks (probcc-random-selector src-copy)))
(⊦ equal?
(picks (probcc-random-selector (make-stdgen-source 3)))
(picks (probcc-random-selector (make-stdgen-source 3))))
(⊦ equal?
(picks (probcc-random-selector
(make-chicken-random-source (string->utf8 "seed"))))
(picks (probcc-random-selector
(make-chicken-random-source (string->utf8 "seed")))))
(⊦ equal? '(b 1.0) ((probcc-random-selector 1) '((a 0.0) (b 1.0))))
(⊦ equal?
'((b 1.0) (b 1.0) (b 1.0))
(let1 (s (probcc-random-selector (make-stdgen-source 7)))
(draws 3 (τ (s '((a 0.0) (b 1.0) (c 0.0)))))))
(⊦ equal?
'(3 7/2)
(let1 (s (probcc-random-selector 17)) (s '((x 1/2) (3 3))))))
((eta 0.001) (memory #(12582912 2666728 1048576)) (stdout "") (stderr ""))
test/selector/raises: pass
A selector given no entries raises (exn hansei) of kind selector, upstream's Choice selection: can't
happen (P:224).
(define (test/selector/raises _)
(⊦raises (exn hansei) ((probcc-random-selector 1) '()))
(⊦ equal?
'selector
(handle-exceptions
e
((condition-property-accessor 'hansei 'kind) e)
((probcc-random-selector 1) '()))))
((eta 0.0) (memory #(12582912 2666760 1048576)) (stdout "") (stderr ""))
test/replay/nested: pass
The records of nested.ml (N:20-287, 346-441), replayed bit for bit:
a sampler inside a model that makes its own random_selector 17 re-seeds the global state, which the outer sampler then draws from (plan D6);
nested rejection and importance sampling with probcc-dist-selector give the distribution of the inner estimates, with the inner failures as empty
tables. Exact masses of 1 are the exact 1 where upstream has 1. (deviation 26), and
the tables are in the order of probcc-explore (deviation 5). The
red-light game (N:346-441) plans with risk: a cheat on a red light fails the world
when a nested 10-sample estimate of the chance of being caught reaches maxrisk.
(define (test/replay/nested _)
(define (testn12)
(let1 (coin (probcc-distribution
`((,(τ (probcc-coin 0.5)) 0.5) (,(τ #t) 0.5))))
(probcc-sample-rejection/a (probcc-random-selector 17) 10 coin)))
(define (testn22)
(let1 (coin-choice (probcc-coin 0.5))
(probcc-sample-rejection
(probcc-random-selector 17)
10
(or (probcc-coin 0.5) coin-choice))))
(define (testnl13/14 inferencer)
(τ (let* ((coin-choice (probcc-letlazy/nesting (τ (probcc-coin 0.5))))
(v1 (inferencer
(τ (let* ((v (probcc-coin 0.5))
(c (coin-choice))
(b (coin-choice)))
(list v b c)))))
(b (coin-choice))
(a (coin-choice)))
(list a b v1))))
(define (tf2-nested n selector)
(probcc-sample-rejection/a selector n tflip2-noshared))
(define (dcoin-and-true n) (or (dcoin-and n) (probcc-impossible)))
(define (nested-rejection n model)
(probcc-sample-rejection/a probcc-dist-selector n model))
(define (nested-importance n model)
(probcc-sample-importance/a probcc-dist-selector n model))
(define police-monitoring
'#(0.10000000000000001
0.20000000000000001
0.10000000000000001
0.20000000000000001
0.20000000000000001))
(define bad-luck 0.10000000000000001)
(define maxrisk 0.20000000000000001)
(define (noisy-or a astrength b bstrength baserate)
(or (and (probcc-coin astrength) (a))
(and (probcc-coin bstrength) (b))
(probcc-coin baserate)))
(define ((being-caught i t))
(or (and (probcc-coin 0.5) (probcc-coin (vector-ref police-monitoring i)))
(and (probcc-coin 0.29999999999999999)
(if (< t 3)
(probcc-coin 0.20000000000000001)
(probcc-coin 0.40000000000000002)))
(probcc-coin bad-luck)))
(define ((game cheatingp))
(let loop ((t 0) (pos 0))
(cond ((>= pos 5) t)
((probcc-coin 0.5)
(cond ((not (probcc-coin cheatingp)) (loop (add1 t) pos))
((probcc-at-least
maxrisk
#t
(nested-importance 10 (being-caught pos t)))
(probcc-impossible))
(else (loop (add1 t) (add1 pos)))))
(else (loop (add1 t) (add1 pos))))))
(define coin-tables
'(((V (((V #t) 1))) 0.5)
((V (((V #t) 0.59999999999999998) ((V #f) 0.40000000000000002))) 0.5)))
(⊦ equal? coin-tables (probcc-reify/exact/a testn12))
(⊦ equal?
coin-tables
(probcc-sample-rejection/a (probcc-random-selector 17) 10 testn12))
(⊦ equal?
'(((V (((V #t) 1.0))) 0.5)
((V (((V #t) 0.59999999999999998) ((V #f) 0.40000000000000002))) 0.5))
(probcc-reify/exact/a testn22))
(⊦ equal?
'(((V (#t #t (((V (#t #t #t)) 0.5) ((V (#f #t #t)) 0.5))))
0.59999999999999998)
((V (#f #f (((V (#t #f #f)) 0.5) ((V (#f #f #f)) 0.5))))
0.40000000000000002))
(probcc-sample-rejection/a
(probcc-random-selector 17)
10
(testnl13/14 probcc-reify/exact/a)))
(⊦ equal?
'(((V #f) 0.76000000000000001) ((V #t) 0.23999999999999999))
(tf2-nested 100 (probcc-random-selector 1)))
(⊦ equal?
'(((V (((V #f) 1.0))) 0.76000000000000001)
((V (((V #t) 1.0))) 0.23999999999999999))
(probcc-sample-rejection
(probcc-random-selector 1)
100
(tf2-nested 1 probcc-dist-selector)))
(⊦ equal?
'(((V (((V #f) 1.0))) 0.41999999999999998)
((V (((V #f) 0.66666666666666663) ((V #t) 0.33333333333333331)))
0.40000000000000002)
((V (((V #t) 0.66666666666666663) ((V #f) 0.33333333333333331)))
0.17999999999999999))
(probcc-sample-rejection
(probcc-random-selector 1)
100
(tf2-nested 3 probcc-dist-selector)))
(⊦ equal?
'(((V (#f
#f
(((V (#f #f #f)) 0.69999999999999996)
((V (#t #f #f)) 0.29999999999999999))))
0.29999999999999999)
((V (#t #t (((V (#t #t #t)) 0.5) ((V (#f #t #t)) 0.5))))
0.20000000000000001)
((V (#t
#t
(((V (#t #t #t)) 0.69999999999999996)
((V (#f #t #t)) 0.29999999999999999))))
0.10000000000000001)
((V (#f
#f
(((V (#t #f #f)) 0.69999999999999996)
((V (#f #f #f)) 0.29999999999999999))))
0.10000000000000001)
((V (#f #f (((V (#t #f #f)) 0.5) ((V (#f #f #f)) 0.5))))
0.10000000000000001)
((V (#f
#f
(((V (#f #f #f)) 0.80000000000000004)
((V (#t #f #f)) 0.20000000000000001))))
0.10000000000000001)
((V (#f
#f
(((V (#f #f #f)) 0.59999999999999998)
((V (#t #f #f)) 0.40000000000000002))))
0.10000000000000001))
(probcc-sample-rejection/a
(probcc-random-selector 17)
10
(testnl13/14 (λ (m) (nested-rejection 10 m)))))
(⊦ equal?
'(((V ()) 0.9476980344) ((V (((V #f) 1.0))) 0.052301965599999993))
(probcc-sample-importance
(probcc-random-selector 17)
100
(nested-rejection 1 (τ (dcoin-and 10)))))
(⊦ equal?
'(((V (((V #f) 0.054999999999999986))) 0.5)
((V (((V #f) 0.050499999999999989))) 0.25)
((V (((V #f) 0.05004999999999999))) 0.125)
((V (((V #f) 0.050004999999999987))) 0.0625)
((V (((V #f) 0.050000499999999989))) 0.03125)
((V (((V #f) 0.05000004999999999))) 0.015625)
((V (((V #f) 0.050000004999999986))) 0.0078125)
((V (((V #f) 0.050000000499999989))) 0.00390625)
((V (((V #f) 0.050000000049999986))) 0.001953125)
((V (((V #f) 0.049999999999999989) ((V #t) 4.9999999999999885e-11)))
0.001953125))
(probcc-sample-importance
(probcc-random-selector 17)
100
(nested-importance 1 (τ (dcoin-and 10)))))
(⊦ equal?
'(((V ()) 0.99999999999999989))
(probcc-sample-importance
(probcc-random-selector 1)
100
(nested-rejection 1 (τ (dcoin-and-true 10)))))
(⊦ equal?
'(((V ()) 0.998046875)
((V (((V #t) 4.9999999999999885e-11))) 0.001953125))
(probcc-sample-importance
(probcc-random-selector 1)
100
(nested-importance 1 (τ (dcoin-and-true 10)))))
(⊦ equal?
'(((V ()) 0.88826171875000004)
((V (((V #t) 9.9999999999999776e-13))) 0.10154296875)
((V (((V #t) 1.9999999999999955e-12))) 0.01017578125)
((V (((V #t) 2.9999999999999933e-12))) 1.953125e-05))
(probcc-sample-importance
(probcc-random-selector 1)
100
(nested-importance 50 (τ (dcoin-and-true 10)))))
(⊦ equal?
'(((V #t) 11/20) ((V #f) 9/20))
(probcc-reify/exact (noisy-or (τ #t) 1/2 (τ #f) 1/2 1/10)))
(⊦ equal?
'(((V 9) 0.10625)
((V 6) 0.09999999999999995)
((V 13) 0.053515624999999997)
((V 10) 0.053124999999999999)
((V 7) 0.049999999999999975)
((V 14) 0.026757812499999999)
((V 11) 0.026562499999999999)
((V 8) 0.024999999999999988)
((V 12) 0.01328125)
((V 15) 0.01259765625)
((V 16) 4.8828124999999949e-05))
(probcc-sample-importance/a
(probcc-random-selector 17)
10
(game 0.20000000000000001)))
(⊦ equal?
'(((V 7) 0.052499999999999998)
((V 5) 0.040000000000000001)
((V 6) 0.025000000000000001)
((V 8) 0.016250000000000001)
((V 9) 0.010625000000000001))
(probcc-sample-importance/a (probcc-random-selector 17) 50 (game 0.5))))
((eta 0.13600000000000001)
(memory #(12582912 2688376 1048576))
(stdout "")
(stderr ""))
(aux hansei lazy): lazy lists, permutations and sorted lists
Correlated memoization and lazy permutations
Lazy lists, shared lazily-evaluated variables and memoized stochastic functions, after Correlated memoization and lazy permutations, slazy.ml and Realistic Example: sorted list:
❝
Sampling of permutations, stochastic bijective functions, sequences of correlated random variables, the all-different constraint -- are related to the same problem. A problem with so many names must be important, and so must be its efficient solution.
The goal is a sequence of N all-different random variables over 0 .. N-1, a stochastic
bijection: the naive f1 of LP:36-41, which draws uniformly_list [0;1] at every application, is not even a function
(test test/correlated); probcc-memo makes it one, and probcc-lazy-permutation (LP:127-229) and probcc-lazy-permutation/ins (LP:266-330) correlate the memoized
values, so that each one is drawn only when asked for, among the values still free.
Laziness and non-determinism
Laziness and non-determinism make a tricky combination (SL:10-93): a variable bound with let is computed in every world, also where it is not needed; a thunk
is not shared; and the built-in delay is shared across worlds, which
gives a wrong distribution (DSL p.21 l.912-922). probcc-letlazy evaluates its body at most once per world, when it is first needed (P:270-330).
Realistic example: sorted list
Harik's problem (H #sorted, SO:1-11) draws N numbers uniformly
from 0 .. M-1 and asks for the probability that they are sorted; with numbers as lazy
lists of bits, the comparison stops at the first differing bit, which is crucial:
❝The exact probability for four numbers of five bits, over 2^20 worlds, and upstream's importance-sampling runs of it are checked in
Harik's problem is ideal for a benchmark: it is easy to explain and to model. The expected result is also known. Getting any useful estimate of the result is rather difficult, especially without extensive fine-tuning of the inference algorithm.
test/sorted.test/llist/kit: pass
The lazy-list kit of SL:325-355, M1:51-101 and M2:113-119, promoted to a library
(deviation 24): '() and pairs (head-thunk . tail) stand for LNil and LCons; probcc-lappend and probcc-lmap are built with probcc-letlazy, so their cells are shared in a
world; probcc-lobserve fails as soon as a head or the length
differs, forcing no more than that.
(define (test/llist/kit _)
(define (to-list l)
(match/first
(l)
(() '())
((,h unquote t) (let1 (x (h)) (cons x (to-list t))))))
(define abc (probcc-list->llist '(a b c)))
(⊦ equal? '() (to-list probcc-lnil))
(⊦ equal? '(a b c) (to-list abc))
(⊦ equal? '(x) (to-list (probcc-lcons (τ 'x) probcc-lnil)))
(⊦ equal? 3 (probcc-llength abc))
(⊦ equal? 0 (probcc-llength probcc-lnil))
(⊦ equal? 2 (probcc-llength/at-most 2 abc))
(⊦ equal? 3 (probcc-llength/at-most 7 abc))
(⊦ equal? 0 (probcc-llength/at-most 0 abc))
(⊦ equal? 'a (probcc-lreflect (probcc-lhead abc)))
(⊦ equal? '(b c) (to-list (probcc-ltail abc)))
(⊦ equal? 'z (probcc-lreflect (τ 'z)))
(⊦ equal?
'(((V (a b c d e)) 1))
(probcc-reify/exact
(to-list (probcc-lappend abc (probcc-list->llist '(d e))))))
(⊦ equal?
'(((V (a b c)) 1))
(probcc-reify/exact (to-list (probcc-lappend probcc-lnil abc))))
(⊦ equal?
'(((V (10 20 30)) 1))
(probcc-reify/exact
(to-list (probcc-lmap (λ (x) (* 10 x)) (probcc-list->llist '(1 2 3))))))
(⊦ equal? '() (probcc-reify/exact (probcc-lhead probcc-lnil)))
(⊦ equal? '() (probcc-reify/exact (probcc-ltail probcc-lnil)))
(⊦ equal?
'(((V ((#t #t) (#t #t))) 1/4)
((V ((#t #f) (#t #f))) 1/4)
((V ((#f #t) (#f #t))) 1/4)
((V ((#f #f) (#f #f))) 1/4))
(probcc-reify/exact
(let* ((l (probcc-lmap
(λ (_) (probcc-coin 1/2))
(probcc-list->llist '(1 2))))
(b (to-list l))
(a (to-list l)))
(list a b))))
(⊦ equal?
'(((V (a b)) 1/2) ((V (a)) 1/2))
(probcc-reify/exact
(let* ((l (probcc-lappend
(probcc-list->llist '(a))
(τ (if (probcc-coin 1/2) '() (cons (τ 'b) probcc-lnil)))))
(b (to-list l))
(a (to-list l)))
(if (equal? a b) a 'unshared))))
(⊦ equal?
`(((V ,(void)) 1))
(probcc-reify/exact (probcc-lobserve abc '(a b c))))
(⊦ equal? '() (probcc-reify/exact (probcc-lobserve abc '(a x c))))
(⊦ equal? '() (probcc-reify/exact (probcc-lobserve abc '(a b))))
(⊦ equal? '() (probcc-reify/exact (probcc-lobserve abc '(a b c d))))
(⊦ equal?
`(((V ,(void)) 1/4))
(probcc-reify/exact
(probcc-lobserve
(probcc-lmap (λ (_) (probcc-coin 1/2)) (probcc-list->llist '(1 2)))
'(#t #f)))))
((eta 0.0) (memory #(8396592 3282944 1048576)) (stdout "") (stderr ""))
test/lazy/sharing: pass
Call-by-need and non-determinism (SL:10-93): a strict x is computed 4 times, also in the world u = 0 where it is not
needed; a thunk is not shared, so the two components may differ (12 calls); a
promise of delay is shared too widely, by the worlds that split
inside its first force, so the first world to finish fixes the
value for the others and the distribution is wrong (tests3, SL:57-70, where OCaml's Lazy is wrong too, differently); probcc-letlazy gives the strict distribution and computes x only where it is needed, 4 times: the world is split by the coin
before x is forced (SL:85-93).
(define (test/lazy/sharing _)
(⊦ equal?
'((((V #f) 1/2)
((V (21 21 #t)) 1/8)
((V (21 21 #f)) 1/8)
((V (11 11 #t)) 1/8)
((V (11 11 #f)) 1/8))
4)
(with-calls (τ (probcc-reify/exact/a (tests1 1/2)))))
(⊦ equal?
'((((V #f) 1/2)
((V (21 21 #t)) 1/16)
((V (21 21 #f)) 1/16)
((V (21 11 #t)) 1/16)
((V (21 11 #f)) 1/16)
((V (11 21 #t)) 1/16)
((V (11 21 #f)) 1/16)
((V (11 11 #t)) 1/16)
((V (11 11 #f)) 1/16))
12)
(with-calls (τ (probcc-reify/exact/a (tests2 1/2)))))
(⊦ equal?
'((((V #f) 1/2) ((V (11 11 #t)) 1/4) ((V (11 11 #f)) 1/4)) 2)
(with-calls (τ (probcc-reify/exact/a (tests3 1/2)))))
(⊦ equal?
'((((V #f) 1/2)
((V (21 21 #t)) 1/8)
((V (21 21 #f)) 1/8)
((V (11 11 #t)) 1/8)
((V (11 11 #f)) 1/8))
4)
(with-calls (τ (probcc-reify/exact/a (tests4 1/2))))))
((eta 0.0) (memory #(8396592 3286960 1048576)) (stdout "") (stderr ""))
test/lazy/variables: pass
Lazy variables in testl1 .. testl54 (SL:99-284): the lazy and the strict models have the same distribution, and the
expensive function is called as often as upstream's prints its line; a variable that
escapes the scope of its letlazy, forced after the world splits,
is computed in each world (8 calls in testl4).
(define (test/lazy/variables _)
(define table
'(((V (2 2)) 3/8)
((V (1 1)) 1/4)
((V (6 6)) 1/8)
((V (4 4)) 1/8)
((V (3 3)) 1/8)))
(define table3
'(((V (2 #t 200 200)) 1/8)
((V (2 #t 20 20)) 1/8)
((V (2 #f 200 200)) 1/8)
((V (2 #f 20 20)) 1/8)
((V (1 #t 100 100)) 1/8)
((V (1 #t 10 10)) 1/8)
((V (1 #f 100 100)) 1/8)
((V (1 #f 10 10)) 1/8)))
(⊦ equal?
(list table 4)
(with-calls (τ (probcc-reify/exact/a (testl1 1/2)))))
(⊦ equal?
(list table 4)
(with-calls (τ (probcc-reify/exact/a (testl2 1/2)))))
(⊦ equal?
(list table3 4)
(with-calls (τ (probcc-reify/exact/a (testl3 1/2)))))
(⊦ equal?
(list table3 8)
(with-calls (τ (probcc-reify/exact/a (testl4 1/2)))))
(⊦ equal?
(list table 4)
(with-calls (τ (probcc-reify/exact/a (testl51 1/2)))))
(⊦ equal?
(list table 4)
(with-calls (τ (probcc-reify/exact/a (testl52 1/2)))))
(⊦ equal?
(list table 4)
(with-calls (τ (probcc-reify/exact/a (testl53 1/2)))))
(⊦ equal?
'((((V (#t 2 2 2 2)) 1/4)
((V (#t 1 1 1 1)) 1/4)
((V (#f 2 2 2 2)) 1/4)
((V (#f 1 1 1 1)) 1/4))
8)
(with-calls (τ (probcc-reify/exact/a (testl54 1/2))))))
((eta 0.001) (memory #(8396592 3289088 1048576)) (stdout "") (stderr ""))
test/llist/at-most: pass
Ken Shan's test (SL:349-367): the length of a lazy list of unbounded length, counting at most 5 cells, explores finitely many worlds.
(define (test/llist/at-most _)
(⊦ equal?
'(((V 0) 1/2)
((V 1) 1/4)
((V 2) 1/8)
((V 3) 1/16)
((V 5) 1/32)
((V 4) 1/32))
(probcc-reify/exact
(probcc-llength/at-most 5 (probcc-letlazy (ll-count 0))))))
((eta 0.0) (memory #(8396592 3289408 1048576)) (stdout "") (stderr ""))
test/llist/pcfg: pass
The probabilistic context-free grammar of IBAL's Music TR (SL:288-478): noun phrases of unbounded length, as lazy lists that are forced only as far as the question needs; the tables are upstream's bit for bit (SL:387, 401, 451-454, 468).
(define (test/llist/pcfg _)
(⊦ equal?
'(((V ("ants" "ants")) 0.59999999999999998)
((V ("flies" "flies")) 0.39999999999999997))
(probcc-reify/exact/a pcfg11))
(⊦ equal?
'(((V ("ants" "ants")) 0.17999999999999999) ((V ("flies" "flies")) 0.12))
(probcc-reify/exact/a pcfg12))
(⊦ equal?
'(((V ("ants" "ants" "ants" "ants")) 0.108)
((V ("flies" "flies" "ants" "ants")) 0.072000000000000008)
((V ("ants" "ants" "flies" "flies")) 0.071999999999999995)
((V ("flies" "flies" "flies" "flies")) 0.048000000000000001))
(probcc-reify/exact/a pcfg13))
(⊦ equal?
'(((V ("flies" "ants")) 0.071999999999999995)
((V ("flies" "flies")) 0.048000000000000001))
(probcc-reify/exact/a pcfg)))
((eta 0.0) (memory #(8396592 3291856 1048576)) (stdout "") (stderr ""))
test/memo/functions: pass
Without memoization, two applications of a stochastic function to the same
argument may give different results in the same world (test1m, SL:661-692); test2m and test3m of hansei.scm show probcc-memo. The sprinkler of SL:731-778 memoizes rain, sprinkler and wet grass as functions
of the day: the grass of a day is wet or not once in a world, and given that it is
wet on day 2, it rained with the weights 0.2838 and 0.322, those of the lawn in the suite hansei.
(define (test/memo/functions _)
(⊦ equal?
'(((V (#t 2 2 3 3)) 1/32)
((V (#t 2 2 3 2)) 1/32)
((V (#t 2 2 2 3)) 1/32)
((V (#t 2 2 2 2)) 1/32)
((V (#t 2 1 3 3)) 1/32)
((V (#t 2 1 3 2)) 1/32)
((V (#t 2 1 2 3)) 1/32)
((V (#t 2 1 2 2)) 1/32)
((V (#t 1 2 3 3)) 1/32)
((V (#t 1 2 3 2)) 1/32)
((V (#t 1 2 2 3)) 1/32)
((V (#t 1 2 2 2)) 1/32)
((V (#t 1 1 3 3)) 1/32)
((V (#t 1 1 3 2)) 1/32)
((V (#t 1 1 2 3)) 1/32)
((V (#t 1 1 2 2)) 1/32)
((V (#f 3 3 2 2)) 1/32)
((V (#f 3 3 2 1)) 1/32)
((V (#f 3 3 1 2)) 1/32)
((V (#f 3 3 1 1)) 1/32)
((V (#f 3 2 2 2)) 1/32)
((V (#f 3 2 2 1)) 1/32)
((V (#f 3 2 1 2)) 1/32)
((V (#f 3 2 1 1)) 1/32)
((V (#f 2 3 2 2)) 1/32)
((V (#f 2 3 2 1)) 1/32)
((V (#f 2 3 1 2)) 1/32)
((V (#f 2 3 1 1)) 1/32)
((V (#f 2 2 2 2)) 1/32)
((V (#f 2 2 2 1)) 1/32)
((V (#f 2 2 1 2)) 1/32)
((V (#f 2 2 1 1)) 1/32))
(probcc-reify/exact (test1m)))
(⊦ equal?
'(((V #t) 1.0))
(probcc-reify/exact/a
(grass-model
(λ (grass-is-wet rain)
(let* ((b (grass-is-wet 2)) (a (grass-is-wet 2))) (equal? a b))))))
(⊦ equal?
'(((V #f) 0.32200000000000001) ((V #t) 0.2838))
(probcc-reify/exact/a
(grass-model
(λ (grass-is-wet rain)
(probcc-observe (λ (_) (grass-is-wet 2)) (τ (rain 2))))))))
((eta 0.0) (memory #(8396592 3293288 1048576)) (stdout "") (stderr ""))
test/lazy/delayed: pass
Lazy versus delayed evaluation (SL:782-860): an expression whose value is not
used is not evaluated, so its failure is not seen; lazy evaluation is sound only for
expressions that cannot fail, or when only the ratios of the final probabilities
matter. In testd1 the lazy x fails in the worlds where it is forced, in testd2 it is never forced.
(define (test/lazy/delayed _)
(⊦ equal?
'((((V (1 1)) 1/2) ((V (2 2)) 1/4)) 4)
(with-calls (τ (probcc-reify/exact/a testd1))))
(⊦ equal? '(((V 2) 1/2) ((V 1) 1/2)) (probcc-reify/exact/a testd2)))
((eta 0.0) (memory #(8396592 3294840 1048576)) (stdout "") (stderr ""))
test/split-at: pass
probcc-split-at (cf. LP:139-153, LPI:32-37) answers the zipper of
a list at an element, the elements before it reversed; a list that is too short gives
an empty tail.
(define (test/split-at _)
(⊦ equal? '(() (1 2 3)) (probcc-split-at 0 '(1 2 3)))
(⊦ equal? '((1) (2 3)) (probcc-split-at 1 '(1 2 3)))
(⊦ equal? '((2 1) (3)) (probcc-split-at 2 '(1 2 3)))
(⊦ equal? '((3 2 1) ()) (probcc-split-at 3 '(1 2 3)))
(⊦ equal? '((3 2 1) ()) (probcc-split-at 5 '(1 2 3)))
(⊦ equal? '(() ()) (probcc-split-at 2 '())))
((eta 0.0) (memory #(8396592 3295688 1048576)) (stdout "") (stderr ""))
test/uncons: pass
probcc-uniform/either/uncons (cf. LP:155-172) chooses an element
uniformly and answers it with the others; a one-element list makes no choice, and the
empty list fails.
(define (test/uncons _)
(⊦ equal?
'(((V (3 (1 2))) 1/3) ((V (2 (1 3))) 1/3) ((V (1 (2 3))) 1/3))
(probcc-reify/exact (probcc-uniform/either/uncons '(1 2 3))))
(⊦ equal?
'(((V (a ())) 1))
(probcc-reify/exact (probcc-uniform/either/uncons '(a))))
(⊦ equal?
#t
(match/first
(probcc-reify/0 (τ (probcc-uniform/either/uncons '(a))))
((((V (a ())) 1)) #t)
(else #f)))
(⊦ equal? '() (probcc-reify/exact (probcc-uniform/either/uncons '()))))
((eta 0.0) (memory #(8396592 3297800 1048576)) (stdout "") (stderr ""))
test/insert-shift: pass
probcc-insert-shift (cf. LP:233-264) inserts into a sorted list
the integer j given as j minus the number of smaller elements.
(define (test/insert-shift _)
(⊦ equal? '(0 (0)) (probcc-insert-shift 0 '()))
(⊦ equal? '(1 (0 1)) (probcc-insert-shift 0 '(0)))
(⊦ equal? '(0 (0 1)) (probcc-insert-shift 0 '(1)))
(⊦ equal? '(2 (1 2)) (probcc-insert-shift 1 '(1)))
(⊦ equal? '(3 (1 3)) (probcc-insert-shift 2 '(1)))
(⊦ equal? '(2 (0 1 2)) (probcc-insert-shift 0 '(0 1)))
(⊦ equal? '(0 (0 1 2)) (probcc-insert-shift 0 '(1 2)))
(⊦ equal? '(3 (1 2 3)) (probcc-insert-shift 1 '(1 2))))
((eta 0.0) (memory #(8396592 3298736 1048576)) (stdout "") (stderr ""))
test/correlated: pass
Stochastic functions (LP:17-107, H #all-different): f1 is not a function, the memoized f2 is; f3 rejects the worlds where f2 is not injective, wasting half of the choice space; f4 draws a whole permutation lazily.
(define (test/correlated _)
(⊦ equal?
(map (λ (l) (list (list 'V l) 1/16))
'((1 1 1 1)
(1 1 1 0)
(1 1 0 1)
(1 1 0 0)
(1 0 1 1)
(1 0 1 0)
(1 0 0 1)
(1 0 0 0)
(0 1 1 1)
(0 1 1 0)
(0 1 0 1)
(0 1 0 0)
(0 0 1 1)
(0 0 1 0)
(0 0 0 1)
(0 0 0 0)))
(probcc-reify/exact (tuple-0101 f01)))
(⊦ equal?
'(((V (1 1 1 1)) 1/4)
((V (1 0 1 0)) 1/4)
((V (0 1 0 1)) 1/4)
((V (0 0 0 0)) 1/4))
(probcc-reify/exact (tuple-0101 (probcc-memo f01))))
(⊦ equal?
'(((V (1 0 1 0)) 1/4) ((V (0 1 0 1)) 1/4))
(probcc-reify/exact
(let* ((f2 (probcc-memo f01))
(f3 (λ (x)
(let* ((b (f2 1)) (a (f2 0)))
(when (eqv? a b) (probcc-impossible))
(f2 x)))))
(tuple-0101 f3))))
(⊦ equal?
'(((V (1 0 1 0)) 1/2) ((V (0 1 0 1)) 1/2))
(probcc-reify/exact
(let* ((perms (probcc-letlazy
(τ (probcc-uniform/either '((0 1) (1 0))))))
(p (perms)))
(tuple-0101 (λ (i) (list-ref p i)))))))
((eta 0.0) (memory #(8396592 3302504 1048576)) (stdout "") (stderr ""))
test/lazy-permutation: pass
probcc-lazy-permutation (cf. LP:127-229) is a list of correlated
random variables, each choosing among the elements its predecessors left: the 24
permutations of four elements have the exact mass 1/24 (upstream's floats sum to
0.999999999999999667, LP:190), and the lazy f5 is an injective function that wastes no choice space.
(define (test/lazy-permutation _)
(define raw (map (λ (l) (list (list 'V l) 1/24)) permutations-of-4))
(define (f5-0101 perms) (tuple-0101 (λ (n) ((list-ref perms n)))))
(⊦ equal?
raw
(probcc-reify/exact
(probcc-map (λ (th) (th)) (probcc-lazy-permutation '(1 2 3 4)))))
(⊦ equal?
(list 1 (map (λ (e) (list (car e) (exact->inexact 1/24))) raw))
(probcc-normalize/total
(probcc-reify/exact
(probcc-map (λ (th) (th)) (probcc-lazy-permutation '(1 2 3 4))))))
(⊦ equal?
'(((V (3 2 3 2)) 1/6)
((V (3 1 3 1)) 1/6)
((V (2 3 2 3)) 1/6)
((V (2 1 2 1)) 1/6)
((V (1 3 1 3)) 1/6)
((V (1 2 1 2)) 1/6))
(probcc-reify/exact (f5-0101 (probcc-lazy-permutation '(1 2 3)))))
(⊦ equal?
3
(probcc-leaves
(probcc-reify/0 (τ ((car (probcc-lazy-permutation '(1 2 3))))))))
(⊦ equal? '() (probcc-lazy-permutation '()))
(⊦ equal?
'(((V (a)) 1))
(probcc-reify/exact
(probcc-map (λ (th) (th)) (probcc-lazy-permutation '(a))))))
((eta 0.0) (memory #(8396592 3308880 1048576)) (stdout "") (stderr ""))
test/lazy-permutation/ins: pass
probcc-lazy-permutation/ins (cf. LP:266-330) builds the same
sequences of all-different variables over 0 .. n-1 by insertion; a count that is
not a non-negative exact integer raises (deviation 33).
(define (test/lazy-permutation/ins _)
(define (f5-0101 perms) (tuple-0101 (λ (n) ((list-ref perms n)))))
(⊦ equal?
(map (λ (l) (list (list 'V (map sub1 l)) 1/24)) permutations-of-4)
(probcc-reify/exact
(probcc-map (λ (th) (th)) (probcc-lazy-permutation/ins 4))))
(⊦ equal?
'(((V (2 1 2 1)) 1/6)
((V (2 0 2 0)) 1/6)
((V (1 2 1 2)) 1/6)
((V (1 0 1 0)) 1/6)
((V (0 2 0 2)) 1/6)
((V (0 1 0 1)) 1/6))
(probcc-reify/exact (f5-0101 (probcc-lazy-permutation/ins 3))))
(⊦ equal? '() (probcc-lazy-permutation/ins 0))
(⊦ equal?
'(((V (0)) 1))
(probcc-reify/exact
(probcc-map (λ (th) (th)) (probcc-lazy-permutation/ins 1))))
(⊦raises (exn hansei) (probcc-lazy-permutation/ins -1))
(⊦raises (exn hansei) (probcc-lazy-permutation/ins 0.5))
(⊦raises (exn hansei) (probcc-lazy-permutation/ins 2.5)))
((eta 0.0) (memory #(8396592 3311816 1048576)) (stdout "") (stderr ""))
test/sorted: pass
Harik's sorted-list problem (SO:50-200): numbers are lazy lists of random bits,
so gte forces only the bits up to the first difference; the
exact probabilities agree with the theory, testp3_theory = 15/64, testp4_theory = 165/2048 and testp5_theory = 6545/131072 (SO:146-171), the last over 2^20
worlds of which only those that gte looks at are generated.
Importance sampling with (probcc-random-selector 17) gives what upstream gives under OCaml today, sample for sample, both for the
four numbers of five bits and for the full model of n numbers of 20 bits, whose
probability of being sorted is 1/n!.
(define (test/sorted _)
(⊦ equal?
(map (λ (l) (list (list 'V l) 1/16))
'((1 1 1 1)
(1 1 1 0)
(1 1 0 1)
(1 1 0 0)
(1 0 1 1)
(1 0 1 0)
(1 0 0 1)
(1 0 0 0)
(0 1 1 1)
(0 1 1 0)
(0 1 0 1)
(0 1 0 0)
(0 0 1 1)
(0 0 1 0)
(0 0 0 1)
(0 0 0 0)))
(probcc-reify/exact (so-to-list (so-bit-list 4))))
(⊦ equal?
(map (λ (l) (list (list 'V (list l l)) 1/8))
'((1 1 1) (1 1 0) (1 0 1) (1 0 0) (0 1 1) (0 1 0) (0 0 1) (0 0 0)))
(probcc-reify/exact
(let* ((l (so-bit-list 3)) (b (so-to-list l)) (a (so-to-list l)))
(list a b))))
(⊦ equal?
'(((V #t) 5/8) ((V #f) 3/8))
(probcc-reify/exact
(let* ((l1 (so-bit-list 2)) (l2 (so-bit-list 2))) (so-gte l1 l2))))
(⊦ equal?
'(((V #t) 9/16) ((V #f) 7/16))
(probcc-reify/exact
(let* ((l1 (so-bit-list 3)) (l2 (so-bit-list 3))) (so-gte l1 l2))))
(⊦ equal?
'(((V #t) 5/8) ((V #f) 3/8))
(probcc-reify/exact (so-is-sorted so-gte (so-num-list 2 2))))
(⊦ equal?
'(((V #t) 9/16) ((V #f) 7/16))
(probcc-reify/exact (so-is-sorted so-gte (so-num-list 2 3))))
(⊦ equal?
'(((V #f) 49/64) ((V #t) 15/64))
(probcc-reify/exact (so-is-sorted so-gte (so-num-list 3 3))))
(⊦ equal?
15/64
(/ (so-sum 0 7 (λ (i) (so-sum i 7 (λ (j) (so-sum j 7 (λ (k) 1))))))
(* 64 8)))
(⊦ equal?
'(((V #f) 1883/2048) ((V #t) 165/2048))
(probcc-reify/exact (so-is-sorted so-gte (so-num-list 4 3))))
(⊦ equal?
165/2048
(/ (so-sum
0
7
(λ (i0)
(so-sum
i0
7
(λ (i) (so-sum i 7 (λ (j) (so-sum j 7 (λ (k) 1))))))))
(* 64 64)))
(⊦ equal?
'(((V #f) 124527/131072) ((V #t) 6545/131072))
(probcc-reify/exact (so-is-sorted so-gte (so-num-list 4 5))))
(⊦ equal?
6545/131072
(/ (so-sum
0
31
(λ (i0)
(so-sum
i0
31
(λ (i) (so-sum i 31 (λ (j) (so-sum j 31 (λ (k) 1))))))))
(* 32 32 32 32)))
(⊦ equal?
'(((V #f) 77429/81920) ((V #t) 4491/81920))
(probcc-sample-importance/a
(probcc-random-selector 17)
100
(τ (so-is-sorted so-gte (so-num-list 4 5)))))
(⊦ equal?
'(((V #f) 155999/163840) ((V #t) 7841/163840))
(probcc-sample-importance/a
(probcc-random-selector 17)
200
(τ (so-is-sorted so-gte (so-num-list 4 5)))))
(⊦ equal?
`(((V ,(void)) 51716307/104857600))
(probcc-sample-importance/a (probcc-random-selector 17) 100 (so-model 2)))
(⊦ equal?
`(((V ,(void)) 55302254749/343597383680))
(probcc-sample-importance/a (probcc-random-selector 17) 100 (so-model 3)))
(⊦ equal?
`(((V ,(void)) 97798819259/2748779069440))
(probcc-sample-importance/a (probcc-random-selector 17) 200 (so-model 4)))
(⊦ equal?
'(0.49320513725280762 0.16095074460899922 0.035579003182283485)
(probcc-map
exact->inexact
(list 51716307/104857600
55302254749/343597383680
97798819259/2748779069440))))
((eta 0.45700000000000002)
(memory #(8396592 3316520 1048576))
(stdout "")
(stderr ""))
test/replay/slazy: pass
The sampling records of slazy.ml (SL:125-211, 431-485, 530,
617-655), replayed bit for bit with probcc-random-selector, with upstream's flonum masses (the digits through a transcription of upstream's uniform, P:156-165, since probcc-uniform's are exact, deviation 7), and the calls of the expensive function counted as
upstream prints them: lazy variables are computed only in the sampled worlds that
need them. With lazy lists, importance sampling finds the exact probability of
all-ones digits and of 20 heads in one sample, where 2000 samples of the strict
lists find nothing (SL:617-655).
(define (test/replay/slazy _)
(define (rejection seed n thunk)
(with-calls
(τ (probcc-sample-rejection/a (probcc-random-selector seed) n thunk))))
(define (importance seed n thunk)
(with-calls
(τ (probcc-sample-importance/a (probcc-random-selector seed) n thunk))))
(⊦ equal?
'((((V (2 2)) 0.40000000000000002)
((V (4 4)) 0.20000000000000001)
((V (3 3)) 0.20000000000000001)
((V (6 6)) 0.10000000000000001)
((V (1 1)) 0.10000000000000001))
10)
(rejection 1 10 (testl1 0.5)))
(⊦ equal?
'((((V (2 2)) 0.5)
((V (1 1)) 0.20000000000000001)
((V (6 6)) 0.10000000000000001)
((V (4 4)) 0.10000000000000001)
((V (3 3)) 0.10000000000000001))
4)
(rejection 1 10 (testl2 0.5)))
(⊦ equal?
'((((V (2 #f 20 20)) 0.20000000000000001)
((V (1 #f 100 100)) 0.20000000000000001)
((V (1 #f 10 10)) 0.20000000000000001)
((V (2 #t 200 200)) 0.10000000000000001)
((V (2 #t 20 20)) 0.10000000000000001)
((V (2 #f 200 200)) 0.10000000000000001)
((V (1 #t 100 100)) 0.10000000000000001))
10)
(rejection 1 10 (testl3 0.5)))
(⊦ equal?
'((((V (2 #t 200 200)) 0.125)
((V (2 #t 20 20)) 0.125)
((V (2 #f 200 200)) 0.125)
((V (2 #f 20 20)) 0.125)
((V (1 #t 100 100)) 0.125)
((V (1 #t 10 10)) 0.125)
((V (1 #f 100 100)) 0.125)
((V (1 #f 10 10)) 0.125))
4)
(importance 1 10 (testl3 0.5)))
(⊦ equal?
'((((V (2 #t 200 200)) 0.20000000000000001)
((V (1 #t 100 100)) 0.20000000000000001)
((V (1 #f 100 100)) 0.20000000000000001)
((V (2 #t 20 20)) 0.10000000000000001)
((V (2 #f 200 200)) 0.10000000000000001)
((V (2 #f 20 20)) 0.10000000000000001)
((V (1 #f 10 10)) 0.10000000000000001))
10)
(rejection 1 10 (testl4 0.5)))
(⊦ equal?
'((((V (2 #t 200 200)) 0.125)
((V (2 #t 20 20)) 0.125)
((V (2 #f 200 200)) 0.125)
((V (2 #f 20 20)) 0.125)
((V (1 #t 100 100)) 0.125)
((V (1 #t 10 10)) 0.125)
((V (1 #f 100 100)) 0.125)
((V (1 #f 10 10)) 0.125))
8)
(importance 1 10 (testl4 0.5)))
(⊦ equal?
'(((V ("ants" "ants")) 0.17999999999999999) ((V ("flies" "flies")) 0.12))
(probcc-sample-importance/a (probcc-random-selector 1) 100 pcfg12))
(⊦ equal?
'(((V ("flies" "ants")) 0.07200000000000005)
((V ("flies" "flies")) 0.048000000000000036))
(probcc-sample-importance/a (probcc-random-selector 1) 100 pcfg))
(⊦ equal?
`(((V ,(void)) 0.12))
(probcc-sample-rejection/a
(probcc-random-selector 17)
100
(full-10tree-obs ml-uniform 1)))
(⊦ equal?
'(((V 1) 0.00010000000000000003))
(probcc-sample-importance/a
(probcc-random-selector 17)
1
(full-10ltree-obs ml-uniform 4)))
(⊦ equal?
'(((V 1) 1/10000))
(probcc-sample-importance/a
(probcc-random-selector 17)
1
(full-10ltree-obs probcc-uniform 4)))
(⊦ equal?
'()
(probcc-sample-importance
(probcc-random-selector 17)
2000
(if (trues 20 (flips 0.5 20)) (void) (probcc-impossible))))
(⊦ equal?
`(((V ,(void)) ,(expt 2.0 -20)))
(probcc-sample-importance
(probcc-random-selector 17)
1
(if (ltrues 20 (lflips 0.5 20)) (void) (probcc-impossible))))
(⊦ equal?
`(((V #t) ,(expt 2.0 -20)))
(probcc-sample-importance
(probcc-random-selector 17)
1
(let1 (ts (lflips 0.5 20))
(if (ltrues 20 ts)
((probcc-lhead (τ ts)))
(probcc-impossible))))))
((eta 0.13800000000000001)
(memory #(8396592 3306592 1048576))
(stdout "")
(stderr ""))
(aux hansei): logic programming and reversible parsers
Logic programming
After HANSEI as a Declarative Logic Programming Language, the logic programs of the page, of its files and of the parser combinators `Reversible' parser combinators.
Re-thinking Prolog
After Re-thinking Prolog (LG l.33-75) and its paper (RT):
❝The library interface the paper relies on (RT Figure 1, l.285-308) is
One lesson is that guessing-- non-determinism-- is fundamental, but should not be the default mode of execution. One should guess, but guess lazily.
dist, fail, reify0, letlazy, flip, uniformly, exact_reify and reify_part, here probcc-distribution, probcc-impossible, probcc-reify/0, probcc-letlazy (and probcc-letlazy/nesting), probcc-coin, probcc-uniform/either, probcc-reify/exact and probcc-reify/part. The cheap non-determinism of RT §3.1, a generator
that guesses whole lists before testing them, finds the splits of a list only within a
bound, and diverges past it.No logic variables? No unification?
After No logic variables? No unification? (LG l.230-327): a logic variable is a generator of its domain, the sharing of the generator expresses its repeated occurrences, and lazy evaluation delays the generation as unification would:
❝
It appears that lazy evaluation can make up for logic variables, rendering Herbrand's original approach practical.
The append relation
After The append relation (LG l.328-624) and logic_vars.ml: the Prolog relation
"append([],L,L).\nappend([H|T],L,[H|R]) :- append(T,L,R)."
"del([H|T],H,T).\ndel([H|T],D,[H|R]) :- del(T,D,R).\nprog3(L,S) :- del(L,true,LRem) -> S=deleted(LRem); S=failed.\n%% ?- prog3(L,S).\n%% L = [true|_G16],\n%% S = deleted(_G16)."
Solving scheduling problems
After Solving scheduling problems (LG l.625-757): state the constraints as evidence and ask for the worlds consistent with them; a world is rejected as soon as one fact contradicts the evidence, so only a part of the 8957952 worlds is generated.
Zebra puzzle
After Zebra puzzle (LG l.758-792): the all-different constraint of the five attributes is a lazy
permutation of (aux hansei lazy), and the puzzle a naive
generate-and-test; its unique solution, at weight 1/199065600000, is in the long suite hansei-long.scm. The U2 bridge puzzle (puzzle.ml) and type inference run backwards (type_inference.ml) are two more programs that guess.
Reversible parser combinators
After `Reversible' parser combinators (LG l.793-996) and parsec.ml, the parser combinators of (aux hansei parsec) work over a non-deterministic stream:
❝
They can emulate not only parsing of a stream (that is, parsing forwards) but also determining a stream for the given parse (parsing backwards).
A stream is a thunk answering '() or a pair (char . stream), and a parser maps a stream to the pair (v . stream) of its value and the rest. Alternation is a uniform
choice, so a parse comes with the probability of the choices that led to it. parsec-run/backward generates streams with parsec-stream-over, lazily, each cell a probcc-letlazy/nesting, and keeps those the parser accepts whole.
Maximal munch, committed choice and nested inference
After Maximal munch, committed choice and nested inference (LG l.997-1230): parsec-many is Prolog's
"many(P,S,Rest) :- call(P,S,Srem) *-> many(P,Srem,Rest) ; S = Rest."
probcc-soft-cut, a nested inference: probcc-first-success over the reified parser. Unlike Prolog's, it
keeps parsing backwards complete, since an inner reification lets the stream's choices
take effect in the outer world:❝
many p (s1 | s2) = many p s1 | many p s2
test/append: pass
The append relation run forwards (LV:39-163, LG #append): a logic variable is a generator of its domain, a-boolean or a-blist, made a variable by probcc-letlazy/nesting: once its value is asked for, it stays
the same in its world. Without it the list [X] is drawn anew each time it is forced, and ?- bool(X), append([X],[f,f],R) has four answers where Prolog
has two. Explored to depth 5 the generators give both kinds of answers that Prolog's
depth-first search misses, lists with an f after the three t's among them.
(define (test/append _)
(⊦ equal?
'(((#t #t #t #f #f) 1))
(probcc-reify/part #f (τ (blist->list (bl-append bl-t3 bl-f2)))))
(⊦ equal?
'((((#t) (#t #f #f)) 1/4)
(((#t) (#f #f #f)) 1/4)
(((#f) (#t #f #f)) 1/4)
(((#f) (#f #f #f)) 1/4))
(probcc-reify/part
#f
(τ (let* ((l (τ (let1 (h (probcc-coin 1/2)) (cons h bl-nil))))
(r (blist->list (bl-append l bl-f2)))
(x (blist->list l)))
(list x r)))))
(⊦ equal?
'((() 1/2)
((#t) 1/8)
((#f) 1/8)
((#t #t) 1/32)
((#t #f) 1/32)
((#f #t) 1/32)
((#f #f) 1/32))
(probcc-reify/part 5 (τ (blist->list (a-blist)))))
(⊦ equal?
'((((#t) (#t #f #f)) 1/2) (((#f) (#f #f #f)) 1/2))
(probcc-reify/part
#f
(τ (let* ((l (probcc-letlazy/nesting
(τ (let1 (h (probcc-coin 1/2)) (cons h bl-nil)))))
(r (blist->list (bl-append l bl-f2)))
(x (blist->list l)))
(list x r)))))
(⊦ equal?
'(((#t #t #t) 1/2)
((#t #t #t #t) 1/8)
((#t #t #t #f) 1/8)
((#t #t #t #t #t) 1/32)
((#t #t #t #t #f) 1/32)
((#t #t #t #f #t) 1/32)
((#t #t #t #f #f) 1/32))
(probcc-reify/part 5 (τ (blist->list (bl-append bl-t3 (a-blist))))))
(⊦ equal?
'(((#f #f) 1/2)
((#t #f #f) 1/8)
((#f #f #f) 1/8)
((#t #t #f #f) 1/32)
((#t #f #f #f) 1/32)
((#f #t #f #f) 1/32)
((#f #f #f #f) 1/32))
(probcc-reify/part 5 (τ (blist->list (bl-append (a-blist) bl-f2))))))
((eta 0.001) (memory #(12582912 3460688 1048576)) (stdout "") (stderr ""))
test/append/backwards: pass
The append relation run backwards (LV:165-268): the prefix, the suffix, all the
splits of [t,t,t,f,f] and its last element are found by
generate-and-test, bl-compare forcing the variables only as far as
the comparison needs, so the search is finite and exhaustive; the weights estimate
the size of the search space (e.g. 2^-12 for each of the six splits, 2^-22 for two
variables equal to the five-element list). With the comparisons swapped, x = y is asked first, over two unbounded lists, and the search
must be bounded (depth 30).
(define (test/append/backwards _)
(define (unappend x r)
(unless (bl-compare r bl-t3f2) (probcc-impossible))
(blist->list x))
(⊦ equal?
'(((#t #f #f) 1/128))
(probcc-reify/part
#f
(τ (let1 (x (a-blist))
(unappend x (bl-append (bl-cons #t (bl-cons #t bl-nil)) x))))))
(⊦ equal?
'(((#t #t #t) 1/128))
(probcc-reify/part
#f
(τ (let1 (x (a-blist))
(unappend x (bl-append x (bl-cons #f (bl-cons #f bl-nil))))))))
(⊦ equal?
(probcc-map
(λ (e) (list e 1/4096))
'(((#t #t #t #f #f) ())
((#t #t #t #f) (#f))
((#t #t #t) (#f #f))
((#t #t) (#t #f #f))
((#t) (#t #t #f #f))
(() (#t #t #t #f #f))))
(probcc-reify/part
#f
(τ (let* ((x (a-blist)) (y (a-blist)))
(unless (bl-compare (bl-append x y) bl-t3f2) (probcc-impossible))
(let* ((ly (blist->list y)) (lx (blist->list x))) (list lx ly))))))
(⊦ equal? '((#f 1/1024)) (probcc-reify/part #f (τ ((bl-last bl-t3f2)))))
(⊦ equal?
'(((#t #t #t #f #f) 1/4194304))
(probcc-reify/part
#f
(τ (let* ((x (a-blist)) (y (a-blist)))
(unless (bl-compare x bl-t3f2) (probcc-impossible))
(unless (bl-compare x y) (probcc-impossible))
(blist->list y)))))
(⊦ equal?
'(((#t #t #t #f #f) 1/4194304))
(probcc-reify/part
30
(τ (let* ((x (a-blist)) (y (a-blist)))
(unless (bl-compare x y) (probcc-impossible))
(unless (bl-compare x bl-t3f2) (probcc-impossible))
(blist->list y))))))
((eta 0.023) (memory #(12582912 3935776 1048576)) (stdout "") (stderr ""))
test/committed-choice: pass
Deletion and committed choice (LV:270-444): bl-del deletes the first #t of a list, failing on a list without
one; run backwards over a-blist it generates the lists that have
one. Prolog's committed choice
"prog3(L,S) :- del(L,true,LRem) -> S=deleted(LRem) ; S=failed."
?- prog3(L,S). but the first (PD:202-217);
with the soft-cut probcc-soft-cut, a nested inference, bl-prog3 keeps them all, S=failed (#f) for the lists without #t, and terminates with no bound
over the lists of three elements.(define (test/committed-choice _)
(define (with-remainder l)
(let* ((r (blist->list (bl-del l #t))) (x (blist->list l))) (list x r)))
(⊦ equal?
'(((#t #t #f #f) 1))
(probcc-reify/part #f (τ (blist->list (bl-del bl-t3f2 #t)))))
(⊦ equal?
'(((#t #t #t #f) 1))
(probcc-reify/part #f (τ (blist->list (bl-del bl-t3f2 #f)))))
(⊦ equal?
'((((#t) ()) 1/4) (((#f) ()) 1/4))
(probcc-reify/part
#f
(τ (with-remainder (τ (let1 (h (probcc-coin 1/2)) (cons h bl-nil)))))))
(⊦ equal?
'((((#t) ()) 1/2))
(probcc-reify/part
#f
(τ (with-remainder
(probcc-letlazy/nesting
(τ (let1 (h (probcc-coin 1/2)) (cons h bl-nil))))))))
(⊦ equal?
'((((#t) ()) 1/8)
(((#t #t) (#t)) 1/32)
(((#t #f) (#f)) 1/32)
(((#f #t) (#f)) 1/32))
(probcc-reify/part 5 (τ (with-remainder (a-blist)))))
(⊦ equal?
'((((#t #t #t) (#t #t)) 1/128)
(((#t #t #f) (#t #f)) 1/128)
(((#t #f #t) (#f #t)) 1/128)
(((#t #f #f) (#f #f)) 1/128)
(((#f #t #t) (#f #t)) 1/128)
(((#f #t #f) (#f #f)) 1/128)
(((#f #f #t) (#f #f)) 1/128))
(probcc-reify/part
#f
(τ (let1 (l (a-blist)) (bl-haslen l 3) (with-remainder l)))))
(⊦ equal?
'((((#f #f)) 1))
(probcc-reify/part
#f
(τ (bl-prog3/value (bl-cons #f (bl-cons #t (bl-cons #f bl-nil)))))))
(⊦ equal?
'(((() #f) 1/2)
(((#t) (())) 1/8)
(((#f) #f) 1/8)
(((#t #t) ((#t))) 1/32)
(((#t #f) ((#f))) 1/32)
(((#f #t) ((#f))) 1/32)
(((#f #f) #f) 1/32)
(((#t #t #t) ((#t #t))) 1/128)
(((#t #t #f) ((#t #f))) 1/128)
(((#t #f #t) ((#f #t))) 1/128)
(((#t #f #f) ((#f #f))) 1/128)
(((#f #t #t) ((#f #t))) 1/128)
(((#f #t #f) ((#f #f))) 1/128)
(((#f #f #t) ((#f #f))) 1/128)
(((#f #f #f) #f) 1/128)
(((#f #f #f #f) #f) 1/512))
(probcc-reify/part
9
(τ (let* ((l (a-blist)) (r (bl-prog3/value l)))
(list (blist->list l) r)))))
(⊦ equal?
'((((#t #t #t) ((#t #t))) 1/128)
(((#t #t #f) ((#t #f))) 1/128)
(((#t #f #t) ((#f #t))) 1/128)
(((#t #f #f) ((#f #f))) 1/128)
(((#f #t #t) ((#f #t))) 1/128)
(((#f #t #f) ((#f #f))) 1/128)
(((#f #f #t) ((#f #f))) 1/128)
(((#f #f #f) #f) 1/128))
(probcc-reify/part
#f
(τ (let1 (l (a-blist))
(bl-haslen l 3)
(let1 (r (bl-prog3/value l)) (list (blist->list l) r)))))))
((eta 0.002) (memory #(12582912 3461752 1048576)) (stdout "") (stderr ""))
test/committed-choice/letlazy: pass
Why LV:19 makes letlazy letlazy_nesting: the soft-cut of bl-prog3 reifies bl-del in an inner world, which starts with an empty memory (N:75-81); with the plain probcc-letlazy the list is drawn anew inside, unrelated to the
outer one, and the answers are wrong, e.g. the empty list has a remainder. The table
is upstream's under the same change.
(define (test/committed-choice/letlazy _)
(⊦ equal?
'(((() (())) 1/16)
((() ((#f))) 1/32)
(((#t) (())) 1/64)
(((#f) (())) 1/64)
((() ((#t))) 1/64)
((() ((#f #f))) 3/256)
(((#t) ((#f))) 1/128)
(((#f) ((#f))) 1/128)
((() ((#f #t))) 1/128)
(((#t #t) (())) 1/256)
(((#t #f) (())) 1/256)
(((#t) ((#t))) 1/256)
(((#f #t) (())) 1/256)
(((#f #f) (())) 1/256)
(((#f) ((#t))) 1/256)
((() ((#t #t))) 1/256)
((() ((#t #f))) 1/256))
(probcc-reify/part
9
(τ (let* ((l (a-blist/with probcc-letlazy)) (r (bl-prog3/value l)))
(list (blist->list l) r))))))
((eta 0.002) (memory #(12582912 3672920 1048576)) (stdout "") (stderr ""))
test/cheap-nondeterminism: pass
Cheap non-determinism (RT §3.1, p.4-5 l.318-392): append run backwards by
generate-and-test over ordinary lists, rt-a-list drawing a
whole list before append looks at it; the six splits of [t,t,t,f,f], at 2^-12 each (the paper's rounded 0.0002), are
found at depth 15 and at the paper's 25, where the generator has made many more
candidates to reject; at depth 35 the search practically diverges (RT l.390-392).
Compare the lazy a-blist above, which needs no bound.
(define (test/cheap-nondeterminism _)
(define (split bound)
(probcc-reify/part
bound
(τ (let* ((x (rt-a-list)) (y (rt-a-list)))
(unless (equal? (append x y) '(#t #t #t #f #f)) (probcc-impossible))
(list x y)))))
(define splits
(probcc-map
(λ (e) (list e 1/4096))
'(((#t #t #t #f #f) ())
((#t #t #t #f) (#f))
((#t #t #t) (#f #f))
((#t #t) (#t #f #f))
((#t) (#t #t #f #f))
(() (#t #t #t #f #f)))))
(⊦ equal? splits (split 15))
(⊦ equal? splits (split 25)))
((eta 0.73599999999999999)
(memory #(50331648 5019768 1048576))
(stdout "")
(stderr ""))
test/daily-schedule: pass
Mikael More's scheduling problem (DS, LG #daily-schedule): the actions of each day are a memoized subset of the four, a logic variable,
and every constraint a fact about a uniformly chosen day. The schedule is unique: on
Monday the supermarket and the shop, on Tuesday the walk and the barber; its weight,
1/8957952 (upstream 1.11632658893461385e-07), estimates the search space, of which only a part is visited, since a world
is rejected as soon as one of its facts contradicts the evidence.
(define (test/daily-schedule _)
(⊦ equal?
(probcc-map
(λ (s) `((V ,s) 1/8))
'((3 2 1) (3 2) (3 1) (3) (2 1) (2) (1) ()))
(probcc-reify/exact (ds-powerset '(1 2 3))))
(⊦ equal?
'(((V #((Supermarket Shop) (Barber Walk) ())) 1/8957952))
(probcc-reify/exact (ds-schedule-model))))
((eta 0.056000000000000001)
(memory #(50331648 3466632 1048576))
(stdout "")
(stderr ""))
test/bridge: pass
The U2 bridge puzzle (PZ:1-185): four people, one flashlight, at most two on
the bridge, 17 minutes. With HANSEI's SimpleNonDet (choose a uniform choice, failing on the empty list) the
exhaustive search finds the two known solutions, each crossing a pair (party forward): Bono and Edge cross, one of them returns, Adam and Larry cross, the other
returns, Bono and Edge cross. The weights, upstream's under OCaml (1/6912 and
1/13824), are those of the choices that led to each trace.
(define (test/bridge _)
(⊦ equal?
'(((V (((Bono Edge) #t)
((Edge) #f)
((Adam Larry) #t)
((Bono) #f)
((Bono Edge) #t)))
1/6912)
((V (((Bono Edge) #t)
((Bono) #f)
((Adam Larry) #t)
((Edge) #f)
((Bono Edge) #t)))
1/13824))
(probcc-reify/exact (pz-loop '() #t 17 (list pz-u2 '())))))
((eta 0.222) (memory #(50331648 3467816 1048576)) (stdout "") (stderr ""))
test/type-inference/types: pass
Types read on demand (TI:81-199): a type is a procedure asked either for its value
or to equal a given value, and a type variable, ti-tp-memo, keeps in a world-local reference (probcc-wref-new, deviation 17) the type it has become, so that it
is the same whenever it is asked. Equating two types determines the first only as
far as the second needs. WorldLocalMem is not shipped upstream; these tables are TI's, and are also what TI gives under
OCaml with wref_new, wref_get and wref_set written over ProbM's memory.
(define (test/type-inference/types _)
(define (tp-and-arrow both?)
(let* ((tp (ti-tp-memo ti-a-tp))
(t2 (ti-tp-memo ti-a-tp))
(t1 (ti-tp-memo ti-a-tp))
(tpa (ti-arr t1 t2)))
(ti-tp-same tp tpa)
(when both? (ti-tp-same tpa tp))
(let* ((b (ti-string-of-tp tpa)) (a (ti-string-of-tp tp))) (list a b))))
(⊦ equal?
'(("int" 1/2) ("int -> int" 1/8))
(probcc-reify/part 3 (τ (ti-string-of-tp ti-a-tp))))
(⊦ equal?
`(((V ,(void)) 1/2))
(probcc-explore 3 (probcc-reify/0 (τ (ti-tp-same ti-a-tp ti-int)))))
(⊦ equal?
`(((V ,(void)) 1/2))
(probcc-explore
3
(probcc-reify/0 (τ (ti-tp-same ti-a-tp (ti-arr ti-int ti-int))))))
(⊦ equal?
'((#t 11/16))
(probcc-reify/part
5
(τ (let1 (tp (ti-tp-memo ti-a-tp))
(let* ((b (ti-string-of-tp tp)) (a (ti-string-of-tp tp)))
(equal? a b))))))
(⊦ equal?
'((("int -> int" "int -> int") 1/8)
(("int -> int -> int" "int -> int -> int") 1/32)
(("(int -> int) -> int" "(int -> int) -> int") 1/32))
(probcc-reify/part 5 (τ (tp-and-arrow #f))))
(⊦ equal?
'((#t 3/16))
(probcc-reify/part
5
(τ (match/first (tp-and-arrow #t) ((,a ,b) (equal? a b))))))
(⊦ equal?
'((#t 1/2))
(probcc-reify/part
5
(τ (let* ((tpa (ti-arr ti-int (ti-tp-memo ti-a-tp)))
(int->int (ti-arr ti-int ti-int))
(tpb (ti-arr (ti-tp-memo ti-a-tp) int->int)))
(ti-tp-same tpa tpb)
(let* ((b (ti-string-of-tp tpb)) (a (ti-string-of-tp tpa)))
(equal? a b))))))
(⊦ equal?
`((,(void) 1/2))
(probcc-reify/part
5
(τ (let* ((tpa (ti-tp-memo ti-a-tp)) (tpb (ti-arr tpa ti-int)))
(ti-tp-same tpa tpb)))))
(⊦ equal?
'()
(probcc-reify/part
#f
(τ (let* ((tpa (ti-tp-memo ti-a-tp)) (tpb (ti-arr tpa ti-int)))
(ti-tp-same tpa tpb)
(ti-tp-same ti-int tpa))))))
((eta 0.002) (memory #(50331648 3466768 1048576)) (stdout "") (stderr ""))
test/type-inference: pass
Type inference for the simply-typed lambda calculus (TI:209-263), forwards: ti-typeof gives every subterm a fresh type variable and equates
types as the rules ask; an ill-typed term has no type, with no bound, and a
polymorphic one has a stream of candidates, its instances, of which a bounded
search shows the first. There is no occurs check, so Lx. x x gets cyclic types, printed up to depth 10.
(define (test/type-inference _)
(⊦ equal?
'(("int -> int" 1/2)
("(int -> int) -> int" 1/8)
("(int -> int -> int) -> int" 1/32)
("((int -> int) -> int) -> int" 1/32))
(ti-infer 4 (ti-lam "f" (ti-num 1))))
(⊦ equal? '() (ti-infer #f (ti-lam "f" ti-vx)))
(⊦ equal?
'(("int -> int" 1/2)
("(int -> int) -> int -> int" 1/8)
("(int -> int -> int) -> int -> int -> int" 1/32)
("((int -> int) -> int) -> (int -> int) -> int" 1/32))
(ti-infer 4 (ti-lam "f" ti-vf)))
(⊦ equal?
'(("int" 1/2))
(ti-infer #f (ti-app (ti-lam "f" ti-vf) (ti-num 1))))
(⊦ equal?
'(("(int -> int) -> int" 1/4)
("(int -> int -> int) -> int -> int" 1/16)
("(int -> int -> int -> int) -> int -> int -> int" 1/64)
("(int -> (int -> int) -> int) -> (int -> int) -> int" 1/64))
(ti-infer 5 (ti-lam "f" (ti-app ti-vf (ti-num 1)))))
(⊦ equal?
'(("((((((((((... -> int) -> int) -> int) -> int) -> int) -> int) -> int) -> int) -> int) -> int) -> int"
1/4)
("((((((((((... -> ...) -> int -> int) -> int -> int) -> int -> int) -> int -> int) -> int -> int) -> int -> int) -> int -> int) -> int -> int) -> int -> int) -> int -> int"
1/16))
(ti-infer 3 (ti-lam "x" (ti-app ti-vx ti-vx))))
(⊦ equal?
(append
'(("(int -> int) -> (int -> int) -> int" 1/16)
("(int -> int -> int) -> (int -> int) -> int -> int" 1/64)
("((int -> int) -> int) -> (int -> int -> int) -> int" 1/64))
(probcc-map
(λ (s) (list s 1/256))
'("(int -> int -> int -> int) -> (int -> int) -> int -> int -> int"
"(int -> (int -> int) -> int) -> (int -> int) -> (int -> int) -> int"
"((int -> int) -> int -> int) -> (int -> int -> int) -> int -> int"
"((int -> int -> int) -> int) -> (int -> int -> int -> int) -> int"
"(((int -> int) -> int) -> int) -> (int -> (int -> int) -> int) -> int"))
(probcc-map
(λ (s) (list s 1/1024))
'("(int -> int -> int -> int -> int) -> (int -> int) -> int -> int -> int -> int"
"(int -> int -> (int -> int) -> int) -> (int -> int) -> int -> (int -> int) -> int"
"(int -> (int -> int) -> int -> int) -> (int -> int) -> (int -> int) -> int -> int"
"(int -> (int -> int -> int) -> int) -> (int -> int) -> (int -> int -> int) -> int"
"(int -> ((int -> int) -> int) -> int) -> (int -> int) -> ((int -> int) -> int) -> int"
"((int -> int) -> int -> int -> int) -> (int -> int -> int) -> int -> int -> int"
"((int -> int) -> (int -> int) -> int) -> (int -> int -> int) -> (int -> int) -> int"
"((int -> int -> int) -> int -> int) -> (int -> int -> int -> int) -> int -> int"
"((int -> int -> int -> int) -> int) -> (int -> int -> int -> int -> int) -> int"
"((int -> (int -> int) -> int) -> int) -> (int -> int -> (int -> int) -> int) -> int"
"(((int -> int) -> int) -> int -> int) -> (int -> (int -> int) -> int) -> int -> int"
"(((int -> int) -> int -> int) -> int) -> (int -> (int -> int) -> int -> int) -> int"
"(((int -> int -> int) -> int) -> int) -> (int -> (int -> int -> int) -> int) -> int"
"((((int -> int) -> int) -> int) -> int) -> (int -> ((int -> int) -> int) -> int) -> int")))
(ti-infer
10
(ti-lam "f" (ti-lam "x" (ti-app ti-vf (ti-app ti-vx (ti-num 1)))))))
(⊦ equal?
'()
(ti-infer
#f
(ti-lam
"f"
(ti-lam
"x"
(ti-app
(ti-app ti-vf (ti-app ti-vx (ti-num 1)))
(ti-app ti-vx (ti-lam "z" ti-vz))))))))
((eta 0.0080000000000000002)
(memory #(50331648 3476424 1048576))
(stdout "")
(stderr ""))
test/type-inference/backwards: pass
Type inference backwards (TI:265-310): ti-a-term generates terms read on demand, and asking for their types, or for the type int -> int, keeps the well-typed ones: a type checker run as a
generator of programs of a given type. The masses are the exact rationals of
upstream's (e.g. 1/500 for 0.00200000000000000048).
(define (test/type-inference/backwards _)
(define (generate bound)
(probcc-reify/part
bound
(τ (let* ((term (ti-a-term))
(tp (ti-string-of-tp (ti-typeof '() term)))
(s (ti-string-of-term term)))
(list s tp)))))
(⊦ equal?
'(("1" 1/10)
("x0" 1/100)
("x1" 9/1000)
("x2" 81/10000)
("x3" 729/100000)
("Lx0. 1" 1/250)
("1 1" 1/250)
("Lx1. 1" 9/2500)
("Lx2. 1" 81/25000)
("x0 1" 1/2500)
("Lx0. x0" 1/2500)
("1 x0" 1/2500)
("x1 1" 9/25000)
("1 x1" 9/25000)
("Lx0. 1 1" 1/6250)
("1 1 1" 1/6250)
("1 (Lx0. 1)" 1/6250)
("1 (1 1)" 1/6250)
("(Lx0. 1) 1" 1/6250))
(probcc-reify/part 5 (τ (ti-string-of-term (ti-a-term)))))
(⊦ equal? '((("1" "int") 1/10) (("Lx0. 1" "int -> int") 1/500)) (generate 4))
(⊦ equal?
'((("1" "int") 1/10)
(("Lx0. 1" "int -> int") 1/500)
(("Lx1. 1" "int -> int") 9/5000))
(generate 5))
(⊦ equal?
'((("1" "int") 1/10)
(("Lx0. 1" "int -> int") 1/500)
(("Lx1. 1" "int -> int") 9/5000)
(("Lx2. 1" "int -> int") 81/50000)
(("Lx3. 1" "int -> int") 729/500000)
(("Lx0. 1" "(int -> int) -> int") 1/2000)
(("Lx1. 1" "(int -> int) -> int") 9/20000)
(("Lx0. x0" "int -> int") 1/5000)
(("(Lx0. 1) 1" "int") 1/12500)
(("(Lx1. 1) 1" "int") 9/125000))
(generate 7))
(⊦ equal?
`((("Lx0. 1" ,(void)) 1/500)
(("Lx1. 1" ,(void)) 9/5000)
(("Lx2. 1" ,(void)) 81/50000)
(("Lx3. 1" ,(void)) 729/500000)
(("Lx4. 1" ,(void)) 6561/5000000)
(("Lx5. 1" ,(void)) 59049/50000000)
(("Lx6. 1" ,(void)) 531441/500000000)
(("Lx0. x0" ,(void)) 1/5000)
(("Lx1. x1" ,(void)) 81/500000)
(("Lx2. x2" ,(void)) 6561/50000000)
(("Lx0. (Lx0. 1) 1" ,(void)) 1/625000)
(("(Lx0. Lx0. 1) 1" ,(void)) 1/625000))
(probcc-reify/part
10
(τ (let* ((term (ti-a-term))
(u (ti-tp-same (ti-arr ti-int ti-int) (ti-typeof '() term)))
(s (ti-string-of-term term)))
(list s u))))))
((eta 0.073999999999999996)
(memory #(50331648 4916688 1048576))
(stdout "")
(stderr ""))
test/parsec/kit: pass
The streams and the primitive parsers (PS:24-144, 278-290): parsec-or is a uniform choice of two parsers, nested to the left
for more, parsec-alt of a list of them; a parse that leaves some
input, or fails, gives no world. Explored to depth 5, parsec-stream-over gives upstream's table of PS:272-273, a fair
coin and a uniform character per cell (deviation 22: the unit value is (void), the alphabet a list).
(define (test/parsec/kit _)
(define abc (parsec-stream-of-string "abc"))
(⊦ equal? "abc" (parsec-string-of-stream abc))
(⊦ equal? "" (parsec-string-of-stream (parsec-stream-of-string "")))
(⊦ equal? #\a (car (abc)))
(⊦ equal? "bc" (parsec-string-of-stream (cdr (abc))))
(⊦ equal? `(((V ,(void)) 1)) (parsec-run/forward parsec-eof ""))
(⊦ equal? '() (parsec-run/forward parsec-eof "a"))
(⊦ equal? `(((V ,(void)) 1)) (parsec-run/forward parsec-empty ""))
(⊦ equal? '(((V x) 1)) (parsec-run/forward (parsec-pure 'x) ""))
(⊦ equal? '() (parsec-run/forward (parsec-pure 'x) "a"))
(⊦ equal? '(((V #\a) 1)) (parsec-run/forward p-a "a"))
(⊦ equal? '() (parsec-run/forward p-a "b"))
(⊦ equal? '() (parsec-run/forward p-a ""))
(⊦ equal? '() (parsec-run/forward p-a "ab"))
(⊦ equal? '(((V #\7) 1)) (parsec-run/forward (parsec-sat char-numeric?) "7"))
(⊦ equal? '(((V 7) 1)) (parsec-run/forward digit "7"))
(⊦ equal?
'(((V (#\a #\b)) 1))
(parsec-run/forward (parsec-ap (parsec-fmap pair p-a) p-b) "ab"))
(⊦ equal? '(((V #\b) 1)) (parsec-run/forward (parsec-right p-a p-b) "ab"))
(⊦ equal? '(((V #\a) 1)) (parsec-run/forward (parsec-left p-a p-b) "ab"))
(⊦ equal?
'(((V #\a) 1/2))
(probcc-reify/exact
(car ((parsec-or p-a p-b) (parsec-stream-of-string "ab")))))
(⊦ equal?
'(((V #\a) 1/4))
(parsec-run/forward (parsec-or (parsec-or p-a p-b) p-b) "a"))
(⊦ equal?
'(((V #\a) 1/3))
(parsec-run/forward (parsec-alt (list p-a p-b p-b)) "a"))
(⊦ equal? '() (parsec-run/forward (parsec-alt '()) ""))
(⊦ equal?
'(((V 3) 1/16))
(parsec-run/forward
(parsec-fix
(λ (p)
(parsec-or
(parsec-fmap (λ (_) 0) parsec-eof)
(parsec-right p-a (parsec-fmap add1 p)))))
"aaa"))
(⊦ equal? `(((V ,(void)) 1)) (probcc-reify/exact (parsec-stream-len abc 3)))
(⊦ equal? '() (probcc-reify/exact (parsec-stream-len abc 2)))
(⊦ equal? '() (probcc-reify/exact (parsec-stream-len abc 4)))
(⊦ equal?
`(((V ,(void)) 1))
(probcc-reify/exact (parsec-stream-at-most-len abc 3)))
(⊦ equal?
`(((V ,(void)) 1))
(probcc-reify/exact (parsec-stream-at-most-len abc 4)))
(⊦ equal? '() (probcc-reify/exact (parsec-stream-at-most-len abc 2)))
(⊦ equal?
`(((V "") 1/2)
((V "b") 1/8)
((V "a") 1/8)
((V "bb") 1/32)
((V "ba") 1/32)
((V "ab") 1/32)
((V "aa") 1/32)
,@(make-list 16 '(C 1/128)))
(shape (probcc-explore
5
(probcc-reify/0
(τ (parsec-string-of-stream
(parsec-stream-over '(#\a #\b)))))))))
((eta 0.001) (memory #(50331648 3482112 1048576)) (stdout "") (stderr ""))
test/parsec/palindromes: pass
Recognizing palindromes forwards (PS:188-251), with the left-nested binary parsec-or and with parsec-alt of the five alternatives; the probabilities are the exact rationals of upstream's
floats, e.g. 1/25 for 0.0399999999999999939 (PS:241).
(define (test/parsec/palindromes _)
(⊦ equal? '() (parsec-run/forward pali "ab"))
(⊦ equal? `(((V ,(void)) 1/64)) (parsec-run/forward pali "aa"))
(⊦ equal? `(((V ,(void)) 1/256)) (parsec-run/forward pali "aaaa"))
(⊦ equal? `(((V ,(void)) 1/32)) (parsec-run/forward pali "aba"))
(⊦ equal? `(((V ,(void)) 1/512)) (parsec-run/forward pali "abaaba"))
(⊦ equal? '() (parsec-run/forward pali/alt "ab"))
(⊦ equal? `(((V ,(void)) 1/25)) (parsec-run/forward pali/alt "aa"))
(⊦ equal? `(((V ,(void)) 1/125)) (parsec-run/forward pali/alt "aaaa"))
(⊦ equal? `(((V ,(void)) 1/25)) (parsec-run/forward pali/alt "aba"))
(⊦ equal? `(((V ,(void)) 1/625)) (parsec-run/forward pali/alt "abaaba")))
((eta 0.0) (memory #(50331648 3473696 1048576)) (stdout "") (stderr ""))
test/parsec/backward: pass
Parsing backwards (PS:253-310): parsec-run/backward generates the streams of parsec-stream-over and keeps those the parser accepts whole; bounded to depth 10 it finds the
palindromes of at most three characters, and with streams of exactly five
characters (parsec-stream-len) the exploration is exhaustive:
the eight palindromes, each at 1/256000 (upstream 3.90625000000000093e-06). The Curry counterpart is findall (\[x1,x2,x3,x4,x5] -> pali [x1,x2,x3,x4,x5] =:= []).
(define (test/parsec/backward _)
(⊦ equal?
`(((V (,(void) "")) 1/10)
((V (,(void) "b")) 1/40)
((V (,(void) "a")) 1/40)
((V (,(void) "bb")) 1/800)
((V (,(void) "aa")) 1/800)
((V (,(void) "bbb")) 1/3200)
((V (,(void) "bab")) 1/3200)
((V (,(void) "aba")) 1/3200)
((V (,(void) "aaa")) 1/3200))
(let1 (pv (parsec-run/backward
10
pali/alt
(τ (parsec-stream-over '(#\a #\b)))))
(filter (λ1-match/first (((V ,_) ,_) #t) (else #f)) pv)))
(⊦ equal?
(probcc-map
(λ (s) `((V (,(void) ,s)) 1/256000))
'("bbbbb" "bbabb" "babab" "baaab" "abbba" "ababa" "aabaa" "aaaaa"))
(parsec-run/backward
#f
pali/alt
(τ (let1 (st (parsec-stream-over '(#\a #\b)))
(parsec-stream-len st 5)
st)))))
((eta 0.002) (memory #(50331648 3483328 1048576)) (stdout "") (stderr ""))
test/parsec/many: pass
The tests of parsec-many (PS:316-395; the table of LG
l.1104-1146): it munches maximally, so it has one parse where a Kleene star has
several, and fails where a shorter munch would let the rest parse, e.g. many (p_char 'a') <* p_char 'a' on "aaa"; many/no-munch, upstream's commented-out many (PS:146-151), gives the parses of every munch instead. Parsed
backwards, parsec-many still generates its whole language.
(define (test/parsec/many _)
(define many-a (parsec-many p-a))
(⊦ equal? '(((V ()) 1)) (parsec-run/forward (parsec-left many-a p-b) "b"))
(⊦ equal? '(((V (#\a #\a #\a)) 1)) (parsec-run/forward many-a "aaa"))
(⊦ equal?
'(((V ((#\a #\a #\a) ())) 1))
(parsec-run/forward (parsec-ap (parsec-fmap pair many-a) many-a) "aaa"))
(⊦ equal? '() (parsec-run/forward (parsec-left many-a p-a) "aaa"))
(⊦ equal?
'(((V (#\a #\b #\a #\b #\b)) 1/32))
(parsec-run/forward (parsec-many (parsec-or p-a p-b)) "ababb"))
(⊦ equal?
'(((V ((#\a #\a #\a))) 1))
(parsec-run/forward (parsec-many (parsec-many1 p-a)) "aaa"))
(⊦ equal?
'(((V ((#\a #\a #\a) (#\b) (#\a) (#\b))) 1/16))
(parsec-run/forward
(parsec-many (parsec-or (parsec-many1 p-a) (parsec-many1 p-b)))
"aaabab"))
(⊦ equal?
'(((V ("aa" "a")) 1/4) ((V ("a" "aa")) 1/4) ((V ("a" "a" "a")) 1/8))
(parsec-run/forward
(parsec-left
(parsec-many
(parsec-or
(parsec-left (parsec-left (parsec-pure "aa") p-a) p-a)
(parsec-left (parsec-pure "a") p-a)))
p-b)
"aaab"))
(⊦ equal?
'()
(parsec-run/forward
(parsec-left (parsec-many (parsec-or (parsec-left p-a p-a) p-a)) p-a)
"aaa"))
(⊦ equal?
'(((V (#\a)) 1/8) ((V (#\a #\a)) 1/32))
(parsec-run/forward
(parsec-left (many/no-munch (parsec-or (parsec-left p-a p-a) p-a)) p-a)
"aaa"))
(⊦ equal?
'(((V (#\a #\a)) 1/8))
(parsec-run/forward (parsec-left (many/no-munch p-a) p-a) "aaa"))
(⊦ equal?
'(((V (#\a #\a #\a)) 1/16))
(parsec-run/forward (many/no-munch p-a) "aaa"))
(⊦ equal?
'(((V ((() #\b) "b")) 1/8)
((V (((#\a) #\b) "ab")) 1/32)
((V (((#\a #\a) #\b) "aab")) 1/128))
(parsec-run/backward
#f
(parsec-ap (parsec-fmap pair many-a) p-b)
(τ (let1 (st (parsec-stream-over '(#\a #\b)))
(parsec-stream-at-most-len st 3)
st))))
(⊦ equal?
'(((V (() "")) 1/2)
((V ((#\a) "a")) 1/8)
((V ((#\a #\a) "aa")) 1/32)
((V ((#\a #\a #\a) "aaa")) 1/128)
(C 1/256)
(C 1/256))
(shape (parsec-run/backward
9
many-a
(τ (parsec-stream-over '(#\a #\b)))))))
((eta 0.001) (memory #(50331648 3492088 1048576)) (stdout "") (stderr ""))
test/parsec/expression: pass
The arithmetic expressions of the Curry example (PS:397-466): num "123" is 123 for sure, expression "10" is 10 at 1/8 (three binary alternatives on the way to num), and "(10+5*2)/4" is 5 at 2^-16 (upstream 1.52587890625e-05).
(define (test/parsec/expression _)
(⊦ equal? '(((V 123) 1)) (parsec-run/forward num "123"))
(⊦ equal? '(((V 10) 1/8)) (parsec-run/forward expression "10"))
(⊦ equal? '(((V 5) 1/65536)) (parsec-run/forward expression "(10+5*2)/4"))
(⊦ equal? '(((V 1) 1/64)) (parsec-run/forward expression "7/4"))
(⊦ equal? '() (parsec-run/forward expression "1+")))
((eta 0.001) (memory #(25165824 3479096 1048576)) (stdout "") (stderr ""))
test/parsec/factor: pass
Generating expressions (PS:514-529): factor explored to depth 14 over the four-character streams of the balanced exp_stream (a digit with probability 3/10) has 17768 entries, 10100
of them values, as under OCaml; among them (9999 "9999") and (9998 "9998") at 81/200000000 (upstream 4.05000000000000305e-07).
(define (test/parsec/factor _)
(define pv
(parsec-run/backward
14
factor
(exp-stream/balanced 7/10 probcc-uniform/either)))
(define leaves (filter (λ1-match/first (((V ,_) ,_) #t) (else #f)) pv))
(⊦ equal? '(17768 10100) (list (length pv) (length leaves)))
(⊦ equal?
'(((V (9999 "9999")) 81/200000000) ((V (9998 "9998")) 81/200000000))
(filter
(λ1-match/first
(((V (,_ ,s)) ,_) (member s '("9999" "9998")))
(else #f))
leaves)))
((eta 0.26500000000000001)
(memory #(25165824 9629000 1048576))
(stdout "")
(stderr ""))
test/parsec/nested: pass
The nested-inference law of LG l.1212-1215: many p (s1 | s2) = many p s1 | many p s2. With s1 a random stream and s2 the stream of "aab", a stream chosen lazily inside the soft-cut of parsec-many (probcc-letlazy/nesting) gives the same tree, to depth 6, as a stream chosen before parsing: the inner
reification lets the stream's choice take effect in the outer world, as under
OCaml.
(define (test/parsec/nested _)
(define (many-over choose)
(probcc-explore
6
(probcc-reify/0
(τ (let* ((s1 (parsec-stream-over '(#\a #\b)))
(s2 (parsec-stream-of-string "aab"))
(x (choose s1 s2)))
(car ((parsec-many p-a) x)))))))
(define left
(many-over
(λ (s1 s2)
(probcc-letlazy/nesting (τ ((probcc-uniform/either (list s1 s2))))))))
(define right (many-over (λ (s1 s2) (probcc-uniform/either (list s1 s2)))))
(⊦ equal?
'(((V (#\a #\a)) 1/2) ((V ()) 3/8) ((V (#\a)) 3/32) (C 1/64) (C 1/64))
(shape left))
(⊦ equal? (shape left) (shape right)))
((eta 0.0) (memory #(25165824 9652640 1048576)) (stdout "") (stderr ""))
test/replay/parsec: pass
The sampling records of PS:492-501 and PS:531-546, importance sampling of 80
expressions of four characters and of 2000 over the more balanced streams, with (probcc-random-selector 17), replayed bit for bit: the characters are drawn with a transcription of upstream's uniform (P:156-165), since probcc-uniform's masses are exact (deviation 7), and the parsers are the module's. With parsec-stream-over, whose masses are exact, the draws are the same
and the two expressions come with the exact 1/2560. The runs of 30000 samples,
PS:504-512 and PS:549-557, are in the long suite.
(define (test/replay/parsec _)
(⊦ equal?
'(((V (7139 "7139")) 0.00039062499999999997)
((V (5164 "5164")) 0.00039062499999999997))
(probcc-sample-importance/a
(probcc-random-selector 17)
80
(expression-model (exp-stream ml-stream-over))))
(⊦ equal?
'(((V (7139 "7139")) 1/2560) ((V (5164 "5164")) 1/2560))
(probcc-sample-importance/a
(probcc-random-selector 17)
80
(expression-model (exp-stream parsec-stream-over))))
(⊦ equal?
(append
'(((V (3403 "3403")) 0.00050000000000000001)
((V (1396 "1396")) 0.00050000000000000001)
((V (55 "55+0")) 2.5000000000000022e-05))
(probcc-map
(λ (n)
`((V (,n ,(string-append "55+" (number->string (- n 55)))))
2.5000000000000001e-05))
'(64 63 62 61 60 59 58 57 56))
'(((V (92 "(92)")) 2.0000000000000005e-05)
((V (66 "(66)")) 2.0000000000000005e-05)))
(probcc-sample-importance/a
(probcc-random-selector 17)
2000
(expression-model
(exp-stream/balanced 0.69999999999999996 ml-uniformly)))))
((eta 0.17199999999999999)
(memory #(25165824 3529944 1048576))
(stdout "")
(stderr ""))
(aux hansei): example models
Introductory examples and regression tests
The example models of Embedded probabilistic domain-specific language HANSEI, one section per example of the page, with upstream's own expected values: the lawn
model of the papers, the regression tests of exactInfM.ml, which
liken HANSEI's syntax to IBAL's, and the early collapse of worlds:
❝
Naively, we would have thought we needed to evaluate 2^10 of possible worlds: each (nested) flip splits the world. But as the result shows, we managed with only 10+1 worlds. The reason is that the algorithm is NOT tail recursive, and && is a _special_ form rather than a function.
The models keep upstream's flonum probabilities, so that their tables are upstream's
bit for bit where the masses are summed in upstream's order; where a model draws with uniform, whose masses are the exact 1/n here (deviation 7), or reflects a table that probcc-explore has sorted by mass (deviation 5), the masses agree to a relative 1e-9.
Colored balls
❝
This is an elegant example of a model with indefinite number of objects whose `real' identity is unobservable. The example is a simplified version of population estimation. HANSEI's implementation of the example makes the case for memoization of stochastic functions.
Blood types
❝
The hidden variables are the presence of genes A and B in the chromosomes inherited from the mother and the father. The observed variable is the blood type, which is a stochastic function of the genotype. A typical query asks for the distribution of blood types for a particular individual, possibly conditioned on the known genotype of some of the ancestors.
Hidden Markov Model
Example: Hidden Markov Model (HMM):
❝
A simple HMM model with a single hidden 8-state variable. Our implementation demonstrates the investigation of time complexity of the inference -- in an interactive session at the OCaml top level -- and changing the model by introducing variable elimination, which speeds up inference exponentially.
Noisy OR
❝
Although the graph does not have directed cycles, it does have many undirected ones. As we work backwards, we may encounter a given node several times. Each time the value of node's random attribute must be the same. We have to use memoization.
Motivic development in classical music
Realistic Example: Motivic development in classical music, the warm-up of music1a.ml:
❝
Crucially this simplified model does not recursively subdivide and transform the sequence. Therefore, the state space remains sufficiently small even for exact inference. The lazy evaluation was essential.
Multitarget tracking
Realistic Example: Multitarget tracking:
❝
We found it crucial to sample the state of an aircraft and the location of a false alarm lazily. If the x-coordinate of the detected aircraft does not match the x-coordinate of any observed blip, we immediately report the contradiction (of the possible world with the evidence). We do not need to compute the y-coordinate of the aircraft, let alone the locations of other aircrafts and false alarms.
The ported models are the AILog one (blip_ailog.ml, memoized), the naive one without laziness (blip_naive.ml) and the lazy one with sequential Monte-Carlo (blip.ml).
Bounded-rational theory of mind
Bounded-Rational Theory of Mind for Conversational Implicature:
❝
The payoff is that the same cognitive machinery that the hearer uses to decide what to do, by interpreting utterances, can also be used by the speaker to decide what to say, by predicting how potential utterances would be interpreted or misinterpreted by the hearer.
The speaker's model reifies the hearer's, exactly or by sampling, inside its own
worlds (incite.ml); the weights of the plans are 1, so the masses
of the tables are expected utilities, not probabilities.
test/paper/grass-cloudy: pass
The extended lawn model of PE:9-28: whether it rained, given that the grass is wet, with a cloudy sky as the common cause of rain and sprinkler; the flip of the wet roof, whose value is not used, splits the worlds all the same. The table and its normalization are upstream's bit for bit.
(define (test/paper/grass-cloudy _)
(define table (probcc-reify/exact/a grass-model/cloudy))
(⊦ equal? '(((V #t) 0.45810000000000001) ((V #f) 0.18899999999999997)) table)
(⊦ equal?
'(((V #t) 0.70792767732962447) ((V #f) 0.29207232267037547))
(probcc-normalize table)))
((eta 0.001) (memory #(12582912 3994768 1048576)) (stdout "") (stderr ""))
test/exact/basic: pass
The basic tests of X:13-63: a shared flip and two independent ones, the alarm
network (whose let … and … draws the earthquake first, D19),
the lawn model and uniform_range 1 8, which is (probcc-uniform/range 1 9) (deviation 8). The tables are
upstream's bit for bit, and the number of leaves of each tree is the number of
worlds upstream's comments report as accepted.
(define (test/exact/basic _)
(⊦ equal? '(((V #t) 0.5) ((V #f) 0.5)) (probcc-reify/exact/a tflip2-shared))
(⊦ equal? 2 (accepted tflip2-shared))
(⊦ equal?
'(((V #f) 0.75) ((V #t) 0.25))
(probcc-reify/exact/a tflip2-noshared))
(⊦ equal? 3 (accepted tflip2-noshared))
(⊦ equal?
'(((V #f) 0.89128000000000007) ((V #t) 0.10872000000000001))
(probcc-reify/exact/a talarm))
(⊦ equal? 8 (accepted talarm))
(⊦ equal?
'(((V #f) 0.3219999999999999) ((V #t) 0.2838))
(probcc-reify/exact/a grass-model))
(⊦ equal? 14 (accepted grass-model))
(⊦ equal?
'(((V 8) 1/8)
((V 7) 1/8)
((V 6) 1/8)
((V 5) 1/8)
((V 4) 1/8)
((V 3) 1/8)
((V 2) 1/8)
((V 1) 1/8))
(probcc-reify/exact (probcc-uniform/range 1 9))))
((eta 0.0) (memory #(12582912 3994840 1048576)) (stdout "") (stderr ""))
test/exact/ibal: pass
The tests of IBAL (X:81-231), whose syntax HANSEI likens: observations through probcc-observe, records and objects as lists, the methods of an
object run when they are called; upstream's tables bit for bit.
(define (test/exact/ibal _)
(⊦ equal?
'(((V (#t #t)) 0.59999999999999998))
(probcc-reify/exact/a tfooibl1))
(⊦ equal?
'(((V (#\c #\e)) 0.47999999999999998)
((V (#\b #\e)) 0.32000000000000006)
((V (#\d #\e)) 0.20000000000000001))
(probcc-reify/exact/a tfooibl2))
(⊦ equal? '(((V #t) 0.25)) (probcc-reify/exact/a tfooibl3))
(⊦ equal?
'(((V #t) 0.021000000000000001) ((V #f) 0.0090000000000000011))
(probcc-reify/exact/a tiblt2))
(⊦ equal? '(((V #t) 0.14000000000000001)) (probcc-reify/exact/a tmusictr41))
(⊦ equal?
'(((V #f) 0.0069999999999999993) ((V #t) 0.0030000000000000001))
(probcc-reify/exact/a tmusictr42)))
((eta 0.0) (memory #(12582912 3996560 1048576)) (stdout "") (stderr ""))
test/exact/flips: pass
Early collapse of worlds (X:236-347): the and of ten flips needs 11 worlds, not
2^10, because and does not evaluate its second operand when the
first is false; the xor of ten flips needs all 1024, unless the xor of the others is
reified and reflected, collapsing its worlds into two, so that each level explores
4.
(define ((flips-xor n))
(let loop ((n n))
(if (= n 1)
(probcc-coin 0.5)
(let* ((r (loop (sub1 n))) (x (probcc-coin 0.5))) (not (eq? x r))))))
(define ((flips-xor* n))
(let loop ((n n))
(if (= n 1)
(probcc-coin 0.5)
(let* ((r (probcc-reflect (probcc-reify/exact (loop (sub1 n)))))
(x (probcc-coin 0.5)))
(not (eq? x r))))))
(define (test/exact/flips _)
(⊦ equal?
'(((V #f) 1023/1024) ((V #t) 1/1024))
(probcc-reify/exact/a (flips-and 1/2 10)))
(⊦ equal?
'(((V #f) 0.9990234375) ((V #t) 0.0009765625))
(probcc-reify/exact/a (flips-and 0.5 10)))
(⊦ equal? 11 (accepted (flips-and 1/2 10)))
(⊦ equal? '(((V #t) 0.5) ((V #f) 0.5)) (probcc-reify/exact/a (flips-xor 10)))
(⊦ equal? 1024 (accepted (flips-xor 10)))
(⊦ equal?
'(((V #t) 0.5) ((V #f) 0.5))
(probcc-reify/exact/a (flips-xor* 10)))
(⊦ equal? 4 (accepted (flips-xor* 10))))
((eta 0.0040000000000000001)
(memory #(12582912 4001992 1048576))
(stdout "")
(stderr ""))
test/paper/lawn: pass
The lawn model of the DSL paper (DSL p.4 l.160-166), normalized as on p.5
l.194-197, whose 14 accepted worlds are those of test/exact/basic; the paper's normalized masses were printed by an earlier version and differ
from today's upstream, which this port gives bit for bit, in the last two digits.
The two coins of the NIPS abstract (NIPS l.24-27, 63), at least one of them true,
with exact halves.
(define (grass-model)
(let* ((rain (probcc-coin 0.29999999999999999))
(sprinkler (probcc-coin 0.5))
(grass-is-wet
(or (and (probcc-coin 0.90000000000000002) rain)
(and (probcc-coin 0.80000000000000004) sprinkler)
(probcc-coin 0.10000000000000001))))
(probcc-when grass-is-wet rain)))
(define (test/paper/lawn _)
(define normalized (probcc-normalize (probcc-reify/exact/a grass-model)))
(⊦= '(((V #f) 0.53152855727963022) ((V #t) 0.46847144272036972)) normalized)
(⊦ equal?
'(((V #f) 0.53152855727963011) ((V #t) 0.46847144272036984))
normalized)
(⊦ equal?
'(((V (#t #t)) 1/4) ((V (#t #f)) 1/4) ((V (#f #t)) 1/4))
(probcc-reify/exact/a two-coins)))
((eta 0.001) (memory #(12582912 4004104 1048576)) (stdout "") (stderr ""))
test/noisy-or: pass
Noisy-or networks (NO): the attribute of a vertex propagated from the roots,
forwards with the graph built in, or backwards from the vertex asked, where a vertex
reached by several paths must have one value per world: memoized with fixm over probcc-memo, or with the evidence
threaded as the memo table, which also holds the observations (tq2, vertex 0 off). The large graph has 20 vertices, 0 to 19 (NO:187 says 19);
every table is upstream's bit for bit.
(define (test/noisy-or _)
(define preds (no-predecessors no-edges))
(define on3 '(((V #t) 0.81441875000000008) ((V #f) 0.18558125000000003)))
(⊦ equal? '(4 2 1) (pmap-find 3 preds))
(⊦ equal? on3 (probcc-reify/exact/a no-model-fwd))
(⊦ equal?
on3
(probcc-reify/exact (no-model/memo no-roots preds 3 (λ (_) #f))))
(⊦ equal? on3 (probcc-reify/exact (no-model no-roots preds 3 pmap-empty)))
(⊦ equal?
'(((V #t) 0.68646750000000001) ((V #f) 0.31353250000000005))
(probcc-reify/exact
(no-model no-roots preds 3 (pmap-add 0 #f pmap-empty))))
(⊦ equal?
'(((V #t) 0.90144617982013708) ((V #f) 0.098553820179875601))
(probcc-reify/exact
(no-model
no-roots/large
(no-predecessors no-edges/large)
15
pmap-empty))))
((eta 0.22600000000000001)
(memory #(12582912 4011768 1048576))
(stdout "")
(stderr ""))
test/bloodtype: pass
Blood types in a pedigree (BT:224-275), exactly, with variable elimination of the
parents' chromosomes: unconditionally, and given that Uwe's chromosomes are both A.
With the reflected tables in upstream's order the masses are upstream's bit for
bit; probcc-variable-elimination reflects them sorted by mass
(deviation 5), so they are summed in another order and agree to a relative
1e-9.
(define (test/bloodtype _)
(define test1exact
'(((V A) 0.32180904495702678)
((V B) 0.32085269277525014)
((V AB) 0.19786055151667717)
((V O) 0.15947771075104664)))
(define test2exact
'(((V A) 0.43868996766624768)
((V AB) 0.23650376976933171)
((V B) 0.21308545732967843)
((V O) 0.1117208052347424)))
(define uwe '(((Uwe GA)) ((Uwe GA))))
(⊦ equal?
test1exact
(probcc-reify/exact/a
(bt-bldtype variable-elimination/upstream '(() ()) 'Linus)))
(⊦ equal?
test2exact
(probcc-reify/exact/a
(bt-bldtype variable-elimination/upstream uwe 'Linus)))
(⊦ (rel-approx? 1.0000000000000001e-09)
test1exact
(probcc-reify/exact/a
(bt-bldtype probcc-variable-elimination '(() ()) 'Linus)))
(⊦ (rel-approx? 1.0000000000000001e-09)
test2exact
(probcc-reify/exact/a
(bt-bldtype probcc-variable-elimination uwe 'Linus))))
((eta 0.23100000000000001)
(memory #(12582912 4004792 1048576))
(stdout "")
(stderr ""))
test/hmm: pass
The hidden Markov model of HM, a random walk on 8 states observed at time 5,
exactly at time 10: without variable elimination (HM:147-167) the table is
upstream's bit for bit (State i is i); with it (HM:274-296), the model runs backwards and the previous state is
reified and reflected at each step, which makes the inference linear in the
number of steps; bit for bit when the reflected tables are in upstream's order, to
a relative 1e-9 with probcc-variable-elimination. The
normalized tables are computed in upstream's order, the order of the sum.
(define (test/hmm _)
(define forward (probcc-reify/exact/a (hm-query/forward 10)))
(define ve
(probcc-reify/exact/a (hm-query/ve variable-elimination/upstream 10)))
(define ve/sorted
(probcc-reify/exact/a (hm-query/ve probcc-variable-elimination 10)))
(⊦ equal?
'(((V 1) 0.10542819045498758)
((V 2) 0.098555483416225692)
((V 3) 0.086272791068737309)
((V 4) 0.070768449664987229)
((V 5) 0.054230446809986121)
((V 6) 0.038726200851251925)
((V 7) 0.026443647103752445)
((V 8) 0.019571040630000598))
forward)
(⊦ equal?
'(((V 8) 0.039142374827817968)
((V 7) 0.052887690865193897)
((V 6) 0.077452982599884354)
((V 5) 0.10846170708279078)
((V 4) 0.14153796086470105)
((V 3) 0.17254687623907095)
((V 2) 0.19711244517581822)
((V 1) 0.21085796234472273))
(probcc-normalize (upstream-order forward)))
(⊦ equal?
'(((V 1) 0.10765110792499999)
((V 2) 0.10024138302499999)
((V 3) 0.08721199276250001)
((V 4) 0.071051930187500001)
((V 5) 0.053946932537500006)
((V 6) 0.037786982462500007)
((V 7) 0.024757760350000001)
((V 8) 0.017348160749999997))
ve)
(⊦ equal?
'(((V 8) 0.034696581724362928)
((V 7) 0.04951589206919052)
((V 6) 0.07557453173398801)
((V 5) 0.10789467428505714)
((V 4) 0.1421049261619462)
((V 3) 0.17442529371470286)
((V 2) 0.2004842696820226)
((V 1) 0.21530383062872968))
(probcc-normalize (upstream-order ve)))
(⊦ (rel-approx? 1.0000000000000001e-09) ve ve/sorted))
((eta 0.628) (memory #(12582912 4009368 1048576)) (stdout "") (stderr ""))
test/hmm/counters: pass
The investigation of complexity of HM:170-279, by counting the calls of the worker at each step: the forward and the backward model grow by a factor of 3 per step; with variable elimination the counts stay flat.
(define (test/hmm/counters _)
(⊦ equal? '#(0 8 22 62 176 502 0 0 0 0 0 0) (hm-count (hm-query/forward 4)))
(⊦ equal? '#(0 8 22 62 176 0 0 0 0 0 0 0) (hm-count (hm-query/backward 4)))
(⊦ equal?
'#(0 8 22 22 22 0 0 0 0 0 0 0)
(hm-count (hm-query/ve probcc-variable-elimination 4)))
(⊦ equal?
'#(0 8 22 22 22 20 20 22 22 22 22 0)
(hm-count (hm-query/ve probcc-variable-elimination 10))))
((eta 0.002) (memory #(12582912 4008488 1048576)) (stdout "") (stderr ""))
test/blip/naive: pass
Airplanes without laziness (BN:183-252): a plane at (3,3) has probability 1/100,
exact here; a blip at (3,3) is a false alarm or the detection of the one plane of
the bounded geometric prior (BN:192: 0.02 +. 0.98 *. 0.00130434782608695684 *. 0.9); one plane cannot
make two blips, two planes make them with probability 1/5000 (either may be at
either blip; upstream 0.000199999999999999901), and two blips come from two planes of the prior with 0.0191897654584221797 · 2/10000; and the two blips of a moving
plane, with variable elimination of the previous state. Without laziness the two
planes are drawn whole, 160000 worlds, which is as far as exact inference goes
(BN:238-245).
(define (test/blip/naive _)
(define ((evidence0 check blips) ps t) (when (= t 0) (check blips ps)))
(define (unit _) (void))
(⊦ equal?
`(((V ,(void)) 1/100))
(probcc-reify/exact
(if (equal? (cadr (bn-plane-state-0 1)) '(3 3))
(void)
(probcc-impossible))))
(⊦ equal?
'(((V 0) 0.85287846481876328)
((V 1) 0.1279317697228145)
((V 2) 0.01918976545842218))
(probcc-reify/exact (probcc-geometric/bounded 2 0.84999999999999998)))
(⊦ (rel-approx? 1.0000000000000001e-09)
`(((V ,(void)) 0.02115043478260931))
(probcc-reify/exact/a
(bn-model 1 0 (evidence0 bn-blips-at-least '((3 3))) unit)))
(⊦ (rel-approx? 1.0000000000000001e-09)
`(((V ,(void)) 0.0013043478260869568))
(probcc-reify/exact/a
(bn-model 1 0 (evidence0 bn-blips-at-least/ideal '((3 3))) unit)))
(⊦ equal?
'()
(probcc-reify/exact/a
(bn-model 1 0 (evidence0 bn-blips-at-least/ideal '((3 3) (3 5))) unit)))
(⊦ equal?
`(((V ,(void)) 1/5000))
(probcc-reify/exact
(let* ((p2 (bn-plane-state-0 2)) (p1 (bn-plane-state-0 1)))
(bn-blips-at-least/ideal '((3 3) (3 5)) (list p1 p2)))))
(⊦ (rel-approx? 1.0000000000000001e-09)
'(((V 2) 3.8379530916844328e-06))
(probcc-reify/exact/a
(bn-model
2
0
(evidence0 bn-blips-at-least/ideal '((3 3) (3 5)))
length)))
(⊦ (rel-approx? 1.0000000000000001e-09)
`(((V ,(void)) 0.00059257469565211177))
(probcc-reify/exact/a (bn-model 1 1 bn-two-blips unit))))
((eta 1.7410000000000001)
(memory #(12582912 4010624 1048576))
(stdout "")
(stderr ""))
test/music: pass
Motivic development, the warm-up (M1:147-191): the motive A B C is split at a
random place and each part transformed by an operation chosen at random; given that
the result has three notes, the first A sharp, the distribution of the second note.
The notes are lazy lists of (aux hansei lazy), forced only as far as the evidence asks; the table, the evidence 0.054 and the
normalized table are upstream's bit for bit, from 29 accepted worlds, and so is
importance sampling of 1000 worlds with (probcc-random-selector 17) (M1:195-201; "1000 samples is quite good, actually").
(define (test/music _)
(define table (probcc-reify/exact/a m1-main))
(define table/upstream (upstream-order table))
(define evidence (total table/upstream))
(⊦ equal?
'(((V C) 0.027000000000000003)
((V Csharp) 0.017999999999999999)
((V B) 0.0074999999999999997)
((V Fsharp) 0.00045000000000000004)
((V F) 0.00045000000000000004)
((V E) 0.00045000000000000004)
((V G) 0.00015000000000000004))
table)
(⊦ equal? 0.053999999999999999 evidence)
(⊦ equal?
'(((V G) 0.0027777777777777788)
((V Fsharp) 0.008333333333333335)
((V F) 0.008333333333333335)
((V E) 0.008333333333333335)
((V Csharp) 0.33333333333333331)
((V C) 0.50000000000000011)
((V B) 0.1388888888888889))
(probcc-map
(λ1-match/first ((,v ,p) (list v (/ p evidence))))
table/upstream))
(⊦ equal? 29 (accepted m1-main))
(⊦ equal?
'(((V C) 0.026973000000000195)
((V Csharp) 0.017982000000000154)
((V B) 0.0078570000000000185)
((V Fsharp) 0.00060749999999999986)
((V F) 0.00060749999999999986)
((V E) 0.00060749999999999986)
((V G) 0.00020249999999999996))
(probcc-sample-importance/a (probcc-random-selector 17) 1000 m1-main)))
((eta 0.014) (memory #(12582912 4080920 1048576)) (stdout "") (stderr ""))
test/colored-balls: pass
Colored balls (CB:32-48), exactly, for three draws that all look blue: the number
of balls in the urn, 1 to 8 a priori, with the color of every ball memoized. With
the exact 4/5 the table is the analytic one, 1/8 · 1/2^n · Σ_colorings (Σ_b q(color b) / n)^3; with upstream's 0.8 it agrees to 1e-6. Upstream only samples this model; the
sampling of ten draws is planned for the long suite.
(define (cb-observed-color q c)
(if (probcc-coin q) c (case c ((Blue) 'Green) ((Green) 'Blue))))
(define ((cb-model-nballs q obs))
(let* ((nballs (probcc-uniform/range 1 9))
(ball-color
(probcc-memo (λ (b) (probcc-uniform/either '(Blue Green))))))
(let loop ((obs obs))
(unless
(null? obs)
(let1 (b (probcc-uniform/range 1 (add1 nballs)))
(unless
(eq? (cb-observed-color q (ball-color b)) (car obs))
(probcc-impossible)))
(loop (cdr obs))))
nballs))
(define (test/colored-balls _)
(⊦ equal?
'(((V 1) 13/400)
((V 2) 77/3200)
((V 3) 17/800)
((V 4) 127/6400)
((V 5) 19/1000)
((V 6) 59/3200)
((V 7) 101/5600)
((V 8) 227/12800))
(probcc-reify/exact/a (cb-model-nballs 4/5 '(Blue Blue Blue))))
(⊦= '(((V 1) 0.032500000000000001)
((V 2) 0.024062500000000001)
((V 3) 0.021250000000000002)
((V 4) 0.01984375)
((V 5) 0.019)
((V 6) 0.018437499999999999)
((V 7) 0.018035714285714287)
((V 8) 0.017734375))
(probcc-reify/exact/a
(cb-model-nballs 0.80000000000000004 '(Blue Blue Blue)))))
((eta 0.090999999999999998)
(memory #(12582912 4013864 1048576))
(stdout "")
(stderr ""))
test/incite/two-bit: pass
The two-bit game (IN:78-207): the hearer interprets the form F_b1 as two people
coming and buys two pieces of food; the speaker, who knows that two come, best says
F_b1 when the hearer infers exactly; if the hearer only draws one sample, what the
speaker says does not matter, with two samples it pays, by 1%, to tell the truth.
The utilities are -0. where nothing is lost, as upstream's; the
tables are upstream's bit for bit, in the Scheme order (the order of the form
variables F_b0, F_b1 is the same in OCaml and by name).
(define (test/incite/two-bit _)
(define (ut1 fv) (eq? fv 'F_b1))
(define (speaker inferencer)
(in-agent-S-model in-forms in-agent-A-model in-utility inferencer))
(⊦ equal?
'(((V (-0.0 (A 2))) 0.25)
((V (-10.0 (A 3))) 0.25)
((V (-10.0 (A 1))) 0.25)
((V (-20.0 (A 0))) 0.25))
(probcc-reify/exact/a (in-agent-A-model ut1)))
(⊦ equal?
'((-0.0 (A 2)))
(in-argmax (probcc-reify/exact/a (in-agent-A-model ut1))))
(⊦ equal?
'(((V (-0.0 (F_b1))) 1.0)
((V (-10.0 (F_b1 F_b0))) 1.0)
((V (-10.0 (F_b0))) 1.0)
((V (-20.0 ())) 1.0))
(probcc-reify/exact/a (speaker in-exact)))
(⊦ equal?
'(((V (-10.0 (F_b1 F_b0))) 0.5)
((V (-10.0 (F_b1))) 0.5)
((V (-10.0 (F_b0))) 0.5)
((V (-10.0 ())) 0.5)
((V (-0.0 (F_b1 F_b0))) 0.25)
((V (-0.0 (F_b1))) 0.25)
((V (-0.0 (F_b0))) 0.25)
((V (-0.0 ())) 0.25)
((V (-20.0 (F_b1 F_b0))) 0.25)
((V (-20.0 (F_b1))) 0.25)
((V (-20.0 (F_b0))) 0.25)
((V (-20.0 ())) 0.25))
(probcc-reify/exact/a (speaker (in-rejection 1))))
(⊦ equal?
'(((V (-10.0 (F_b1 F_b0))) 0.5078125)
((V (-10.0 (F_b1))) 0.5)
((V (-10.0 (F_b0))) 0.5)
((V (-10.0 ())) 0.4921875)
((V (-0.0 (F_b1))) 0.26171875)
((V (-20.0 ())) 0.26171875)
((V (-0.0 (F_b1 F_b0))) 0.25390625)
((V (-0.0 (F_b0))) 0.25)
((V (-20.0 (F_b0))) 0.25)
((V (-0.0 ())) 0.24609375)
((V (-20.0 (F_b1 F_b0))) 0.23828125)
((V (-20.0 (F_b1))) 0.23828125))
(probcc-reify/exact/a (speaker (in-rejection 2)))))
((eta 0.0050000000000000001)
(memory #(12582912 4019944 1048576))
(stdout "")
(stderr ""))
test/incite/some-all: pass
The some/all game (IN:210-426): the forms are the sixteen subsets of all, some,
not all and none, and a conjunction is understood only when the hearer also guesses
the conjunction state. With the symmetric utility the speaker best says some; with
the asymmetric one, where buying too little costs ten times more, not all, or some
and not all, when the hearer infers exactly, and some, or nothing, when the hearer
draws two samples. OCaml orders the form variables by declaration, pmap-compare by
name (deviation 6), so the tables and maxima keyed by forms are compared as sets.
IN:416-422 sample the speaker with random_selector: the records do not replay (plan D20), but what upstream gives under OCaml today
does, bit for bit.
(define (test/incite/some-all _)
(define (ut1 fv) (eq? fv 'F_some))
(define (ut2 fv) (and (memq fv '(F_some F_notall)) #t))
(define agent-A (in-agent-A-model/some-all in-utility))
(define agent-A/asymmetric (in-agent-A-model/some-all in-utility/asymmetric))
(define (speaker agent-A utility inferencer)
(in-agent-S-model in-all-forms agent-A utility inferencer))
(⊦ equal?
'(((V (-10.0 (A 2))) 0.5)
((V (-0.0 (A 3))) 0.25)
((V (-0.0 (A 2))) 0.25)
((V (-0.0 (A 1))) 0.25)
((V (-10.0 (A 3))) 0.25)
((V (-10.0 (A 1))) 0.25)
((V (-10.0 (A 0))) 0.25)
((V (-20.0 (A 3))) 0.25)
((V (-20.0 (A 1))) 0.25)
((V (-20.0 (A 0))) 0.25)
((V (-30.0 (A 0))) 0.25))
(probcc-reify/exact/a (agent-A ut1)))
(⊦ equal?
'((-6.666666666666667 (A 2)))
(in-argmax (probcc-reify/exact/a (agent-A ut1))))
(⊦ equal?
'(((V (-0.0 (A 2))) 0.125)
((V (-0.0 (A 1))) 0.125)
((V (-10.0 (A 3))) 0.125)
((V (-10.0 (A 2))) 0.125)
((V (-10.0 (A 1))) 0.125)
((V (-10.0 (A 0))) 0.125)
((V (-20.0 (A 3))) 0.125)
((V (-20.0 (A 0))) 0.125))
(probcc-reify/exact/a (agent-A ut2)))
(⊦ equal?
'((-5.0 (A 2)) (-5.0 (A 1)))
(in-argmax (probcc-reify/exact/a (agent-A ut2))))
(let1 (table (probcc-reify/exact/a (speaker agent-A in-utility in-exact)))
(⊦ equal?
(canon '(((V (-0.0 (F_some F_none))) 0.25)
((V (-0.0 (F_some F_notall F_none))) 0.25)
((V (-0.0 (F_some F_notall))) 0.5)
((V (-0.0 (F_some))) 1.0)
((V (-0.0 (F_all F_none))) 0.25)
((V (-0.0 (F_all F_notall F_none))) 0.25)
((V (-0.0 (F_all F_notall))) 0.25)
((V (-0.0 (F_all F_some F_none))) 0.25)
((V (-0.0 (F_all F_some F_notall F_none))) 0.25)
((V (-0.0 (F_all F_some F_notall))) 0.25)
((V (-0.0 ())) 0.5)
((V (-10.0 (F_notall))) 1.0)
((V (-10.0 (F_some F_none))) 0.5)
((V (-10.0 (F_some F_notall F_none))) 0.5)
((V (-10.0 (F_some F_notall))) 0.5)
((V (-10.0 (F_all F_none))) 0.5)
((V (-10.0 (F_all F_notall F_none))) 0.5)
((V (-10.0 (F_all F_notall))) 0.5)
((V (-10.0 (F_all F_some F_none))) 0.5)
((V (-10.0 (F_all F_some F_notall F_none))) 0.5)
((V (-10.0 (F_all F_some F_notall))) 0.5)
((V (-10.0 (F_all F_some))) 1.0)
((V (-10.0 (F_all))) 1.0)
((V (-10.0 ())) 0.5)
((V (-20.0 (F_none))) 1.0)
((V (-20.0 (F_notall F_none))) 1.0)
((V (-20.0 (F_some F_none))) 0.25)
((V (-20.0 (F_some F_notall F_none))) 0.25)
((V (-20.0 (F_all F_none))) 0.25)
((V (-20.0 (F_all F_notall F_none))) 0.25)
((V (-20.0 (F_all F_notall))) 0.25)
((V (-20.0 (F_all F_some F_none))) 0.25)
((V (-20.0 (F_all F_some F_notall F_none))) 0.25)
((V (-20.0 (F_all F_some F_notall))) 0.25)))
(canon table))
(⊦ equal? '((-0.0 (F_some))) (in-argmax table)))
(⊦ equal?
'((-8.90625 (F_some)))
(in-argmax
(probcc-reify/exact/a (speaker agent-A in-utility (in-rejection 2)))))
(⊦ equal?
'(((V (-0.0 (A 3))) 0.25)
((V (-0.0 (A 2))) 0.25)
((V (-0.0 (A 1))) 0.25)
((V (-10.0 (A 3))) 0.25)
((V (-10.0 (A 2))) 0.25)
((V (-20.0 (A 3))) 0.25)
((V (-100.0 (A 2))) 0.25)
((V (-100.0 (A 1))) 0.25)
((V (-100.0 (A 0))) 0.25)
((V (-200.0 (A 1))) 0.25)
((V (-200.0 (A 0))) 0.25)
((V (-300.0 (A 0))) 0.25))
(probcc-reify/exact/a (agent-A/asymmetric ut1)))
(⊦ equal?
'((-10.0 (A 3)))
(in-argmax (probcc-reify/exact/a (agent-A/asymmetric ut1))))
(⊦ equal?
'(((V (-0.0 (A 2))) 0.125)
((V (-0.0 (A 1))) 0.125)
((V (-10.0 (A 3))) 0.125)
((V (-10.0 (A 2))) 0.125)
((V (-20.0 (A 3))) 0.125)
((V (-100.0 (A 1))) 0.125)
((V (-100.0 (A 0))) 0.125)
((V (-200.0 (A 0))) 0.125))
(probcc-reify/exact/a (agent-A/asymmetric ut2)))
(⊦ equal?
'((-5.0 (A 2)))
(in-argmax (probcc-reify/exact/a (agent-A/asymmetric ut2))))
(⊦ equal?
(canon '((-0.0 (F_notall)) (-0.0 (F_some F_notall))))
(canon (in-argmax
(probcc-reify/exact/a
(speaker agent-A/asymmetric in-utility/asymmetric in-exact)))))
(⊦ equal?
(canon '((-60.078125 (F_some)) (-60.078125 ())))
(canon (in-argmax
(probcc-reify/exact/a
(speaker
agent-A/asymmetric
in-utility/asymmetric
(in-rejection 2))))))
(⊦ equal?
'((-0.0 (F_notall)))
(in-argmax
(probcc-reify/exact/a
(speaker
agent-A/asymmetric
in-utility/asymmetric
(in-importance 2)))))
(let1 (speaker
(speaker agent-A/asymmetric in-utility/asymmetric (in-rejection 10)))
(⊦ equal?
'((-25.694444444444443 (F_some)))
(in-argmax
(probcc-sample-importance/a
(probcc-random-selector 17)
100
speaker)))
(⊦ equal?
'((-19.003472222222221 (F_some)))
(in-argmax
(probcc-sample-importance/a
(probcc-random-selector 17)
1000
speaker)))))
((eta 0.40600000000000003)
(memory #(12582912 4031208 1048576))
(stdout "")
(stderr ""))
test/blip/ailog: pass
The AILog model of the blips (BA:181-216), where every prior and every step of the dynamics is memoized, so that a fact observed once stays the same: the position of a plane asked twice, at one and at two times, and the blips at one and two times. The memo tables live in the world-local memory, so the model is instantiated inside the reification, as upstream's functor is.
(define (test/blip/ailog _)
(define (query q)
(probcc-reify/exact
(match/first (ba-blip-model) ((,pos ,anyblip1 ,_) (q pos anyblip1)))))
(⊦ equal?
`(((V ,(void)) 1/100))
(query (λ (pos _)
(if (and (equal? (pos '(1 0)) '(3 3))
(equal? (pos '(1 0)) '(3 3)))
(void)
(probcc-impossible)))))
(⊦ equal?
'()
(query (λ (pos _)
(if (and (equal? (pos '(1 0)) '(3 3))
(equal? (pos '(1 0)) '(3 4)))
(void)
(probcc-impossible)))))
(⊦ (rel-approx? 1.0000000000000001e-09)
`(((V ,(void)) 0.0015))
(query (λ (pos _)
(if (and (equal? (pos '(1 0)) '(3 3))
(equal? (pos '(1 1)) '(3 4)))
(void)
(probcc-impossible)))))
(⊦ (rel-approx? 1.0000000000000001e-09)
'(((V #f) 0.97118000000000115) ((V #t) 0.028820000000000002))
(query (λ (_ anyblip1) (anyblip1 '((3 3) 0)))))
(⊦ (rel-approx? 1.0000000000000001e-09)
'(((V #f) 0.97118000000000115) ((V #t) 0.028820000000000002))
(query (λ (_ anyblip1)
(and (anyblip1 '((3 3) 0)) (anyblip1 '((3 3) 0))))))
(⊦ (rel-approx? 1.0000000000000001e-09)
'(((V #f) 0.99595959400000089) ((V #t) 0.0019196860000000003))
(query (λ (_ anyblip1)
(and (anyblip1 '((3 3) 0)) (anyblip1 '((3 4) 1)))))))
((eta 0.02) (memory #(12582912 4020648 1048576)) (stdout "") (stderr ""))
test/blip/lazy: pass
Lazy planes and the whole screen observed (BL:263-428): a coordinate is drawn
only when an observation compares it, so two and three planes are within reach of
exact inference: two blips may come from two planes, or from three when two are
in one cell, and three blips from three planes with probability 6 · 0.01^3 · 0.00287020304028914744. The motion of one plane over
two steps is inferred in one query, and step by step, by reflecting the exact
state of the first step into the second (sequential Monte-Carlo, BL:419-428), with
the same result; so are three moving planes over three steps (BL:430-493), whose
number of worlds stays 384, since one motion only is consistent with the blips.
(define (test/blip/lazy _)
(define ((evidence0 blips) ps t) (when (= t 0) (bl-blips-ideal blips ps)))
(define (count blips n)
(probcc-reify/exact
(car (bl-step (bl-parm-n n) bl-state0 0 (evidence0 blips)))))
(⊦ equal?
`(((V ,(void)) 1/100))
(probcc-reify/exact
(match/first
(bl-lplane-state-0 1)
((,_ (,x ,y) ,_)
(if (and (= (x) 3) (= (y) 3)) (void) (probcc-impossible))))))
(⊦ (rel-approx? 1.0000000000000001e-09)
`(((V ,(void)) 0.0013043478260869568))
(probcc-reify/exact
(begin
(bl-step (bl-parm-n 1) bl-state0 0 (evidence0 '((3 3))))
(void))))
(⊦ equal?
`(((V ,(void)) 1/5000))
(probcc-reify/exact
(let* ((p2 (probcc-letlazy (τ (bl-lplane-state-0 2))))
(p1 (probcc-letlazy (τ (bl-lplane-state-0 1)))))
(bl-blips-ideal '((3 3) (3 5)) (list p1 p2)))))
(⊦ equal?
'()
(probcc-reify/exact
(begin
(bl-step (bl-parm-n 1) bl-state0 0 (evidence0 '((3 3) (3 5))))
(void))))
(⊦ (rel-approx? 1.0000000000000001e-09)
'(((V 2) 3.8379530916844328e-06))
(count '((3 3) (3 5)) 2))
(⊦ (rel-approx? 1.0000000000000001e-09)
'(((V 2) 3.8269373870522e-06) ((V 3) 1.7221218241734858e-08))
(count '((3 3) (3 5)) 3))
(⊦ (rel-approx? 1.0000000000000001e-09)
'(((V 3) 1.7221218241734858e-08))
(count '((3 3) (3 5) (3 7)) 3))
(⊦ (rel-approx? 1.0000000000000001e-09)
`(((V ,(void)) 0.00019565217391304352))
(probcc-reify/exact
(let1 (st (bl-step (bl-parm-n 1) bl-state0 0 bl-two-blips-ideal))
(bl-step bl-parm-no (bl-lazy-state st) 1 bl-two-blips-ideal)
(void))))
(let1 (q2b-s0
(probcc-reify/exact
(bl-step (bl-parm-n 1) bl-state0 0 bl-two-blips-ideal)))
(⊦ (rel-approx? 1.0000000000000001e-09)
`(((V ,(void)) 0.00019565217391304352))
(probcc-reify/exact
(let1 (st (bl-lazy-state (probcc-reflect q2b-s0)))
(bl-step bl-parm-no st 1 bl-two-blips-ideal)
(void)))))
(let* ((q3b-s0
(probcc-reify/exact
(bl-step (bl-parm-n 3) bl-state0 0 blip333-ideal)))
(q3b-s1
(probcc-reify/exact
(bl-step
bl-parm-no
(bl-lazy-state (probcc-reflect q3b-s0))
1
blip333-ideal)))
(q3b-s2
(probcc-reify/exact
(bl-step
bl-parm-no
(bl-lazy-state (probcc-reflect q3b-s1))
2
blip333-ideal)))
(nplanes (λ (st) (length (cadr st)))))
(⊦ equal? '(384 384 384) (probcc-map length (list q3b-s0 q3b-s1 q3b-s2)))
(⊦ (rel-approx? 1.0000000000000001e-09)
'(((V 3) 1.7221218241734858e-08))
(probcc-map-dist nplanes q3b-s0))
(⊦ (rel-approx? 1.0000000000000001e-09)
'(((V 3) 5.8121611565855201e-11))
(probcc-map-dist nplanes q3b-s1))
(⊦ (rel-approx? 1.0000000000000001e-09)
'(((V 3) 2.0760830241309679e-12))
(probcc-map-dist nplanes q3b-s2))))
((eta 0.26600000000000001)
(memory #(12582912 4533888 1048576))
(stdout "")
(stderr ""))