diff --git a/ChangeLog.md b/ChangeLog.md new file mode 100644 --- /dev/null +++ b/ChangeLog.md diff --git a/LICENSE b/LICENSE new file mode 100644 --- /dev/null +++ b/LICENSE @@ -0,0 +1,30 @@ +Copyright (c) 2021, Niki Vazou + +All rights reserved. + +Redistribution and use in source and binary forms, with or without +modification, are permitted provided that the following conditions are met: + + * Redistributions of source code must retain the above copyright + notice, this list of conditions and the following disclaimer. + + * Redistributions in binary form must reproduce the above + copyright notice, this list of conditions and the following + disclaimer in the documentation and/or other materials provided + with the distribution. + + * Neither the name of Niki Vazou nor the names of other + contributors may be used to endorse or promote products derived + from this software without specific prior written permission. + +THIS SOFTWARE IS PROVIDED BY THE COPYRIGHT HOLDERS AND CONTRIBUTORS +"AS IS" AND ANY EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT +LIMITED TO, THE IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR +A PARTICULAR PURPOSE ARE DISCLAIMED. IN NO EVENT SHALL THE COPYRIGHT +OWNER OR CONTRIBUTORS BE LIABLE FOR ANY DIRECT, INDIRECT, INCIDENTAL, +SPECIAL, EXEMPLARY, OR CONSEQUENTIAL DAMAGES (INCLUDING, BUT NOT +LIMITED TO, PROCUREMENT OF SUBSTITUTE GOODS OR SERVICES; LOSS OF USE, +DATA, OR PROFITS; OR BUSINESS INTERRUPTION) HOWEVER CAUSED AND ON ANY +THEORY OF LIABILITY, WHETHER IN CONTRACT, STRICT LIABILITY, OR TORT +(INCLUDING NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY OUT OF THE USE +OF THIS SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF SUCH DAMAGE. diff --git a/README.md b/README.md new file mode 100644 --- /dev/null +++ b/README.md @@ -0,0 +1,86 @@ +# safe-coupling +Library for relational verification of probabilistic algorithms. + +Supports two proving methods: + - Upper bound _Kantorovich distance_ between two distributions + - Establish a _boolean relation_ on samples from two distributions (this is stronger) + +Includes two larger examples of verification: + - Stability of stochastic gradient descent (src/SGD) using Kantorovich distance + - Convergence of temporal difference learning (src/TD0) using boolean relations + +## A smaller example (src/Bins/Bins.hs) + +This function recursively counts how many times the ball hit the bin after n attempted throws: + + bins :: Double -> Nat -> PrM Nat + bins _ 0 = ppure 0 + bins p n = liftA2 (+) (bernoulli p) (bins p (n - 1)) + +Throws succeed with probability _p_ which is simulated by `bernoulli p`. The function returns a distribution over natural numbers. When comparing results of two throwers with respective chances of success _p_ and _q > p_, we expect the second thrower to score notably better with the increase of _n_. Formally, we can show that Kantorovich distance between `bins p n` and `bins q n` is upper bounded by _(q - p)·n_. + +## Proof (src/Bins/Theorem.hs) + +The proof uses four definitions from the library: + * In the first case, no throws were made. Axiom `pureDist` allows deriving Kantorovich distance between pure expressions. In our case, _0_ and _0_. + * In the second case, axiom `liftA2Dist` derives Kantorovich distance between the inductive cases. Numeric arguments specify the expected bound in format _a·x + b·y + c_ where _x_ and _y_ are bounds for the second and third arguments of `liftA2` respectively. As the last argument, the axiom requires proof of linearity of plus. It is empty since it can be automatically constructed by an SMT-solver. + * Axiom `bernoulliDist` upper bounds the distance between calls to `bernoulli` with _q - p_ — this is our _x_. The second upper bound _y_ is provided by a recursive call to our theorem. + * A function `distInt` is used to measure the distance between arguments of `liftA2`. In this case, all of them provide integer values. A pre-defined distance between _n_ and _m_ is _|n - m|_ but this allows customization. + +``` +{-@ binsDist :: p:Prob -> {q:Prob|p <= q} -> n:Nat + -> {dist (kant distInt) (bins p n) (bins q n) <= n * (q - p)} / [n] @-} +binsDist :: Prob -> Prob -> Nat -> () +binsDist p q 0 = pureDist distInt 0 0 +binsDist p q n += liftA2Dist d d d 1 (q - p) 1 ((n - 1) * (q - p)) 0 + (+) (bernoulli p) (bins p (n - 1)) + (+) (bernoulli q) (bins q (n - 1)) + (bernoulliDist d p q) + (binsDist p q (n - 1)) + (\_ _ _ _ -> ()) +where + d = distInt +``` + +This concludes the mechanized proof of the boundary _(q-p)·n_. + +## Installation +1. Install stack https://docs.haskellstack.org/en/stable/install_and_upgrade/ + +2. Compile the library and case studies + + $ cd safe-coupling + $ stack install --fast + ... + Registering library for safe-coupling-0.1.0.0.. + + +3. Run unit tests on executable case studies + + $ stack test + ... + test/Spec.hs + Spec + Bins + mockbins 1 it: OK + mockbins 2 it: OK + bins 1 it: OK + bins 2 it: OK (0.02s) + exp dist mockbins: OK (0.12s) + SGD + sgd: OK + TD0 + td0 base: OK + td0 simple: OK + + All 8 tests passed (0.12s) + + safe-coupling> Test suite safe-coupling-test passed + Completed 2 action(s). + + +In case of errors try + + $ stack clean + diff --git a/Setup.hs b/Setup.hs new file mode 100644 --- /dev/null +++ b/Setup.hs @@ -0,0 +1,2 @@ +import Distribution.Simple +main = defaultMain diff --git a/safe-coupling.cabal b/safe-coupling.cabal new file mode 100644 --- /dev/null +++ b/safe-coupling.cabal @@ -0,0 +1,91 @@ +name: safe-coupling +category: Formal Methods +version: 0.1.0.0 +synopsis: Relational proof system for probabilistic algorithms +description: Relational proof system for probabilistic algorithms. Supports two proving methods: upper bound Kantorovich distance between two distributions and establish a boolean relation on samples from two distributions (the latter is stronger). +license: BSD3 +license-file: LICENSE +author: Lisa Vasilenko, Niki Vazou +maintainer: Lisa Vasilenko +homepage: https://github.com/nikivazou/safe-coupling +bug-reports: https://github.com/nikivazou/safe-coupling/issues +copyright: 2020-21 Lisa Vasilenko & Niki Vazou, IMDEA Software Institute +build-type: Simple +extra-source-files: ChangeLog.md, README.md +cabal-version: >=1.10 + +source-repository head + type: git + location: https://github.com/nikivazou/safe-coupling + +library + exposed-modules: + Data.Dist + , Data.List + , Data.Derivative + + , Monad.PrM + , Monad.PrM.Predicates + , Monad.PrM.Laws + + , Monad.PrM.Relational.TCB.Spec + , Monad.PrM.Relational.TCB.EDist + , Monad.PrM.Relational.Theorems + + , Misc.ProofCombinators + +-- Toy Examples + , Examples.ExpDist + +-- Bins Example + , Bins.Bins + , Bins.Theorem + +-- Bins Example using Applicatives + , ApplicativeBins.Bins + , ApplicativeBins.Theorem + +-- TD Case Study + , TD.TD0 + , TD.Lemmata.Relational.Update + , TD.Lemmata.Relational.Sample + , TD.Lemmata.Relational.Act + , TD.Lemmata.Relational.Iterate + , TD.Theorem + +-- SGD Case Study + , SGD.SGD + , SGD.Theorem + + build-depends: + liquid-base >= 4.14.0 && < 4.16 + , liquidhaskell >= 0.8.10 && < 0.9 + , liquid-containers >= 0.6.2 && < 0.7 + , liquid-prelude >= 0.8.10 && < 0.9 + , probability >= 0.2.7 && < 0.3 + hs-source-dirs: src + default-language: Haskell2010 + ghc-options: -fplugin=LiquidHaskell + +test-suite safe-coupling-test + type: exitcode-stdio-1.0 + main-is: Spec.hs + other-modules: + Spec.TD0 + Spec.SGD + Spec.Bins + Spec.Utils + hs-source-dirs: + test + ghc-options: -threaded -rtsopts -with-rtsopts=-N + build-tool-depends: + tasty-discover:tasty-discover + build-depends: + liquid-base + , safe-coupling + , HUnit >= 1.6.1 && < 1.7 + , tasty >= 1.2.3 && < 1.5 + , tasty-hunit >= 0.10.0 && < 1.11 + , probability >= 0.2.7 && < 0.3 + , sort >= 1.0.0 && < 1.1 + default-language: Haskell2010 diff --git a/src/ApplicativeBins/Bins.hs b/src/ApplicativeBins/Bins.hs new file mode 100644 --- /dev/null +++ b/src/ApplicativeBins/Bins.hs @@ -0,0 +1,19 @@ +{-@ LIQUID "--reflection" @-} + +module ApplicativeBins.Bins where + +import Monad.PrM +import Data.Dist + +{-@ type PDouble = {v:Double | 0 <= v } @-} + +{-@ reflect bins @-} +{-@ bins :: p:Prob -> n:PDouble -> PrM {d:Double | 0 <= d && d <= n } / [n] @-} +bins :: Double -> Double -> PrM Double +bins _ n | n < 1.0 = ppure 0 +bins p n = liftA2 plus (bins p (n - 1)) (bernoulli p) + +{-@ reflect plus @-} +{-@ plus :: x:Double -> y:Double -> {d:Double | d = x + y } @-} +plus :: Double -> Double -> Double +plus x y = x + y diff --git a/src/ApplicativeBins/Theorem.hs b/src/ApplicativeBins/Theorem.hs new file mode 100644 --- /dev/null +++ b/src/ApplicativeBins/Theorem.hs @@ -0,0 +1,35 @@ +{-@ LIQUID "--reflection" @-} +{-@ LIQUID "--ple" @-} + +module ApplicativeBins.Theorem where + +import Monad.PrM +import Data.Dist + +import Monad.PrM.Relational.TCB.EDist +import ApplicativeBins.Bins + +import Language.Haskell.Liquid.ProofCombinators +import Misc.ProofCombinators + +{-@ binsDist :: p:Prob -> {q:Prob|p <= q} -> n:PDouble + -> {dist (kant distDouble) (bins p n) (bins q n) <= n * (q - p)} / [n] @-} +binsDist :: Prob -> Prob -> Double -> () +binsDist p q n | n < 1.0 + = pureDist distDouble 0 0 + ? assert (0 <= n) + ? assert (0 <= (q - p)) + ? assert (dist (kant distDouble) (bins p n) (bins q n) <= n * (q - p)) +binsDist p q n + = liftA2Dist d d d 1 ((n - 1) * (q - p)) 1 (q - p) 0 + plus (bins p (n - 1)) (bernoulli p) + plus (bins q (n - 1)) (bernoulli q) + (binsDist p q (n - 1)) + (bernoulliDist d p q) + plusDist + where d = distDouble + +{-@ plusDist :: x1:Double -> y1:Double -> x2:Double -> y2:Double + -> {distD (plus x1 y1) (plus x2 y2) <= distD x1 x2 + distD y1 y2} @-} +plusDist :: Double -> Double -> Double -> Double -> () +plusDist _ _ _ _ = () diff --git a/src/Bins/Bins.hs b/src/Bins/Bins.hs new file mode 100644 --- /dev/null +++ b/src/Bins/Bins.hs @@ -0,0 +1,41 @@ +{-@ LIQUID "--reflection" @-} + +module Bins.Bins where + +import Monad.PrM +import Data.Dist +import Data.List + +import Prelude hiding ( map + , max + , repeat + , foldr + , fmap + , mapM + , iterate + , uncurry + ) + +{-@ type NBool = {v:Int | 0 <= v && v <= 1} @-} +type NBool = Int +{-@ type NDouble = {v:Double | 0 <= v && v <= 1} @-} +type NDouble = Double +{-@ type PDouble = {v:Double | 0 <= v } @-} + +{-@ reflect bins @-} +{-@ bins :: p:Prob -> n:PDouble -> PrM {d:Double | 0 <= d && d <= n } / [n] @-} +bins :: Double -> Double -> PrM Double +bins _ n | n < 1.0 = ppure 0 +bins p n = bind (bins p (n - 1)) (addBernoulli p (n - 1)) + +-- bins = liftA2 (+) (bins p (n - 1)) (bernoulli p) + +{-@ reflect addBernoulli @-} +{-@ addBernoulli :: Prob -> n:PDouble -> {d:Double | 0 <= d && d <= n } -> PrM {d:Double | 0 <= d && d <= n + 1 } @-} +addBernoulli :: Double -> Double -> Double -> PrM Double +addBernoulli p n x = bind (bernoulli p) (ppure . plus x) + +{-@ reflect plus @-} +{-@ plus :: x:Double -> y:Double -> {d:Double | d = x + y } @-} +plus :: Double -> Double -> Double +plus x y = x + y diff --git a/src/Bins/Theorem.hs b/src/Bins/Theorem.hs new file mode 100644 --- /dev/null +++ b/src/Bins/Theorem.hs @@ -0,0 +1,207 @@ +{-@ LIQUID "--reflection" @-} +{-@ LIQUID "--ple" @-} + +module Bins.Theorem where + +import Monad.PrM +import Monad.PrM.Laws +import Data.Dist +import Data.List + +import Prelude hiding ( flip ) + +import Monad.PrM.Predicates +import Monad.PrM.Relational.TCB.Spec +import Monad.PrM.Relational.TCB.EDist +import Monad.PrM.Relational.Theorems +import Bins.Bins + +import Language.Haskell.Liquid.ProofCombinators +import Misc.ProofCombinators + +{-@ relationalinccond :: x1:Double -> {x2:Double|x1 <= x2} -> y1:Double -> {y2:Double|leDoubleP y1 y2} + -> {lift leDoubleP ((ppure . (plus x1)) (y1)) ((ppure . (plus x2)) (y2))} @-} +relationalinccond :: Double -> Double -> Double -> Double -> () +relationalinccond x1 x2 y1 y2 = pureSpec leDoubleP + (plus y1 x1) + (plus y2 x2) + () + +{-@ relationalbinsrec :: p:Prob -> {q:Prob|leDoubleP p q} -> n:Double -> x1:Double -> {x2:Double| x1 <= x2} + -> {lift leDoubleP (addBernoulli p n x1) (addBernoulli q n x2)} / [n, 1] @-} +relationalbinsrec :: Double -> Double -> Double -> Double -> Double -> () +relationalbinsrec p q n x1 x2 + = bindSpec leDoubleP leDoubleP + (bernoulli p) (ppure . plus x1) + (bernoulli q) (ppure . plus x2) + (bernoulliSpec p q) + (relationalinccond x1 x2) + +{-@ binsSpec :: p:Prob -> {q:Prob|leDoubleP p q} -> n:Double + -> {lift leDoubleP (bins p n) (bins q n)} / [n, 0] @-} +binsSpec :: Double -> Double -> Double -> () +binsSpec p q n | n < 1 + = pureSpec leDoubleP 0 0 () +binsSpec p q n + = bindSpec leDoubleP leDoubleP + (bins p (n - 1)) (addBernoulli p (n - 1)) + (bins q (n - 1)) (addBernoulli q (n - 1)) + (binsSpec p q (n - 1)) + (relationalbinsrec p q (n - 1)) + +{-@ plusDist :: y:Double -> x1:Double -> x2:Double + -> {distD (plus y x1) (plus y x2) = distD x1 x2} @-} +plusDist :: Double -> Double -> Double -> () +plusDist _ _ _ = () + +{-@ addBernoulliDist :: p:Prob -> {q:Prob|p <= q} -> n:PDouble -> {y:PDouble|y <= n} + -> {dist (kant distDouble) (addBernoulli p n y) (addBernoulli q n y) <= q - p} @-} +addBernoulliDist :: Prob -> Prob -> Double -> Double -> () +addBernoulliDist p q n y + = dist (kant distDouble) (addBernoulli p n y) (addBernoulli q n y) + ? fmapDist distDouble distDouble + 0 + (plus y) (bernoulli p) + (plus y) (bernoulli q) + (plusDist y) + =<= dist (kant distDouble) (bernoulli p) (bernoulli q) + ? (bernoulliDist distDouble p q) + ? assert (distD 1.0 0.0 == 1.0) + ? assert (distD 1 0 * (q-p) == q-p) + =<= distD 1 0 * (q -p) + =<= q - p + *** QED + +{-@ binsDistL :: p:Prob -> {q:Prob|p <= q} -> {n:PDouble|1 <= n} + -> {dist (kant distDouble) (bins p n) (bins' p q n) <= q - p} @-} +binsDistL :: Prob -> Prob -> Double -> () +binsDistL p q n + = bindDistEq distDouble + (q - p) + (addBernoulli p (n - 1)) (bins p (n - 1)) + (addBernoulli q (n - 1)) (bins p (n - 1)) + (addBernoulliDist p q (n - 1)) + +{-@ addBinsDist :: p:Prob -> {q:Prob|p <= q} -> n:PDouble -> x:Double + -> {dist (kant distDouble) + (seqBind (bins p n) (flip (pure2 plus)) x) + (seqBind (bins q n) (flip (pure2 plus)) x) + <= n * (q - p)} / [n, 2] @-} +addBinsDist :: Prob -> Prob -> Double -> Double -> () +addBinsDist p q n x + = dist (kant distDouble) + (seqBind (bins p n) (flip (pure2 plus)) x) + (seqBind (bins q n) (flip (pure2 plus)) x) + === dist (kant distDouble) + (bind (bins p n) (flip (pure2 plus) x)) + (bind (bins q n) (flip (pure2 plus) x)) + ? flipPlus x + === dist (kant distDouble) + (bind (bins p n) (ppure . (plus x))) + (bind (bins q n) (ppure . (plus x))) + ? fmapDist distDouble distDouble + 0 + (plus x) (bins p n) + (plus x) (bins q n) + (plusDist x) + =<= dist (kant distDouble) (bins p n) (bins q n) + ? binsDist p q n + =<= n * (q - p) + *** QED + +{-@ reflect pure2 @-} +pure2 :: (a -> b -> c) -> a -> b -> PrM c +pure2 f a b = ppure (f a b) + +{-@ addBernoulliEq :: n:{Double | 0 <= n - 1 } -> p:Prob -> q:Prob + -> {addBernoulli q (n - 1) == seqBind (bernoulli q) (pure2 plus)} @-} +addBernoulliEq :: Double -> Double -> Double -> () +addBernoulliEq n p q + = extDouble (addBernoulli q (n - 1)) (seqBind (bernoulli q) (pure2 plus)) + (addBernoulliEq' n p q) + +{-@ addBernoulliEq' :: n:{Double | 0 <= n - 1 } -> p:Prob -> q:Prob + -> x:{Double | 0 <= x && x <= n - 1 } + -> {addBernoulli q (n - 1) x == seqBind (bernoulli q) (pure2 plus) x} @-} +addBernoulliEq' :: Double -> Double -> Double -> Double -> () +addBernoulliEq' n p q x + = addBernoulli q (n - 1) x + === bind (bernoulli q) (ppure . plus x) + ? extDouble (ppure . plus x) (pure2 plus x) ( + \z -> (ppure . plus x) z === pure2 plus x z *** QED + ) + === bind (bernoulli q) (pure2 plus x) + === bind (bernoulli q) (pure2 plus x) + === seqBind (bernoulli q) (pure2 plus) x + *** QED + +{-@ binsDistR :: p:Prob -> {q:Prob|p <= q} -> {n:PDouble|1 <= n} + -> {dist (kant distDouble) (bins' p q n) (bins q n) <= (n - 1) * (q - p)} + / [n, 0] @-} +binsDistR ::Prob -> Prob -> Double -> () +binsDistR p q n + = dist (kant d) (bins' p q n) (bins q n) + ? addBernoulliEq n p q + ? assert (bins' p q n == bind (bins p (n - 1)) (seqBind (bernoulli q) (pure2 plus))) + ? assert (bins q n == bind (bins q (n - 1)) (seqBind (bernoulli q) (pure2 plus))) + === dist (kant d) + (bind (bins p (n - 1)) (seqBind (bernoulli q) (pure2 plus))) + (bind (bins q (n - 1)) (seqBind (bernoulli q) (pure2 plus))) + ? commutative (bins p (n - 1)) (bernoulli q) (pure2 plus) + ? commutative (bins q (n - 1)) (bernoulli q) (pure2 plus) + === dist (kant d) + (bind (bernoulli q) (seqBind (bins p (n - 1)) (flip (pure2 plus)))) + (bind (bernoulli q) (seqBind (bins q (n - 1)) (flip (pure2 plus)))) + ? bindDistEq d + ((n - 1) * (q - p)) + (seqBind (bins p (n - 1)) (flip (pure2 plus))) (bernoulli q) + (seqBind (bins q (n - 1)) (flip (pure2 plus))) (bernoulli q) + (addBinsDist p q (n - 1)) + =<= (n - 1) * (q - p) + *** QED + where d = distDouble + +{-@ binsDist :: p:Prob -> {q:Prob|p <= q} -> n:PDouble + -> {dist (kant distDouble) (bins p n) (bins q n) <= n * (q - p)} + / [n, 1] @-} +binsDist :: Prob -> Prob -> Double -> () +binsDist p q n | n < 1.0 + = pureDist distDouble 0 0 + ? assert (0 <= n) + ? assert (0 <= (q - p)) + ? assert (dist (kant distDouble) (bins p n) (bins q n) <= n * (q - p)) +binsDist p q n + = dist (kant d) (bins p n) (bins q n) + ? triangularIneq (kant d) (bins p n) (bins' p q n) (bins q n) + ? assert (dist (kant d) (bins p n) (bins q n) + <= dist (kant d) (bins p n) (bins' p q n) + + dist (kant d) (bins' p q n) (bins q n)) + =<= dist (kant d) (bins p n) (bins' p q n) + + dist (kant d) (bins' p q n) (bins q n) + ? binsDistL p q n + =<= q - p + + dist (kant d) (bins' p q n) (bins q n) + ? binsDistR p q n + =<= q - p + (n - 1) * (q - p) + =<= n * (q - p) + *** QED + where d = distDouble + +{-@ reflect bins' @-} +{-@ bins' :: Prob -> Prob -> n:Double -> PrM Double @-} +bins' :: Double -> Double -> Double -> PrM Double +bins' _ q n | n < 1.0 = ppure 0 +bins' p q n = bind (bins p (n - 1)) (addBernoulli q (n - 1)) + +{-@ flipPlus :: x:Double -> {(flip (pure2 plus) x) == (ppure . (plus x))} @-} +flipPlus :: Double -> () +flipPlus x = extDouble (flip (pure2 plus) x) (ppure . (plus x)) (flipPlus' x) + +{-@ flipPlus' :: x:Double -> y:Double -> {(flip (pure2 plus) x y) == (ppure . (plus x)) (y)} @-} +flipPlus' :: Double -> Double -> () +flipPlus' _ _ = () + +{-@ assume extDouble :: f:(a -> b) -> g:(a -> b) + -> (x:a -> {v:() | f x == g x}) -> {f == g } @-} +extDouble :: (a -> b) -> (a -> b) -> (a -> ()) -> () +extDouble _ _ _ = () diff --git a/src/Data/Derivative.hs b/src/Data/Derivative.hs new file mode 100644 --- /dev/null +++ b/src/Data/Derivative.hs @@ -0,0 +1,5 @@ +module Data.Derivative where + +-- TODO: find implementation, e.g. Numeric.AD +grad :: (Double -> Double) -> (Double -> Double) +grad f x = 2 * x + 2 diff --git a/src/Data/Dist.hs b/src/Data/Dist.hs new file mode 100644 --- /dev/null +++ b/src/Data/Dist.hs @@ -0,0 +1,120 @@ +----------------------------------------------------------------- +-- | Distance as a desugared type class ------------------------- +----------------------------------------------------------------- + +{-@ LIQUID "--reflection" @-} +{-@ LIQUID "--ple-local" @-} + +module Data.Dist where + +import Prelude hiding (max) +import Language.Haskell.Liquid.ProofCombinators +import Misc.ProofCombinators +import Data.List + +----------------------------------------------------------------- +-- | class Dist a ----------------------------------------------- +----------------------------------------------------------------- +data Dist a = Dist { + dist :: a -> a -> Double + , distEq :: a -> () + , triangularIneq :: a -> a -> a -> () + , symmetry :: a -> a -> () + } + + +{-@ data Dist a = Dist { + dist :: a -> a -> {v:Double | 0.0 <= v } + , distEq :: a:a -> {dist a a == 0} + , triangularIneq :: x:a -> y:a -> z:a -> {dist x z <= dist x y + dist y z} + , symmetry :: a:a -> b:a -> {dist a b = dist b a} + } @-} + +-- TODO: define this +-- distFun :: Dist b -> Dist (a -> b) + +----------------------------------------------------------------- +-- | instance Dist Double --------------------------------------- +----------------------------------------------------------------- + +{-@ reflect distDouble@-} +distDouble :: Dist Double +distDouble = Dist distD distEqD triangularIneqD symmetryD + +{-@ ple distEqD @-} +{-@ reflect distEqD @-} +distEqD :: Double -> () +{-@ distEqD :: x:Double -> {distD x x == 0 } @-} +distEqD _ = () + +{-@ ple triangularIneqD @-} +{-@ reflect triangularIneqD @-} +{-@ triangularIneqD :: a:Double -> b:Double -> c:Double -> { distD a c <= distD a b + distD b c} @-} +triangularIneqD :: Double -> Double -> Double -> () +triangularIneqD _ _ _ = () + +{-@ ple symmetryD @-} +{-@ reflect symmetryD @-} +{-@ symmetryD :: a:Double -> b:Double -> {distD a b = distD b a} @-} +symmetryD :: Double -> Double -> () +symmetryD _ _ = () + +{-@ reflect distD @-} +{-@ distD :: Double -> Double -> {d:Double | 0.0 <= d } @-} +distD :: Double -> Double -> Double +distD x y = if x <= y then y - x else x - y + +----------------------------------------------------------------- +-- | instance Dist a => Dist (List a) --------------------------- +----------------------------------------------------------------- +-- Note the proof obligations hold, but this is not a real metric +-- since the two lists should have the same len +-- The following cannot type check +-- listDist :: Dist a -> Dist (List a) +-- listDist d = Dist (distList d) (distListEq d) (distListTri d) (distListSym d) + +{-@ type ListEq a XS = {ys:List a | llen ys == llen XS } @-} +{-@ reflect distList @-} +{-@ distList :: Dist a -> x:List a -> y:ListEq a {x} + -> {d:Double | 0 <= d } @-} +distList :: Dist a -> List a -> List a -> Double +distList d Nil _ = 0 +distList d _ Nil = 0 +distList d (Cons x xs) (Cons y ys) = max (dist d x y) (distList d xs ys) + +{-@ ple distListEq @-} +{-@ distListEq :: d:Dist a -> x:List a -> { distList d x x == 0 } @-} +distListEq :: Dist a -> List a -> () +distListEq d Nil = () +distListEq d (Cons x xs) = distEq d x ? distListEq d xs + +{-@ ple distListSym @-} +{-@ distListSym :: d:Dist a -> x:List a -> y:ListEq a {x} -> { distList d x y == distList d y x } @-} +distListSym :: Dist a -> List a -> List a -> () +distListSym d Nil _ = () +distListSym d _ Nil = () +distListSym d (Cons x xs) (Cons y ys) = symmetry d x y ? distListSym d xs ys + + +{-@ ple distListTri @-} +{-@ distListTri :: d:Dist a -> x:List a -> y:ListEq a {x} -> z:ListEq a {x} + -> { distList d x z <= distList d x y + distList d y z } @-} +distListTri :: Dist a -> List a -> List a -> List a -> () +distListTri d x@Nil y z = assert (distList d x z <= distList d x y + distList d y z) +distListTri d x y z@Nil = assert (distList d x z <= distList d x y + distList d y z) +distListTri d (Cons x xs) (Cons y ys) (Cons z zs) + = triangularIneq d x y z ? distListTri d xs ys zs + +----------------------------------------------------------------- +-- | Linearity on Doubles +-- | Does not type check forall a, so cannot just get axiomatized +----------------------------------------------------------------- + +{-@ ple linearity @-} +{-@ linearity :: k:{Double | 0 <= k } -> l:Double -> a:Double -> b:Double + -> { distD (k * a + l) (k * b + l) = k * distD a b} @-} +linearity :: Double -> Double -> Double -> Double -> () +linearity k l a b + | a <= b = assert (k * a + l <= k * b + l) + | otherwise = assert (distD (k * a + l) (k * b + l) == k * distD a b) + ? assert (k * a + l >= k * b + l) diff --git a/src/Data/List.hs b/src/Data/List.hs new file mode 100644 --- /dev/null +++ b/src/Data/List.hs @@ -0,0 +1,90 @@ +{-@ LIQUID "--reflection" @-} + +module Data.List where + +import Prelude hiding ( map + , max + , zipWith + , all + , foldr + ) + + +{-@ type SameLen L = {v:_|llen v = llen L} @-} +{-@ type ListN N = {v:_|llen v = N} @-} + +data List a = Nil | Cons a (List a) + deriving (Eq, Show) + +{-@ reflect consDouble @-} +{-@ consDouble :: Double -> xs:List Double -> {v:List Double | llen v == llen xs + 1 } @-} +consDouble :: Double -> List Double -> List Double +consDouble = Cons + +{-@ measure llen @-} +{-@ llen :: List a -> Nat @-} +llen :: List a -> Int +llen Nil = 0 +llen (Cons _ xs) = 1 + llen xs + +{-@ type Idx V = {i:Int | 0 <= i && i < llen V} @-} + +{-@ reflect at @-} +{-@ at :: xs:List a -> Idx xs -> a @-} +at :: List a -> Int -> a +at (Cons x _) i | i == 0 = x +at (Cons _ xs) i = at xs (i - 1) + +{-@ reflect range @-} +{-@ range :: i:Nat -> len:Nat -> {v:List {j:Nat|j < i + len}|llen v = len} / [len] @-} +range :: Int -> Int -> List Int +range _ 0 = Nil +range i len = Cons i (range (i + 1) (len - 1)) + +{-@ reflect map @-} +{-@ map :: (a -> b) -> xs:List a -> {ys:List b|llen ys = llen xs} @-} +map :: (a -> b) -> List a -> List b +map f Nil = Nil +map f (Cons x xs) = Cons (f x) (map f xs) + +zipWith :: (a -> b -> c) -> List a -> List b -> List c +zipWith _ Nil _ = Nil +zipWith _ _ Nil = Nil +zipWith f (Cons x xs) (Cons x' xs') = Cons (f x x') (zipWith f xs xs') + +all :: List Bool -> Bool +all Nil = True +all (Cons x xs) = x && all xs + +{-@ reflect max @-} +max :: Double -> Double -> Double +max a b = if a < b then b else a + +{-@ reflect pow @-} +{-@ pow :: {v:Double|v >= 0} -> i:Nat -> {v:Double|v >= 0} / [i] @-} +pow :: Double -> Int -> Double +pow x 0 = 1 +pow x i = x * pow x (i - 1) + +{-@ reflect ap @-} +ap :: List (a -> b) -> List a -> List b +ap _ Nil = Nil +ap Nil _ = Nil +ap (Cons f fs) (Cons x xs) = Cons (f x) (ap fs xs) + +{-@ reflect zip3With @-} +zip3With :: (a -> b -> c -> d) -> List a -> List b -> List c -> List d +zip3With _ Nil _ _ = Nil +zip3With _ _ Nil _ = Nil +zip3With _ _ _ Nil = Nil +zip3With f (Cons a as) (Cons b bs) (Cons c cs) = Cons (f a b c) (zip3With f as bs cs) + +{-@ zip4 :: as:List a -> {bs:List b|llen bs = llen as} -> {cs:List c|llen cs = llen as} -> {ds:List d|llen ds = llen as} -> List (a, b, c, d) @-} +zip4 :: List a -> List b -> List c -> List d -> List (a, b, c, d) +zip4 Nil Nil Nil Nil = Nil +zip4 (Cons a as) (Cons b bs) (Cons c cs) (Cons d ds) = Cons (a, b, c, d) (zip4 as bs cs ds) + +{-@ reflect foldr @-} +foldr :: (a -> b -> b) -> b -> List a -> b +foldr _ z Nil = z +foldr f z (Cons x xs) = f x (foldr f z xs) diff --git a/src/Examples/ExpDist.hs b/src/Examples/ExpDist.hs new file mode 100644 --- /dev/null +++ b/src/Examples/ExpDist.hs @@ -0,0 +1,40 @@ +{-@ LIQUID "--reflection" @-} +{-@ LIQUID "--ple" @-} + +module Examples.ExpDist where + +import Monad.PrM +import Data.Dist +import Data.List + +import Monad.PrM.Relational.TCB.EDist +import Misc.ProofCombinators + +import Prelude hiding ( map + , max + , repeat + , foldr + , fmap + , mapM + , iterate + , uncurry + ) + +{-@ relationalu :: d:Dist a -> xs:[a] -> { dist (kant d) (unif xs) (unif xs) == 0} @-} +relationalu :: Dist a -> [a] -> () +relationalu d xs = unifDist d xs xs + + +-- Attention: In the Haskell code you need to write 4.0 instead of just 4 to avoid implicit conversion +{-@ exDistPure :: () -> {dist (kant distDouble) (ppure 4.0) (ppure 2.0) <= 2.0 } @-} +exDistPure :: () -> () +exDistPure _ = pureDist distDouble 4.0 2.0 + +{-@ ex2DistPure :: p:Prob -> xs:{[Double] | 0 < len xs } + -> {dist (kant distDouble) (choice p (ppure 4.0) (unif xs)) (choice p (ppure 2.0) (unif xs)) <= p * 2.0 } @-} +ex2DistPure :: Prob -> [Double] -> () +ex2DistPure p xs + = relationalu distDouble xs `const` + exDistPure () `const` + choiceDist distDouble p (ppure 4.0) (unif xs) p (ppure 2.0) (unif xs) + diff --git a/src/Misc/ProofCombinators.hs b/src/Misc/ProofCombinators.hs new file mode 100644 --- /dev/null +++ b/src/Misc/ProofCombinators.hs @@ -0,0 +1,10 @@ +module Misc.ProofCombinators where + + +{-@ assert :: {b:Bool | b} -> {v:() | b } @-} +assert :: Bool -> () +assert _ = () + +{-@ assume assume :: b:Bool -> {v:() | b } @-} +assume :: Bool -> () +assume _ = () diff --git a/src/Monad/PrM.hs b/src/Monad/PrM.hs new file mode 100644 --- /dev/null +++ b/src/Monad/PrM.hs @@ -0,0 +1,110 @@ +----------------------------------------------------------------- +-- | Implementation of the PrM monad as a wrapper ------------ +-- | This module includes only executable code ------------ +----------------------------------------------------------------- + +{-@ LIQUID "--reflection" @-} +module Monad.PrM where + +import Data.Dist +import Data.List hiding (all) + +import Prelude hiding (max, mapM) +import Numeric.Probability.Distribution hiding (Cons, cons) + +{-@ type Prob = {v:Double| 0 <= v && v <= 1} @-} +type Prob = Double + +type PrM a = T Prob a + +{-@ measure Monad.PrM.bind :: PrM a -> (a -> PrM b) -> PrM b @-} +{-@ assume bind :: x1:PrM a -> x2:(a -> PrM b) -> {v:PrM b | v = bind x1 x2 } @-} +bind :: PrM a -> (a -> PrM b) -> PrM b +bind = (>>=) + +{-@ measure Monad.PrM.ppure :: a -> PrM a @-} +{-@ ppure :: x:a -> {v:PrM a | v = Monad.PrM.ppure x } @-} +ppure :: a -> PrM a +ppure = pure + +{-@ measure Monad.PrM.liftA2 :: (a -> b -> c) -> PrM a -> PrM b -> PrM c @-} +{-@ liftA2 :: f:(a -> b -> c) -> x:PrM a -> y:PrM b -> {v:PrM c | v = Monad.PrM.liftA2 f x y} @-} +liftA2 :: (a -> b -> c) -> PrM a -> PrM b -> PrM c +liftA2 f a b = do x <- a + y <- b + ppure (f x y) + +{-@ reflect fmap @-} +fmap :: (a -> b) -> PrM a -> PrM b +fmap f a = bind a (ppure . f) + +{-@ measure Monad.PrM.choice :: Prob -> PrM a -> PrM a -> PrM a @-} +{-@ assume choice :: x1:Prob -> x2:PrM a -> x3:PrM a -> {v:PrM a | v == choice x1 x2 x3 } @-} +choice :: Prob -> PrM a -> PrM a -> PrM a +choice p x y = cond (fromFreqs [(True, p), (False, 1 - p)]) x y + +{-@ measure Monad.PrM.bernoulli :: Prob -> PrM Double @-} +{-@ assume bernoulli :: p:Prob -> {v:PrM {n:Double | 0 <= n && n <= 1}| v == bernoulli p } @-} +bernoulli :: Prob -> PrM Double +bernoulli p = fromFreqs [(1, p), (0, 1 - p)] + +{-@ reflect unif @-} +{-@ unif :: {xs:[a]|0 < len xs} -> PrM a @-} +unif :: [a] -> PrM a +unif [a] = ppure a +unif (x:xs) = choice (1 `mydiv` fromIntegral (len xs + 1)) (ppure x) (unif xs) + +{-@ measure Monad.PrM.lift :: (a -> b -> Bool) -> PrM a -> PrM b -> Bool @-} +{-@ assume lift :: p1:(a -> b -> Bool) -> x1:PrM a -> x2:PrM b + -> {v:Bool | v == Monad.PrM.lift p1 x1 x2 } @-} +lift :: (a -> b -> Bool) -> PrM a -> PrM b -> Bool +lift p e1 e2 = and (fst <$> (decons act)) + where act = do x <- e1 + y <- e2 + return (p x y) + +----------------------------------------------------------------- +-- | mapM: Standard monadic mapM instantiated for LH limitations +----------------------------------------------------------------- + +{-@ reflect mapM @-} +{-@ mapM :: (a -> PrM Double) -> xs:List a -> PrM ({ys:List Double| llen ys = llen xs }) @-} +mapM :: (a -> PrM Double) -> List a -> PrM (List Double) +mapM _ Nil = ppureDouble Nil +mapM f (Cons x xs) = bind (f x) (cons (llen xs) (mapM f xs)) + +----------------------------------------------------------------- +-- | Helper Definitions for Reflection +----------------------------------------------------------------- + +{-@ reflect len @-} +len :: [a] -> Int +len [] = 0 +len (_:xs) = 1 + len xs + +{-@ reflect mydiv @-} +{-@ mydiv :: Double -> {i:Double | i /= 0 } -> Double @-} +mydiv :: Double -> Double -> Double +mydiv x y = x / y + +{-@ reflect ppureDouble @-} +{-@ ppureDouble :: xs:List Double -> PrM ({v:List Double | llen v == llen xs}) @-} +ppureDouble :: List Double -> PrM (List Double) +ppureDouble x = ppure x + +{-@ reflect cons @-} +{-@ cons :: n:Nat -> PrM ({xs:List Double | llen xs == n}) -> Double -> PrM ({v:List Double | llen v = n + 1}) @-} +cons :: Int -> PrM (List Double) -> Double -> PrM (List Double) +cons n xs x = bind xs (ppure `o` (consDouble x)) + +{-@ reflect o @-} +o :: (b -> c) -> (a -> b) -> a -> c +o g f x = g (f x) + +{-@ reflect seqBind @-} +seqBind :: PrM b -> (a -> b -> PrM c) -> a -> PrM c +seqBind u f x = bind u (f x) + +{-@ reflect flip @-} +flip :: (a -> b -> c) -> b -> a -> c +flip f x y = f y x diff --git a/src/Monad/PrM/Laws.hs b/src/Monad/PrM/Laws.hs new file mode 100644 --- /dev/null +++ b/src/Monad/PrM/Laws.hs @@ -0,0 +1,18 @@ +----------------------------------------------------------------- +-- | Monad Laws for the PrM monad ----------------------------- +----------------------------------------------------------------- + +{-@ LIQUID "--reflection" @-} +module Monad.PrM.Laws where + +import Monad.PrM + +{-@ assume leftId :: x:a -> f:(a -> PrM b) -> { bind (ppure x) f = f x } @-} +leftId :: a -> (a -> PrM b) -> () +leftId _ _ = () + +{-@ assume commutative :: e:PrM a -> u:PrM b -> f:(a -> b -> PrM c) + -> {bind e (seqBind u f) + = bind u (seqBind e (flip f))} @-} +commutative :: PrM a -> PrM b -> (a -> b -> PrM c) -> () +commutative _ _ _ = () diff --git a/src/Monad/PrM/Predicates.hs b/src/Monad/PrM/Predicates.hs new file mode 100644 --- /dev/null +++ b/src/Monad/PrM/Predicates.hs @@ -0,0 +1,55 @@ +----------------------------------------------------------------- +-- | Reflected Predicates (required for lifting) ---------------- +----------------------------------------------------------------- + +{-@ LIQUID "--reflection" @-} +{-@ LIQUID "--ple-local" @-} + +module Monad.PrM.Predicates where + +import Data.Dist +import Data.List + +{-@ reflect trueP @-} +trueP :: a -> a -> Bool +trueP _ _ = True + +{-@ reflect bounded @-} +{-@ bounded :: Double -> x:List Double -> ListEq Double {x} -> Bool @-} +bounded :: Double -> List Double -> List Double -> Bool +bounded m v1 v2 = distList distDouble v1 v2 <= m && llen v1 == llen v2 + +{-@ reflect boundedD @-} +boundedD :: Dist a -> Double -> a -> a -> Bool +boundedD d m v1 v2 = dist d v1 v2 <= m + +{-@ reflect bounded' @-} +bounded' :: Double -> Double -> Double -> Bool +bounded' m x1 x2 = distD x1 x2 <= m + +{-@ reflect eqP @-} +eqP :: Eq a => a -> a -> Bool +eqP = (==) + + +{-@ reflect leDoubleP @-} +{-@ leDoubleP :: x:Double -> y:Double -> {v:Bool|v <=> (x <= y)} @-} +leDoubleP :: Double -> Double -> Bool +leDoubleP x y = x <= y + +{-@ reflect impP @-} +{-@ impP :: x:Bool -> y:Bool -> {v:Bool|v <=> (x => y)} @-} +impP :: Bool -> Bool -> Bool +impP True False = False +impP _ _ = True + +{-@ reflect leIntP @-} +{-@ leIntP :: x:Int -> y:Int -> {v:Bool|v <=> (x <= y)} @-} +leIntP :: Int -> Int -> Bool +leIntP x y = x <= y + +-- Properties on Predicates +{-@ ple boundedNil @-} +{-@ boundedNil :: {m:_|0 <= m} -> {bounded m Nil Nil} @-} +boundedNil :: Double -> () +boundedNil _ = () diff --git a/src/Monad/PrM/Relational/TCB/EDist.hs b/src/Monad/PrM/Relational/TCB/EDist.hs new file mode 100644 --- /dev/null +++ b/src/Monad/PrM/Relational/TCB/EDist.hs @@ -0,0 +1,77 @@ +----------------------------------------------------------------- +-- | Expected Distance Specifications for PrM Primitives ------ +----------------------------------------------------------------- + +{-@ LIQUID "--reflection" @-} + +module Monad.PrM.Relational.TCB.EDist where + +import Data.Dist +import Data.List +import Monad.PrM +import Monad.PrM.Relational.TCB.Spec + +{-@ measure Monad.PrM.Relational.TCB.EDist.kant :: Dist a -> Dist (PrM a) @-} +{-@ assume kant :: d:Dist a -> {dd:Dist (PrM a) | dd = Monad.PrM.Relational.TCB.EDist.kant d } @-} +kant :: Dist a -> Dist (PrM a) +kant = undefined + +{-@ reflect edist @-} +{-@ edist :: Dist a -> PrM a -> PrM a -> {v:Double | 0 <= v } @-} +edist :: Dist a -> PrM a -> PrM a -> Double +edist d = dist (kant d) + +{-@ assume pureDist :: d:Dist a -> x1:a -> x2:a + -> { dist (kant d) (ppure x1) (ppure x2) = dist d x1 x2} @-} +pureDist :: Dist a -> a -> a -> () +pureDist _ _ _ = () + +{-@ assume bindDist :: d:Dist b -> m:Double -> p:(a -> a -> Bool) + -> f1:(a -> PrM b) -> e1:PrM a + -> f2:(a -> PrM b) -> e2:{PrM a | lift p e1 e2} + -> lemma:(x1:a -> {x2:a| p x1 x2 } + -> { dist (kant d) (f1 x1) (f2 x2) <= m}) + -> { dist (kant d) (bind e1 f1) (bind e2 f2) <= m } @-} +bindDist :: Dist b -> Double -> (a -> a -> Bool) -> (a -> PrM b) -> PrM a -> (a -> PrM b) -> PrM a -> (a -> a -> ()) -> () +bindDist _ _ _ _ _ _ _ _ = () + +{-@ assume fmapDist :: da:Dist a -> db:Dist b + -> m:Double + -> f1:(a -> b) -> e1:PrM a + -> f2:(a -> b) -> e2:PrM a + -> (x1:a -> x2:a -> { dist db (f1 x1) (f2 x2) <= dist da x1 x2 + m}) + -> { dist (kant db) (fmap f1 e1) (fmap f2 e2) <= dist (kant da) e1 e2 + m } @-} +fmapDist :: Dist a -> Dist b -> Double -> (a -> b) -> PrM a -> (a -> b) -> PrM a -> (a -> a -> ()) -> () +fmapDist _ _ _ _ _ _ _ _ = () + +{-@ assume liftA2Dist :: da:Dist a -> db:Dist b -> dc:Dist c + -> ma:Double -> ka:Double -> mb:Double -> kb:Double -> m:Double + -> f1:(a -> b -> c) -> e1:PrM a -> u1:PrM b + -> f2:(a -> b -> c) -> e2:PrM a -> u2:PrM b + -> {_:()|dist (kant da) e1 e2 <= ka} + -> {_:()|dist (kant db) u1 u2 <= kb} + -> (x1:a -> y1:b -> x2:a -> y2:b + -> {dist dc (f1 x1 y1) (f2 x2 y2) <= ma * dist da x1 x2 + mb * dist db y1 y2 + m}) + -> {dist (kant dc) (liftA2 f1 e1 u1) (liftA2 f2 e2 u2) <= ma * ka + mb * kb + m} @-} +liftA2Dist :: Dist a -> Dist b -> Dist c -> Double -> Double -> Double -> Double -> Double + -> (a -> b -> c) -> PrM a -> PrM b -> (a -> b -> c) -> PrM a -> PrM b + -> () -> () -> (a -> b -> a -> b -> ()) + -> () +liftA2Dist _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ = () + +{-@ assume unifDist :: d:Dist a -> xsl:[a] -> xsr:{[a] | xsl == xsr} + -> { dist (kant d) (unif xsl) (unif xsr) == 0 } @-} +unifDist :: Dist a -> [a] -> [a] -> () +unifDist _ _ _ = () + +{-@ assume choiceDist :: d:Dist a -> p:Prob -> e1:PrM a -> e1':PrM a + -> q:{Prob | p = q } -> e2:PrM a -> e2':PrM a + -> { dist (kant d) (choice p e1 e1') (choice q e2 e2') <= p * (dist (kant d) e1 e2) + (1.0 - p) * (dist (kant d) e1' e2')} @-} +choiceDist :: Dist a -> Prob -> PrM a -> PrM a -> Prob -> PrM a -> PrM a -> () +choiceDist _ _ _ _ _ _ _ = () + +{-@ assume bernoulliDist :: d:Dist Double -> p:Prob -> {q:Prob | p <= q} + -> {dist (kant d) (bernoulli p) (bernoulli q) <= dist d 1 0 * (q - p)} @-} +bernoulliDist :: Dist Double -> Prob -> Prob -> () +bernoulliDist d p q = () + where _ = dist (kant d) (bernoulli p) (bernoulli q) <= dist d 1 0 * (q - p) diff --git a/src/Monad/PrM/Relational/TCB/Spec.hs b/src/Monad/PrM/Relational/TCB/Spec.hs new file mode 100644 --- /dev/null +++ b/src/Monad/PrM/Relational/TCB/Spec.hs @@ -0,0 +1,43 @@ +----------------------------------------------------------------- +-- | Relational Specifications for PrM Primitives ------------- +----------------------------------------------------------------- + +{-@ LIQUID "--reflection" @-} + +module Monad.PrM.Relational.TCB.Spec where + +import Monad.PrM +import Monad.PrM.Predicates + +{-@ assume pureSpec :: p:(a -> b -> Bool) + -> x1:a -> x2:b -> {_:_|p x1 x2} + -> {lift p (ppure x1) (ppure x2)} @-} +pureSpec :: (a -> b -> Bool) -> a -> b -> () -> () +pureSpec _ _ _ _ = () + +{-@ assume bindSpec :: p:(b -> b -> Bool) -> q:(a -> a -> Bool) + -> e1:PrM a -> f1:(a -> PrM b) + -> e2:PrM a -> f2:(a -> PrM b) + -> {_:()|lift q e1 e2} + -> (x1:a -> {x2:a|q x1 x2} -> {lift p (f1 x1) (f2 x2)}) + -> {lift p (bind e1 f1) (bind e2 f2)} @-} +bindSpec :: (b -> b -> Bool) -> (a -> a -> Bool) -> + PrM a -> (a -> PrM b) -> PrM a -> (a -> PrM b) -> + () -> + (a -> a -> ()) -> + () +bindSpec _ _ _ _ _ _ _ _ = () + +{-@ assume bernoulliSpec :: p:Prob -> {q:Prob| leDoubleP p q} + -> {lift leDoubleP (bernoulli p) (bernoulli q)} @-} +bernoulliSpec :: Prob -> Prob -> () +bernoulliSpec _ _ = () + +{-@ assume liftSpec :: e:PrM a -> {lift eqP e e} @-} +liftSpec :: PrM a -> () +liftSpec _ = () + + +{-@ assume liftTrue :: e1:PrM a -> e2:PrM a -> {lift trueP e1 e2} @-} +liftTrue :: PrM a -> PrM a -> () +liftTrue _ _ = () diff --git a/src/Monad/PrM/Relational/Theorems.hs b/src/Monad/PrM/Relational/Theorems.hs new file mode 100644 --- /dev/null +++ b/src/Monad/PrM/Relational/Theorems.hs @@ -0,0 +1,101 @@ +----------------------------------------------------------------- +-- | Proved Theorems for Relational Properties: mapMSpec ------ +----------------------------------------------------------------- + +{-@ LIQUID "--reflection" @-} +{-@ LIQUID "--ple" @-} + +module Monad.PrM.Relational.Theorems where + +import Monad.PrM +import Data.Dist +import Data.List +import Prelude hiding (max, mapM) + +import Monad.PrM.Relational.TCB.Spec +import Monad.PrM.Relational.TCB.EDist +import Monad.PrM.Predicates + +import Language.Haskell.Liquid.ProofCombinators +import Misc.ProofCombinators + + +------------------------------------------------------- +-- | bindDistEq when the bind args are BijCoupling --- +------------------------------------------------------- + +{-@ predicate BijCoupling X Y = X = Y @-} +{-@ bindDistEq :: d:Dist b -> m:Double + -> f1:(a -> PrM b) -> e1:PrM a + -> f2:(a -> PrM b) -> e2:{PrM a | BijCoupling e1 e2 } + -> (x:a -> { dist (kant d) (f1 x) (f2 x) <= m}) + -> { dist (kant d) (bind e1 f1) (bind e2 f2) <= m } @-} +bindDistEq :: Eq a => Dist b -> Double -> (a -> PrM b) -> PrM a -> (a -> PrM b) -> PrM a -> (a -> ()) -> () +bindDistEq d m f1 e1 f2 e2 lemma = + bindDist d m eqP f1 e1 f2 (e2 `const` liftSpec e2) + (makeTwoArg d m f1 f2 lemma) + +{-@ makeTwoArg :: d:Dist b -> m:Double -> f1:(a -> PrM b) -> f2:(a -> PrM b) + -> (x:a -> {v:() | dist (kant d) (f1 x) (f2 x) <= m}) + -> (x:a -> y:{a | eqP x y} -> { dist (kant d) (f1 x) (f2 y) <= m}) @-} +makeTwoArg :: Dist b -> Double -> (a -> PrM b) -> (a -> PrM b) -> (a -> ()) + -> (a -> a -> ()) +makeTwoArg d m f1 f2 lemma x y = lemma x + +------------------------------------------------------- +-- | mapM Spec ---------------------------------------- +------------------------------------------------------- + +{-@ mapMSpec :: {m:_|0 <= m} + -> f1:(a -> PrM Double) -> f2:(a -> PrM Double) + -> is:List a + -> (i:a -> {lift (bounded' m) (f1 i) (f2 i)}) + -> {lift (bounded m) (mapM f1 is) (mapM f2 is)} / [llen is, 0] @-} +mapMSpec :: Double -> (a -> PrM Double) -> (a -> PrM Double) -> List a + -> (a -> ()) + -> () +mapMSpec m f1 f2 is@Nil lemma + = pureSpec (bounded m) Nil Nil (boundedNil m) +mapMSpec m f1 f2 (Cons i is) lemma + = bindSpec (bounded m) (bounded' m) + (f1 i) (cons (llen is) (mapM f1 is)) + (f2 i) (cons (llen is) (mapM f2 is)) + (lemma i) + (consBindLemma m f1 f2 is lemma) + +{-@ consLemma :: m:_ -> r1:_ -> rs1:_ -> {r2:_|bounded' m r1 r2} -> {rs2:_|llen rs1 = llen rs2 && bounded m rs1 rs2} + -> {bounded m (Cons r1 rs1) (Cons r2 rs2)} @-} +consLemma :: Double -> Double -> List Double -> Double -> List Double -> () +consLemma m r1 rs1 r2 rs2 = () + +{-@ consBindLemma :: {m:_|0 <= m} -> f1:_ -> f2:_ -> is:_ + -> (i:a -> {lift (bounded' m) (f1 i) (f2 i)}) + -> r1:_ + -> {r2:_|bounded' m r1 r2} + -> {lift (bounded m) + ((cons (llen is) (mapM f1 is)) (r1)) + ((cons (llen is) (mapM f2 is)) (r2))} / [llen is, 1] @-} +consBindLemma :: Double -> (a -> PrM Double) -> (a -> PrM Double) + -> List a + -> (a -> ()) + -> Double -> Double + -> () +consBindLemma m f1 f2 is lemma r1 r2 + = bindSpec (bounded m) (bounded m) + (mapM f1 is) (ppure `o` (consDouble r1)) + (mapM f2 is) (ppure `o` (consDouble r2)) + (mapMSpec m f1 f2 is lemma) + (pureLemma m r1 r2 f1 f2 is) + +{-@ pureLemma :: {m:_|0 <= m} + -> r1:_ -> {r2:_|bounded' m r1 r2} + -> f1:_ -> f2:_ -> is:_ + -> rs1:_ -> rs2:{_|bounded m rs1 rs2} + -> {lift (bounded m) (o ppure (consDouble r1) rs1) + (o ppure (consDouble r2) rs2)} @-} +pureLemma :: Double -> Double -> Double -> (a -> PrM Double) -> (a -> PrM Double) + -> List a -> List Double -> List Double -> () +pureLemma m r1 r2 f1 f2 is rs1 rs2 = pureSpec (bounded m) + (Cons r1 rs1) (Cons r2 rs2) + (consLemma m r1 rs1 r2 rs2) + diff --git a/src/SGD/SGD.hs b/src/SGD/SGD.hs new file mode 100644 --- /dev/null +++ b/src/SGD/SGD.hs @@ -0,0 +1,90 @@ +{-@ LIQUID "--reflection" @-} + +module SGD.SGD where + +import Prelude hiding ( head, tail, sum) +import Monad.PrM +import Data.Dist +import Data.Derivative + +{-@ type StepSize = {v:Double | 0.0 <= v } @-} +type StepSize = Double +{-@ data StepSizes = SSEmp | SS StepSize StepSizes @-} +data StepSizes = SSEmp | SS StepSize StepSizes +type DataPoint = (Double, Double) +type Weight = Double +type LossFunction = DataPoint -> Weight -> Double + +type Set a = [a] +{-@ type DataSet = {v:Set DataPoint| 1 < lend v && 1 < len v } @-} +type DataSet = Set DataPoint +type DataPrM = PrM DataPoint + + +{-@ reflect sgd @-} +{-@ sgd :: zs:{DataSet | 1 < len zs && 1 < lend zs } -> Weight -> ss:StepSizes -> LossFunction + -> PrM Weight / [ sslen ss, 0 ] @-} +sgd :: DataSet -> Weight -> StepSizes -> LossFunction -> PrM Weight +sgd _ w0 SSEmp _ = ppure w0 +sgd zs w0 (SS α a) f = + choice (one / lend zs) + (bind uhead (sgdRecUpd zs w0 α a f)) + (bind utail (sgdRecUpd zs w0 α a f)) + where + uhead = ppure (head zs) + utail = unif (tail zs) + + +{-@ reflect sgdRecUpd @-} +{-@ sgdRecUpd :: zs:{DataSet | 1 < len zs && 1 < lend zs } -> Weight -> StepSize -> ss:StepSizes -> LossFunction + -> DataPoint -> PrM Weight / [ sslen ss, 1 ] @-} +sgdRecUpd :: DataSet -> Weight -> StepSize -> StepSizes -> LossFunction -> DataPoint -> PrM Weight +sgdRecUpd zs w0 α a f z = bind (sgd zs w0 a f) (pureUpdate z α f) + +{-@ reflect pureUpdate @-} +{-@ pureUpdate :: DataPoint -> StepSize -> LossFunction -> Weight -> PrM Weight @-} +pureUpdate :: DataPoint -> StepSize -> LossFunction -> Weight -> PrM Weight +pureUpdate zs a f = ppure . update zs a f + + +{-@ measure SGD.SGD.update :: DataPoint -> StepSize -> LossFunction -> Weight -> Weight @-} +{-@ update :: x1:DataPoint -> x2:StepSize -> x3:LossFunction -> x4:Weight + -> {v:Weight | v = SGD.SGD.update x1 x2 x3 x4 } @-} +update :: DataPoint -> StepSize -> LossFunction -> Weight -> Weight +update z α f w = w - α * (grad (f z) w) + + +------------------------------------------------------------------------------- +-- | Helper Definitions ------------------------------------------------------- +------------------------------------------------------------------------------- + + +{-@ measure lend @-} +{-@ lend :: xs:[a] -> {v:Double| 0.0 <= v } @-} +lend :: [a] -> Double +lend [] = 0 +lend (_ : xs) = 1 + lend xs + + +{-@ reflect one @-} +{-@ one :: {v:Double| v = 1.0 } @-} +one :: Double +one = 1 + + + +{-@ reflect head @-} +{-@ head :: {xs:[a] | len xs > 0 } -> a @-} +head :: [a] -> a +head (z : _) = z + +{-@ reflect tail @-} +{-@ tail :: {xs:[a] | len xs > 0 } -> {v:[a] | len v == len xs - 1 && lend v == lend xs - 1 } @-} +tail :: [a] -> [a] +tail (_ : zs) = zs + +{-@ measure sslen @-} +sslen :: StepSizes -> Int +{-@ sslen :: StepSizes -> Nat @-} +sslen SSEmp = 0 +sslen (SS _ ss) = 1 + sslen ss diff --git a/src/SGD/Theorem.hs b/src/SGD/Theorem.hs new file mode 100644 --- /dev/null +++ b/src/SGD/Theorem.hs @@ -0,0 +1,196 @@ +{-@ LIQUID "--reflection" @-} +{-@ LIQUID "--fast" @-} +{-@ LIQUID "--ple-local" @-} + +module SGD.Theorem where + +import Prelude hiding ( head + , tail + , sum + , fmap + ) +import Language.Haskell.Liquid.ProofCombinators + +import Misc.ProofCombinators + +import Monad.PrM +import Monad.PrM.Laws +import Monad.PrM.Relational.TCB.EDist +import Monad.PrM.Relational.Theorems (bindDistEq) +import Data.Dist +import SGD.SGD + + +{-@ measure SGD.Theorem.lip :: Double @-} +{-@ assume lip :: {v:Double|SGD.Theorem.lip = v && v >= 0 } @-} +lip :: Double +lip = 10 + + +{-@ assume relationalupdatep :: d:Dist Double -> z1:DataPoint -> α1:StepSize -> f1:LossFunction + -> z2:DataPoint -> {α2:StepSize|α1 = α2} -> {f2:LossFunction|f1 = f2} + -> ws1:Weight -> ws2:Weight -> + {dist d (update z1 α1 f1 ws1) (update z2 α2 f2 ws2) = dist d ws1 ws2 + 2.0 * lip * α1 } @-} +relationalupdatep :: Dist Double -> DataPoint -> StepSize -> LossFunction -> DataPoint -> StepSize -> LossFunction -> Weight -> Weight -> () +relationalupdatep _ _ _ _ _ _ _ _ _ = () + + +{-@ assume relationalupdateq :: d:Dist Double -> z1:DataPoint -> α1:StepSize -> f1:LossFunction + -> {z2:DataPoint|z1 = z2} -> {α2:StepSize|α1 = α2} -> {f2:LossFunction|f1 = f2} + -> ws1:Weight -> ws2:Weight -> + {dist d (update z1 α1 f1 ws1) (update z2 α2 f2 ws2) = dist d ws1 ws2} @-} +relationalupdateq :: Dist Double -> DataPoint -> StepSize -> LossFunction -> DataPoint -> StepSize -> LossFunction -> Weight -> Weight -> () +relationalupdateq = undefined + + +{-@ reflect sum @-} +{-@ sum :: StepSizes -> {v:StepSize | 0.0 <= v } @-} +sum :: StepSizes -> Double +sum SSEmp = 0 +sum (SS a as) = a + sum as + +{-@ reflect estab @-} +{-@ estab :: DataSet -> StepSizes -> {v:Double | 0.0 <= v} @-} +estab :: DataSet -> StepSizes -> Double +estab zs as = 2.0 * lip / (lend zs) * sum as + +{-@ ple estabEmp @-} +estabEmp :: DataSet -> () +{-@ estabEmp :: zs:DataSet -> {estab zs SSEmp == 0.0} @-} +estabEmp zs = + estab zs SSEmp + === 2.0 / (lend zs) * sum SSEmp + *** QED + +{-@ ple estabconsR @-} +{-@ measure Theorem.estabconsR :: DataSet -> StepSize -> StepSizes -> () @-} +{-@ estabconsR :: zs:{DataSet | lend zs /= 0} -> x:StepSize -> xs:StepSizes + -> { estab zs (SS x xs) == 2.0 * lip * x * (one / lend zs) + estab zs xs } @-} +estabconsR :: DataSet -> StepSize -> StepSizes -> () +estabconsR zs x xs + = estab zs (SS x xs) + === 2.0 * lip / (lend zs) * sum (SS x xs) + === 2.0 * lip * x * (one / lend zs) + estab zs xs + *** QED + +{-@ ple thm @-} +{-@ thm :: d:Dist Double -> zs1:DataSet -> ws1:Weight -> α1:StepSizes -> f1:LossFunction -> + zs2:{DataSet | lend zs1 == lend zs2 && tail zs1 = tail zs2} -> + ws2:Weight -> {α2:StepSizes| α2 = α1} -> {f2:LossFunction|f1 = f2} -> + { dist (kant d) (sgd zs1 ws1 α1 f1) (sgd zs2 ws2 α2 f2) <= dist d ws1 ws2 + estab zs1 α1} / [sslen α1, 0]@-} +thm :: Dist Double -> DataSet -> Weight -> StepSizes -> LossFunction -> DataSet -> Weight -> StepSizes -> LossFunction -> () +thm d zs1 ws1 α1@SSEmp f1 zs2 ws2 α2@SSEmp f2 = + dist (kant d) (sgd zs1 ws1 α1 f1) (sgd zs2 ws2 α2 f2) + === dist (kant d) (ppure ws1) (ppure ws2) + ? pureDist d ws1 ws2 + === dist d ws1 ws2 + ? estabEmp zs1 + === dist d ws1 ws2 + estab zs1 α1 + *** QED + +thm d zs1 ws1 as1@(SS α1 a1) f1 zs2 ws2 as2@(SS α2 a2) f2 = + dist (kant d) (sgd zs1 ws1 as1 f1) (sgd zs2 ws2 as2 f2) + === dist (kant d) + (choice (one / lend zs1) (bind uhead1 sgdRec1) (bind utail1 sgdRec1)) + (choice (one / lend zs2) (bind uhead2 sgdRec2) (bind utail2 sgdRec2)) + ? choiceDist d (one / lend zs1) (bind uhead1 sgdRec1) (bind utail1 sgdRec1) + (one / lend zs2) (bind uhead2 sgdRec2) (bind utail2 sgdRec2) + + =<= (one / lend zs1) * (dist (kant d) (bind uhead1 sgdRec1) (bind uhead2 sgdRec2)) + + (1 - (one / lend zs1)) * (dist (kant d) (bind utail1 sgdRec1) (bind utail2 sgdRec2)) + ? leftId (head zs1) sgdRec1 + ? leftId (head zs2) sgdRec2 + + =<= (one / lend zs1) * (dist (kant d) (sgdRec1 (head zs1)) (sgdRec2 (head zs2))) + + (1 - (one / lend zs1)) * (dist (kant d) (bind utail1 sgdRec1) (bind utail2 sgdRec2)) + + =<= (one / lend zs1) * (dist (kant d) (bind (sgd zs1 ws1 a1 f1) pureUpd1) + (bind (sgd zs2 ws2 a2 f2) pureUpd2)) + + (1 - (one / lend zs1)) * (dist (kant d) (bind utail1 sgdRec1) (bind utail2 sgdRec2)) + ? pureUpdateEq (head zs1) α1 f1 + ? pureUpdateEq (head zs2) α2 f2 + + =<= (one / lend zs1) * (dist (kant d) (bind (sgd zs1 ws1 a1 f1) (ppure . update (head zs1) α1 f1 )) + (bind (sgd zs2 ws2 a2 f2) (ppure . update (head zs2) α2 f2 ))) + + (1 - (one / lend zs1)) * (dist (kant d) (bind utail1 sgdRec1) (bind utail2 sgdRec2)) + === (one / lend zs1) * (dist (kant d) (fmap (update (head zs1) α1 f1) (sgd zs1 ws1 a1 f1)) + (fmap (update (head zs2) α2 f2 ) (sgd zs2 ws2 a2 f2))) + + (1 - (one / lend zs1)) * (dist (kant d) (bind utail1 sgdRec1) (bind utail2 sgdRec2)) + + ? fmapDist d d (2 * lip * α1) (update (head zs1) α1 f1) (sgd zs1 ws1 a1 f1) + (update (head zs2) α2 f2) (sgd zs2 ws2 a2 f2) + (relationalupdatep d (head zs1) α1 f1 (head zs2) α2 f2) + + =<= (one / lend zs1) * (dist (kant d) (sgd zs1 ws1 a1 f1) + (sgd zs2 ws2 a2 f2) + (2.0 * lip * α1)) + + (1 - (one / lend zs1)) * (dist (kant d) (bind utail1 sgdRec1) (bind utail2 sgdRec2)) + + ? thm d zs1 ws1 a1 f1 zs2 ws2 a2 f2 + ? assert (dist (kant d) (sgd zs1 ws1 a1 f1) (sgd zs2 ws2 a2 f2) <= dist d ws1 ws2 + estab zs1 a1) + =<= (one / lend zs1) * (dist d ws1 ws2 + estab zs1 a1 + (2.0 * lip * α1)) + + (1 - (one / lend zs1)) * (dist (kant d) (bind utail1 sgdRec1) (bind utail2 sgdRec2)) + ? bindDistEq d (dist d ws1 ws2 + estab zs1 a1) sgdRec1 utail1 sgdRec2 utail2 + (lemma d zs1 ws1 α1 a1 f1 zs2 ws2 α2 a2 f2) + ? assert (dist (kant d) (bind utail1 sgdRec1) (bind utail2 sgdRec2) <= dist d ws1 ws2 + estab zs1 a1) + ? assert (0 <= (1 - (one / lend zs1))) + ? multHelper ((one / lend zs1) * (dist d ws1 ws2 + estab zs1 a1 + (2.0 * lip * α1))) (1 - (one / lend zs1)) + (dist (kant d) (bind utail1 sgdRec1) (bind utail2 sgdRec2)) + (dist d ws1 ws2 + estab zs1 a1) + =<= (one / lend zs1) * (dist d ws1 ws2 + estab zs1 a1 + (2.0 * lip * α1)) + + (1 - (one / lend zs1)) * (dist d ws1 ws2 + estab zs1 a1) + + =<= (one / lend zs1) * (dist d ws1 ws2 + estab zs1 a1 + (2.0 * lip * α1)) + + (1 - (one / lend zs1)) * (dist d ws1 ws2 + estab zs1 a1) + + =<= dist d ws1 ws2 + 2.0 * lip * α1 * (one / lend zs1) + estab zs1 a1 + ? estabconsR zs1 α1 a1 + + =<= dist d ws1 ws2 + estab zs1 (SS α1 a1) + =<= dist d ws1 ws2 + estab zs1 as1 + *** QED + where + pureUpd1 = pureUpdate (head zs1) α1 f1 + pureUpd2 = pureUpdate (head zs2) α2 f2 + sgdRec1 = sgdRecUpd zs1 ws1 α1 a1 f1 + sgdRec2 = sgdRecUpd zs2 ws2 α2 a2 f2 + uhead1 = ppure (head zs1) + utail1 = unif (tail zs1) + uhead2 = ppure (head zs2) + utail2 = unif (tail zs2) +thm d zs1 ws1 _ f1 zs2 ws2 _ f2 = () + +{-@ multHelper :: a:Double -> b:{Double | 0 <= b} -> c:Double -> d:{Double | c <= d } + -> { a + b * c <= a + b * d } @-} +multHelper :: Double -> Double -> Double -> Double -> () +multHelper _ _ _ _ = () + + + +{-@ lemma :: d:Dist Double -> zs1:DataSet -> ws1:Weight -> α1:StepSize -> a1:StepSizes -> f1:LossFunction -> + zs2:{DataSet | lend zs1 == lend zs2 && tail zs1 = tail zs2} -> + ws2:Weight -> α2:{StepSize | α1 = α2} -> {a2:StepSizes| a2 = a1} -> f2:{LossFunction|f1 = f2} -> + z:DataPoint -> + {dist (kant d) (sgdRecUpd zs1 ws1 α1 a1 f1 z) (sgdRecUpd zs2 ws2 α2 a2 f2 z) <= dist d ws1 ws2 + estab zs1 a1} / [sslen a1, 1] @-} +lemma :: Dist Double -> DataSet -> Weight -> StepSize -> StepSizes -> LossFunction -> DataSet -> Weight -> StepSize -> StepSizes -> LossFunction -> DataPoint -> () +lemma d zs1 ws1 α1 a1 f1 zs2 ws2 α2 a2 f2 z = + dist (kant d) (sgdRecUpd zs1 ws1 α1 a1 f1 z) (sgdRecUpd zs2 ws2 α2 a2 f2 z) + === dist (kant d) (bind (sgd zs1 ws1 a1 f1) (pureUpdate z α1 f1)) + (bind (sgd zs2 ws2 a2 f2) (pureUpdate z α2 f2)) + ? pureUpdateEq z α1 f1 + ? pureUpdateEq z α2 f2 + === dist (kant d) (bind (sgd zs1 ws1 a1 f1) (ppure . update z α1 f1)) + (bind (sgd zs2 ws2 a2 f2) (ppure . update z α2 f2)) + === dist (kant d) (fmap (update z α1 f1) (sgd zs1 ws1 a1 f1)) + (fmap (update z α2 f2) (sgd zs2 ws2 a2 f2)) + ? fmapDist d d 0 (update z α1 f1) (sgd zs1 ws1 a1 f1) + (update z α2 f2) (sgd zs2 ws2 a2 f2) + (relationalupdateq d z α1 f1 z α2 f2) + =<= dist (kant d) (sgd zs1 ws1 a1 f1) (sgd zs2 ws2 a2 f2) + ? thm d zs1 ws1 a1 f1 zs2 ws2 a2 f2 + =<= dist d ws1 ws2 + estab zs1 a1 + *** QED + +{-@ assume pureUpdateEq :: zs:DataPoint -> a:StepSize -> f:LossFunction + -> {pureUpdate zs a f == ppure . update zs a f} @-} +pureUpdateEq :: DataPoint -> StepSize -> LossFunction -> () +pureUpdateEq zs a f = () diff --git a/src/TD/Lemmata/Relational/Act.hs b/src/TD/Lemmata/Relational/Act.hs new file mode 100644 --- /dev/null +++ b/src/TD/Lemmata/Relational/Act.hs @@ -0,0 +1,30 @@ +{-@ LIQUID "--reflection" @-} +{-@ LIQUID "--ple" @-} + +module TD.Lemmata.Relational.Act where + +import Monad.PrM +import Data.Dist +import Data.List +import Prelude hiding (max) + +import Monad.PrM.Predicates +import Monad.PrM.Relational.Theorems (mapMSpec) + +import TD.Lemmata.Relational.Sample + +import TD.TD0 +import Language.Haskell.Liquid.ProofCombinators + + + +{-@ relationalact :: l:Nat -> t:TransitionOf l -> m:{_|0 <= m} -> v1:{_ | llen v1 == l} + -> {v2:_|llen v2 = l} + -> {bounded m v1 v2 => lift (bounded (k * m)) (act l t v1) (act l t v2)} @-} +relationalact :: Int -> Transition -> Double -> ValueFunction -> ValueFunction -> () +relationalact _ t m v1 v2 | bounded m v1 v2 + = mapMSpec (k * m) + (sample v1 t) (sample v2 t) + (range 0 (llen v1)) + (relationalsample m (llen v1) t v1 v2) +relationalact _ _ _ _ _ = () diff --git a/src/TD/Lemmata/Relational/Iterate.hs b/src/TD/Lemmata/Relational/Iterate.hs new file mode 100644 --- /dev/null +++ b/src/TD/Lemmata/Relational/Iterate.hs @@ -0,0 +1,39 @@ +{-@ LIQUID "--reflection" @-} +{-@ LIQUID "--ple" @-} + +module TD.Lemmata.Relational.Iterate where + +import Monad.PrM +import Data.Dist +import Data.List +import Prelude hiding (iterate) + +import Monad.PrM.Relational.TCB.Spec +import Monad.PrM.Predicates + +import TD.TD0 +import Language.Haskell.Liquid.ProofCombinators +import Misc.ProofCombinators + + +{-@ relationaliterate :: m:{_|0 <= m} -> {k:_|k >= 0} -> n:Nat -> l:Nat + -> f:(v:ListN l -> PrM (ListN l)) + -> (m:{_|0 <= m} -> y1:{List Double|llen y1 = l} -> y2:{List Double|llen y2 = l} -> {bounded m y1 y2 => lift (bounded (k * m)) (f y1) (f y2)}) + -> x1:ListN l -> x2:ListN l + -> {bounded m x1 x2 => lift (bounded (pow k n * m)) ((iterate n (llen x1) f) (x1)) + ((iterate n (llen x2) f) (x2))} / [n] @-} +relationaliterate :: Double -> Double -> Int -> Int + -> (List Double -> PrM (List Double)) + -> (Double -> List Double -> List Double -> ()) + -> List Double -> List Double + -> () +relationaliterate m k 0 _ _ _ x1 x2 | bounded m x1 x2 + = pureSpec (bounded (pow k 0 * m)) x1 x2 () +relationaliterate m k n l f lemma x1 x2 | bounded m x1 x2 + = assert (pow k (n-1) * (k * m) == pow k n * m) ? + bindSpec (bounded (pow k n * m)) (bounded (k * m)) + (f x1) (iterate (n - 1) (llen x1) f) + (f x2) (iterate (n - 1) (llen x2) f) + (lemma m x1 x2) + (relationaliterate (k * m) k (n - 1) l f lemma) +relationaliterate m k n l f lemma x1 x2 = () diff --git a/src/TD/Lemmata/Relational/Sample.hs b/src/TD/Lemmata/Relational/Sample.hs new file mode 100644 --- /dev/null +++ b/src/TD/Lemmata/Relational/Sample.hs @@ -0,0 +1,87 @@ +{-@ LIQUID "--reflection" @-} +{-@ LIQUID "--ple" @-} + +module TD.Lemmata.Relational.Sample where + +import Monad.PrM +import Data.Dist +import Data.List +import Prelude hiding (max, uncurry) + +import Monad.PrM.Relational.TCB.Spec +import Monad.PrM.Predicates + +import TD.Lemmata.Relational.Update + +import TD.TD0 +import Language.Haskell.Liquid.ProofCombinators +import Misc.ProofCombinators + +{-@ listLemma :: v1:_ -> v2:SameLen v1 -> i:StateOf v1 + -> {distD (at v1 i) (at v2 i) <= distList distDouble v1 v2} @-} +listLemma :: ValueFunction -> ValueFunction -> State -> () +listLemma Nil v2 i = () +listLemma v1 Nil i = () +listLemma (Cons x xs) (Cons y ys) 0 = () +listLemma (Cons x xs) (Cons y ys) i = listLemma xs ys (i - 1) + +{-@ maxLemma :: v1:_ -> v2:SameLen v1 -> i:StateOf v1 -> j:StateOf v1 -> + {max (distD (at v1 i) (at v2 i)) (distD (at v1 j) (at v2 j)) <= distList distDouble v1 v2 } @-} +maxLemma :: ValueFunction -> ValueFunction -> State -> State -> () +maxLemma v1 v2 i j + = max (distD (at v1 i) (at v2 i)) (distD (at v1 j) (at v2 j)) + ? listLemma v1 v2 i + =<= max (distList distDouble v1 v2) (distD (at v1 j) (at v2 j)) + ? listLemma v1 v2 j + =<= max (distList distDouble v1 v2) (distList distDouble v1 v2) + =<= distList distDouble v1 v2 + *** QED + +{-@ updateLemma :: v1:_ -> v2:SameLen v1 -> i:StateOf v1 -> j:StateOf v1 -> r:_ -> + {distD (update v1 i j r) (update v2 i j r) <= k * distList distDouble v1 v2} @-} +updateLemma :: ValueFunction -> ValueFunction -> State -> State -> Reward -> () +updateLemma v1 v2 i j r + = distD (update v1 i j r) (update v2 i j r) + ? relationalupdate v1 v2 i j r + =<= k * max (distD (at v1 i) (at v2 i)) (distD (at v1 j) (at v2 j)) + ? maxLemma v1 v2 i j + =<= k * distList distDouble v1 v2 + *** QED + +{-@ uncurryLemma :: {m:_|0 <= m} -> v1:_ -> v2:SameLen v1 -> {_:_|bounded m v1 v2} -> i:StateOf v1 + -> t1:(StateOf v1, Reward) -> {t2:(StateOf v1, Reward)|t1 = t2} + -> {bounded' (k * m) (uncurry (update v1 i) t1) (uncurry (update v2 i) t2)} @-} +uncurryLemma :: Double -> ValueFunction -> ValueFunction -> () -> State + -> (State, Reward) -> (State, Reward) + -> () +uncurryLemma m v1 v2 b i t1@(j1, r1) t2@(j2, r2) + = distD (uncurry (update v1 i) t1) (uncurry (update v2 i) t2) + === distD (uncurry (update v1 i) t1) (uncurry (update v2 i) t1) + === distD (update v1 i j1 r1) (update v2 i j1 r1) + ? updateLemma v1 v2 i j1 r1 + =<= k * distList distDouble v1 v2 + ? b + =<= k * m + *** QED + +{-@ pureUpdateLemma :: {m:_|0 <= m} -> v1:_ -> v2:SameLen v1 -> {_:_|bounded m v1 v2} -> i:StateOf v1 + -> t1:(StateOf v1, Reward) -> {t2:(StateOf v1, Reward)|eqP t1 t2} + -> {lift (bounded' (k * m)) ((o ppure (uncurry (update v1 i))) (t1)) ((o ppure (uncurry (update v2 i))) (t2))} @-} +pureUpdateLemma :: Double -> ValueFunction -> ValueFunction -> () -> State -> (State, Reward) -> (State, Reward) -> () +pureUpdateLemma m v1 v2 b i t1@(j1, r1) t2@(j2, r2) = + pureSpec (bounded' (k * m)) + (uncurry (update v1 i) t1) (uncurry (update v2 i) t2) + (uncurryLemma m v1 v2 b i t1 t2) + + +{-@ relationalsample :: {m:_|0 <= m} -> n:Nat -> t:TransitionOf n -> {v1:_|llen t = llen v1 && llen v1 == n } -> {v2:_|llen t = llen v2} + -> i:StateOf t + -> {bounded m v1 v2 => lift (bounded' (k * m)) (sample v1 t i) (sample v2 t i)} @-} +relationalsample :: Double -> Int -> Transition -> ValueFunction -> ValueFunction -> State -> () +relationalsample m n t v1 v2 i | bounded m v1 v2 + = bindSpec (bounded' (k * m)) eqP + (t `at` i) (ppure `o` (uncurry (update v1 i))) + (t `at` i) (ppure `o` (uncurry (update v2 i))) + (liftSpec (t `at` i)) + (pureUpdateLemma m v1 v2 () i) +relationalsample _ _ _ _ _ _ = () diff --git a/src/TD/Lemmata/Relational/Update.hs b/src/TD/Lemmata/Relational/Update.hs new file mode 100644 --- /dev/null +++ b/src/TD/Lemmata/Relational/Update.hs @@ -0,0 +1,39 @@ +{-@ LIQUID "--reflection" @-} +{-@ LIQUID "--fast" @-} +{-@ LIQUID "--ple" @-} + +module TD.Lemmata.Relational.Update where + +import Monad.PrM +import Data.Dist +import Data.List +import Prelude hiding (max) + +import TD.TD0 +import Language.Haskell.Liquid.ProofCombinators + +{-@ relationalupdate :: v1:_ -> v2:SameLen v1 -> i:StateOf v1 -> j:StateOf v1 -> r:_ -> + {distD (update v1 i j r) (update v2 i j r) + <= k * max (distD (at v1 i) (at v2 i)) (distD (at v1 j) (at v2 j))} @-} +relationalupdate :: ValueFunction -> ValueFunction -> State -> State -> Reward -> () +relationalupdate v1 v2 i j r + = distD (update v1 i j r) (update v2 i j r) + === distD ((1 - α) * v1 `at` i + α * (r + γ * v1 `at` j)) + ((1 - α) * v2 `at` i + α * (r + γ * v2 `at` j)) + ? triangularIneq distDouble + ((1 - α) * v1 `at` i + α * (r + γ * v1 `at` j)) + ((1 - α) * v2 `at` i + α * (r + γ * v1 `at` j)) + ((1 - α) * v2 `at` i + α * (r + γ * v2 `at` j)) + =<= distD ((1 - α) * v1 `at` i + α * (r + γ * v1 `at` j)) + ((1 - α) * v2 `at` i + α * (r + γ * v1 `at` j)) + + distD ((1 - α) * v2 `at` i + α * (r + γ * v1 `at` j)) + ((1 - α) * v2 `at` i + α * (r + γ * v2 `at` j)) + ? linearity (1 - α) (α * (r + γ * v1 `at` j)) (v1 `at` i) (v2 `at` i) + =<= (1 - α) * distD (v1 `at` i) (v2 `at` i) + + distD ((1 - α) * v2 `at` i + α * (r + γ * v1 `at` j)) + ((1 - α) * v2 `at` i + α * (r + γ * v2 `at` j)) + ? linearity (α * γ) ((1 - α) * v2 `at` i + α * r) (v1 `at` j) (v2 `at` j) + =<= (1 - α) * distD (v1 `at` i) (v2 `at` i) + + α * γ * distD (v1 `at` j) (v2 `at` j) + =<= k * max (distD (v1 `at` i) (v2 `at` i)) (distD (v1 `at` j) (v2 `at` j)) + *** QED diff --git a/src/TD/TD0.hs b/src/TD/TD0.hs new file mode 100644 --- /dev/null +++ b/src/TD/TD0.hs @@ -0,0 +1,76 @@ +{-@ LIQUID "--reflection" @-} + +module TD.TD0 where + +-- import Monad.Implemented.PrM +import Monad.PrM +import Data.Dist +import Data.List + +import Prelude hiding ( map + , max + , repeat + , foldr + , fmap + , mapM + , iterate + , uncurry + ) + +{-@ type StateOf V = Idx V @-} +type State = Int +type Action = Int +type Reward = Double + +{-@ type TransitionOf N = {v:List (PrM ({i:State|0 <= i && i < N}, Reward))| llen v = N} @-} +type Transition = List (PrM (State, Reward)) +type ValueFunction = List Reward +type PrMValueFunction = PrM (List Reward) + +lq_required :: List Int -> () +lq_required _ = () + +{-@ reflect td0 @-} +{-@ td0 :: Nat -> v:ValueFunction -> TransitionOf (llen v) -> PrMValueFunction @-} +td0 :: Int -> ValueFunction -> Transition -> PrMValueFunction +td0 n v t = iterate n (llen v) (act (llen v) t) v + + + +{-@ reflect iterate @-} +{-@ iterate :: n:Nat -> l:Nat -> (v:{ValueFunction | llen v == l} -> PrM ({v':ValueFunction|llen v' = llen v})) -> + v:{ValueFunction | llen v == l} -> PrM ({v':ValueFunction|llen v' = llen v}) @-} +iterate :: Int -> Int -> (ValueFunction -> PrMValueFunction) -> ValueFunction -> PrMValueFunction +iterate n l _ x | n <= 0 = ppure x +iterate n l f x = bind (f x) (iterate (n - 1) l f) + + + +{-@ reflect act @-} +{-@ act :: n:Nat -> TransitionOf n -> v:{ValueFunction|llen v == n} + -> PrM {v':ValueFunction|llen v' = llen v} @-} +act :: Int -> Transition -> ValueFunction -> PrMValueFunction +act n t v = mapM (sample v t) (range 0 (llen v)) + +{-@ reflect uncurry @-} +uncurry :: (a -> b -> c) -> (a, b) -> c +uncurry f (a, b) = f a b + +{-@ reflect sample @-} +{-@ sample :: v:ValueFunction -> TransitionOf (llen v) -> StateOf v -> PrM Reward @-} +sample :: ValueFunction -> Transition -> State -> PrM Reward +sample v t i = bind (t `at` i) (ppure `o` (uncurry (update v i))) + +{-@ reflect γ @-} +{-@ reflect α @-} +{-@ reflect k @-} +γ, α, k :: Double +γ = 0.2 +α = 0.5 +k = 1 - α + α * γ + +{-@ reflect update @-} +{-@ update :: v:ValueFunction -> StateOf v -> StateOf v -> Reward -> Reward @-} +update :: ValueFunction -> State -> State -> Reward -> Reward +update v i j r = (1 - α) * (v `at` i) + α * (r + γ * v `at` j) + diff --git a/src/TD/Theorem.hs b/src/TD/Theorem.hs new file mode 100644 --- /dev/null +++ b/src/TD/Theorem.hs @@ -0,0 +1,25 @@ +{-@ LIQUID "--reflection" @-} +{-@ LIQUID "--ple" @-} + +module TD.Theorem where + +import Monad.PrM +import Data.Dist +import Data.List + +import Monad.PrM.Predicates + + +import TD.Lemmata.Relational.Act +import TD.Lemmata.Relational.Iterate + +import TD.TD0 +import Language.Haskell.Liquid.ProofCombinators +import Misc.ProofCombinators + + +{-@ relationaltd0 :: n:Nat -> l:Nat -> t:TransitionOf l -> {v1:_|llen v1 = l} -> v2:SameLen v1 -> + {lift (bounded (pow k n * (distList distDouble v1 v2))) (td0 n v1 t) (td0 n v2 t)} @-} +relationaltd0 :: Int -> Int -> Transition -> ValueFunction -> ValueFunction -> () +relationaltd0 n l t v1 v2 + = relationaliterate (distList distDouble v1 v2) k n l (act l t) (relationalact l t) v1 v2 diff --git a/test/Spec.hs b/test/Spec.hs new file mode 100644 --- /dev/null +++ b/test/Spec.hs @@ -0,0 +1,1 @@ +{-# OPTIONS_GHC -F -pgmF tasty-discover -optF --tree-display #-} diff --git a/test/Spec/Bins.hs b/test/Spec/Bins.hs new file mode 100644 --- /dev/null +++ b/test/Spec/Bins.hs @@ -0,0 +1,82 @@ +module Spec.Bins where + +import Test.HUnit ( assertEqual + , (@?) + , (@?=) + , Assertion + ) +import Data.Sort ( sort ) +import Numeric.Probability.Distribution + hiding ( map ) +import Spec.Utils + +import Monad.PrM hiding ( fmap ) +import Bins.Bins + +p, q :: Double +p = 0.5 +q = 0.625 + +coupling :: Double -> Double -> PrM (Bool, Bool) +coupling p q = + fromFreqs [((True, True), p), ((False, True), q - p), ((False, False), 1 - q)] + +mockbins :: PrM (Bool, Bool) -> Int -> PrM (Int, Int) +mockbins _ 0 = return (0, 0) +mockbins c n = do + (xl, xr) <- c + (yl, yr) <- mockbins c (n - 1) + return (yl + toInt xl, yr + toInt xr) + where toInt x = if x then 1 else 0 + +binsIter1 = sort [ ((1, 1), p) + , ((0, 1), q - p) + , ((0, 0), 1 - q) + ] + +binsIter2 = sort [ ((2, 2), p ^ 2) + , ((1, 2), 2 * p * (q - p)) + , ((1, 1), 2 * p * (1 - q)) + , ((0, 2), (q - p) ^ 2) + , ((0, 1), 2 * (q - p) * (1 - q)) + , ((0, 0), (1 - q) ^ 2) + ] + +unit_mockbins_1_it :: Assertion +unit_mockbins_1_it = + bins @?= binsIter1 + where + bins = clean $ decons $ mockbins (coupling p q) 1 + +unit_mockbins_2_it :: Assertion +unit_mockbins_2_it = + bins @?= binsIter2 + where + bins = clean $ decons $ mockbins (coupling p q) 2 + +unit_bins_1_it :: Assertion +unit_bins_1_it = do + resl @?= clean (map (\((a, _), p) -> (fromIntegral a, p)) binsIter1) + resr @?= clean (map (\((_, b), p) -> (fromIntegral b, p)) binsIter1) + where + resl = clean $ decons $ bins p 1 + resr = clean $ decons $ bins q 1 + +unit_bins_2_it :: Assertion +unit_bins_2_it = do + resl @?= clean (map (\((a, _), p) -> (fromIntegral a, p)) binsIter2) + resr @?= clean (map (\((_, b), p) -> (fromIntegral b, p)) binsIter2) + where + resl = clean $ decons $ bins p 2 + resr = clean $ decons $ bins q 2 + +unit_exp_dist_mockbins :: Assertion +unit_exp_dist_mockbins = + expDist == fromIntegral n * (q - p) + @? "want: E[dist (bins p n) (bins q n)] <= n * (q - p), got: " ++ show expDist + where + n = 10 + bins = mockbins (coupling p q) n + dist = fmap (\(a, b) -> fromIntegral (b - a)) bins + expDist = expected dist + diff --git a/test/Spec/SGD.hs b/test/Spec/SGD.hs new file mode 100644 --- /dev/null +++ b/test/Spec/SGD.hs @@ -0,0 +1,26 @@ +module Spec.SGD where + +import Test.HUnit ( assertEqual + , (@?) + , (@?=) + , Assertion + ) +import Numeric.Probability.Distribution + ( decons ) +import Spec.Utils + +import SGD.SGD + +{-@ loss :: DataPoint -> {ws:[Weight]|len ws = 1} -> Dbl @-} +loss :: DataPoint -> Weight -> Double +loss (x, y) w = (y - x + w) ^ 2 + +dp :: DataPoint +dp = (0, 1) + +ss :: StepSizes +ss = SS 0.5 (SS 0.5 (SS 0.5 (SS 0.5 SSEmp))) + +unit_sgd :: Assertion +unit_sgd = w @?= (-1) + where [(w, 1)] = clean $ decons $ sgd (replicate 4 dp) 1 ss loss diff --git a/test/Spec/TD0.hs b/test/Spec/TD0.hs new file mode 100644 --- /dev/null +++ b/test/Spec/TD0.hs @@ -0,0 +1,46 @@ +{-# LANGUAGE PackageImports #-} +{-@ LIQUID "--reflection" @-} + +module Spec.TD0 where + +import Test.HUnit ( assertEqual + , (@?) + , (@?=) + , Assertion + ) +import TD.TD0 ( td0 + , ValueFunction + , Transition + ) + +import Monad.PrM +import Data.Dist +import "safe-coupling" Data.List + +import Prelude hiding ( map + , repeat + , foldr + , fmap + , mapM + ) + +import Numeric.Probability.Distribution + ( decons ) + +v0 :: ValueFunction +v0 = Cons 1.0 (Cons (-1.0) Nil) + +t :: Transition +t = Cons (ppure (0, 0)) (Cons (ppure (0, 0)) Nil) + +unit_td0_base :: Assertion +unit_td0_base = + v @?= v0 + where [(v, 1)] = decons $ td0 0 v0 t + +unit_td0_simple :: Assertion +unit_td0_simple = + v @?= Cons 0.36 (Cons (-0.14) Nil) + where + [(v, 1)] = decons $ td0 2 v0 t + diff --git a/test/Spec/Utils.hs b/test/Spec/Utils.hs new file mode 100644 --- /dev/null +++ b/test/Spec/Utils.hs @@ -0,0 +1,10 @@ +module Spec.Utils where + +import Data.Sort ( sort ) + +clean :: (Ord a, Eq a) => [(a, Double)] -> [(a, Double)] +clean = foldr append [] . sort + +append :: Eq a => (a, Double) -> [(a, Double)] -> [(a, Double)] +append (v', p') ((v, p):xs) | v == v' = (v, p + p'):xs +append x xs = x:xs