(* Probabilistic embedded domain specific language *) (* Tests of exact inference and the comparison of the syntax with that of IBAL*) open Ptypes;; open ProbM;; (* ------------------------------------------------------------------------ * Basic tests *) (* illustration of sharing/nonsharing: flip twice and `and' the results *) let tflip2_shared () = (* sharing of flips *) let v = flip 0.5 in v && v;; let [(0.5, V true); (0.5, V false)] = exact_reify tflip2_shared;; (* reify: done; 2 accepted 0 rejected 0 left *) let tflip2_noshared () = (* independent flips *) flip 0.5 && flip 0.5;; let [(0.25, V true); (0.75, V false)] = exact_reify tflip2_noshared;; (* Result: reify: done; 3 accepted 0 rejected 0 left *) let talarm () = let earthquake = flip 0.01 and burglary = flip 0.1 in if earthquake then if burglary then flip 0.99 else flip 0.2 else if burglary then flip 0.98 else flip 0.01 ;; let () = assert ( exact_reify talarm = [(0.108720000000000011, V true); (0.891280000000000072, V false)]);; (* reify: done; 8 accepted 0 rejected 0 left *) (* The famous sprinkle example: given that the grass is wet on a given day, did it rain or did the sprinkler come on? This is the basic model for one single day, just for the sake of introducing our library. More general models are considered in slazy.ml *) let grass_model () = let rain = flip 0.3 and sprinkler = flip 0.5 in let grass_is_wet = (flip 0.9 && rain) || (flip 0.8 && sprinkler) || flip 0.1 in if not grass_is_wet then fail (); rain;; let [(0.2838, V true); (0.321999999999999897, V false)] = exact_reify grass_model;; (* 14 worlds are examined *) (* Tests of library distributions *) let [(0.125, V 8); (0.125, V 7); (0.125, V 6); (0.125, V 5); (0.125, V 4); (0.125, V 3); (0.125, V 2); (0.125, V 1)] = exact_reify (fun () -> uniform_range 1 8);; (* Tests of the geometric distribution: exact inference to a specific depth. *) let [(0.00286875000000000107, V 3); (0.0191250000000000031, V 2); (0.1275, V 1); (0.85, V 0); (1.13906250000000092e-05, C _); (6.45468750000000426e-05, C _); (0.000430312500000000248, C _)] = Inference.explore (Some 4) (reify0 (fun () -> geometric 0.85));; let [(1.4523050597099037e-06, V 7); (9.68203373139935684e-06, V 6); (6.45468915426623699e-05, V 5); (0.000430312610284415691, V 4); (0.00286875073522943743, V 3); (0.0191250049015295812, V 2); (0.127500032676863856, V 1); (0.850000217845759, V 0)] = exact_reify (fun () -> geometric_bounded 7 0.85);; (* ------------------------------------------------------------------------ * Tests from IBAL, to compare the syntax of IBAL with our EDSL *) (* foo.ibl, program 1 // P(e) = 0.6; P(TT) = 1.0 obs { y : _, z : true } in let x = dist [ 0.6 : true, 0.4 : false ] in let y = x in let z = x in ({ y = y, z = z }) *) let tfooibl1 () = observe (fun (_,z) -> z = true) (fun () -> let x = flip 0.6 in let y = x in let z = x in (y,z));; let () = assert ( (exact_reify tfooibl1) = [(0.6, V (true,true))]);; (* reify: done; 1 accepted 1 rejected 0 left val rtfooibl1 : unit pV = PV [(0.6, V ())] *) (* foo.ibl, program 2 // P(e) = 1; P(d,e) = 0.2; P(c,e) = 0.48; P(b,e) = 0.32 { fst = let x = 'a in let f(x,y) = if dist [ 0.2 : true, 0.8 : false ] then x else y in f('d, dist [ 0.4 : 'b, 0.6 : 'c ]), snd = 'e } *) (* Note how similar the syntaxes are. The results match, too *) let tfooibl2 () = (let x = 'a' in let f (x,y) = if dist [ (0.2,true); (0.8,false)] then x else y in f('d', dist [ (0.4,'b'); (0.6,'c') ]), 'e');; let () = assert ( (exact_reify tfooibl2) = [(0.2, V ('d', 'e')); (0.48, V ('c', 'e')); (0.320000000000000062, V ('b', 'e'))]);; (* reify: done; 4 accepted 0 rejected 0 left *) (* foo.ibl, program 3 //P(e) = 0.25; P(T) = 1.0 obs true in ({ f = if dist [ 0.1 : true, 0.9 : false ] then dist [ 0.7 : 'a, 0.3 : 'b ] else dist [ 0.2 : 'a, 0.8 : 'b ], g = 'c }).f == 'a *) let tfooibl3 () = observe (fun x -> x = true) (fun () -> (object method mf = if dist [ (0.1, true); (0.9, false) ] then dist [ (0.7, 'a'); (0.3, 'b') ] else dist [ (0.2, 'a'); (0.8, 'b') ] method mg = 'c' end)#mf = 'a');; let () = assert ( (exact_reify tfooibl3) = [(0.25, V true)]);; (* reify: done; 2 accepted 2 rejected 0 left *) (* test2.ibl let x = obs { z : 'a } in dist [ 0.01 : { z = 'a, w = 'b }, 0.02 : { z = 'a, w = 'c }, 0.97 : { z = 'd, w = 'e } ] in let y = if x.w == 'b then dist [ 0.9 : true, 0.1 : false ] else if x.w == 'c then dist [ 0.6 : true, 0.4 : false ] else dist [ 0.2 : true, 0.8 : false ] in y *) let tiblt2 () = let x = observe (fun r -> r#mz = 'a') (fun () -> dist [ (0.01, object method mz = 'a' method mw = 'b' end); (0.02, object method mz = 'a' method mw = 'c' end); (0.97, object method mz = 'd' method mw = 'e' end) ]) in let y = if x#mw = 'b' then dist [ (0.9, true); (0.1, false) ] else if x#mw = 'c' then dist [ (0.6, true); (0.4, false) ] else dist [ (0.2, true); (0.8, false) ] in y;; let () = assert ( (exact_reify tiblt2) = [(0.021, V true); (0.00900000000000000105, V false)]);; (* reify: done; 4 accepted 1 rejected 0 left *) (* test from the Music TR, Chap 4 let f() = dist [ 0.1 : true, 0.9 : false ] in let g() = dist [ 0.3 : true, 0.7 : false ] in observe true in dist [ 0.8 : f(), 0.2 : g() ] *) let tmusictr41 () = let f () = flip 0.1 in let g () = flip 0.3 in observe (fun x -> x = true) (fun () -> dist [ (0.8, f); (0.2, g) ] ());; let () = assert ( (exact_reify tmusictr41) = [(0.14, V true)]);; (* reify: done; 2 accepted 2 rejected 0 left *) (* test from the Music TR, Chap 4 let f() = dist [ 0.01 : true, 0.99 : false ] in (observe { p : 'a } in if f() then { p = 'a, q = dist [ 0.3 : true, 0.7 : false ] } else { p = 'b, q = true }).q *) let tmusictr42 () = let f () = flip 0.01 in (observe (fun r -> r#mp = 'a') (fun () -> if f () then object method mp = 'a' method mq = flip 0.3 end else object method mp = 'b' method mq = true end))#mq;; let () = assert ( (exact_reify tmusictr42) = [(0.003, V true); (0.00699999999999999928, V false)]);; (* reify: done; 2 accepted 1 rejected 0 left *) (* ------------------------------------------------------------------------ * Reflection and reification: how to express early world collapse and optimize computations *) (* The exact_reify is not too naive, which is can be seen from the result of tflip2_noshared already. The following test makes the optimality (sophistication) clearer. We evaluate AND of ten flips *) let flips_and n () = let rec loop n = if n = 1 then flip 0.5 else flip 0.5 && loop (n-1) in loop n;; let rtflips_10and = exact_reify (flips_and 10);; (* reify: done; 11 accepted 0 rejected 0 left val rtflips_10and : bool pV = PV [(0.0009765625, V true); (0.9990234375, V false)] *) (* 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. That is, when flip 0.5 yields false, we no longer need to evaluate loop (n-1). We finish right away. It is imperative that we write [flip 0.5 && loop (n-1)] rather than [loop (n-1) && flip 0.5]. Such an optimality (exponential reduction) came from the fact that we started to collapse worlds early. Ostensibly, the algorithm collapses worlds only when we finished evaluating it -- when it is too late. Ideally, we should see if two worlds in progress are the same -- but that requires us to (extensionally!) compare continuations -- which seems impossible. In the case above, the comparsion of the continuations and collapse happened implicitly, due to the special _control_ structure of the algorithm: (false && k) yields false regarding of the continuation. Thus evaluation of that _control structure_ embodies the collapse of all worlds represented by k. But the same trick won't work if we want to compute XOR of ten flips. XOR really needs the values of its two arguments in order to proceed. One, neatest, solution to the early collapse is suggested by Ken. He essentially showed the way to _intensionally_ compare two continuations. Indeed, in our case, the continuations result from splitting the current world. One continuation receives 'true' and another receives 'false'. So, the two worlds differ *only* in the value of the variable that receives the result of the flip. If that variable is garbage collected, then the difference between worlds disappears and we can collapse the worlds in progress. Again, the _intension_ of the world (continuation) is the value of its flips and derived quantities. If we can combine the value of the flips and forget about the original flips, two worlds become intensionally equal and so can be collapsed into one. It seems however that solution, in some circumstances, may have a dual. Dual in a sense like CBV is dual to CBN. Variable is a continuation, and garbage collection is done by `region management'. The solution is again _non_ tail recursive. The solution critically relies on our ability to reify and reflect the probilistic monadic computations. Here we truly need delimited continuations. *) let xor (x:bool) y = x <> y;; (* First, the naive evaluation of the xor of 10 flips *) let flips_xor n () = let rec loop n = if n = 1 then flip 0.5 else xor (flip 0.5) (loop (n-1)) in loop n;; let rtflips_10xor = exact_reify (flips_xor 10);; (* reify: done; 1024 accepted 0 rejected 0 left val rtflips_10xor : bool pV = PV [(0.5, V true); (0.5, V false)] *) (* The above was naive: we indeed had to evaluate 512+512 possible worlds. Here is an optimal xor *) let flips_xor' n () = let rec loop n = if n = 1 then flip 0.5 else let r = reflect (exact_reify (fun () -> loop (n-1))) in xor (flip 0.5) r in loop n;; let [(0.5, V true); (0.5, V false)] = exact_reify (flips_xor' 10);; (* reify: done; 2 accepted 0 rejected 0 left reify: done; 4 accepted 0 rejected 0 left reify: done; 4 accepted 0 rejected 0 left reify: done; 4 accepted 0 rejected 0 left reify: done; 4 accepted 0 rejected 0 left reify: done; 4 accepted 0 rejected 0 left reify: done; 4 accepted 0 rejected 0 left reify: done; 4 accepted 0 rejected 0 left reify: done; 4 accepted 0 rejected 0 left reify: done; 4 accepted 0 rejected 0 left *) (* That is obviously better: we needed only 38 threads to do the job. The 'let' operator is important, for `sharing': we needed to compute the join in the original thread, so we can share the result among the childrens, after the flip. *)