backprop 0.0.3.0 → 0.1.0.0

raw patch · 25 files changed

+3539/−6156 lines, 25 filesdep +lensdep +primitivedep −addep −generics-sopdep −microlens-mtldep ~basebinary-addedPVP ok

version bump matches the API change (PVP)

Dependencies added: lens, primitive

Dependencies removed: ad, generics-sop, microlens-mtl, microlens-th, mtl, profunctors, tagged, transformers-base

Dependency ranges changed: base

API changes (from Hackage documentation)

Files

• Build.hs +17/−12
• CHANGELOG.md +29/−1
• LICENSE +1/−1
• README.md +209/−133
• backprop.cabal +78/−97
• bench/MNISTBench.hs +126/−106
• renders/backprop-mnist.md +353/−214
• renders/backprop-mnist.pdf binary
• renders/backprop-neural-test.md +0/−447
• renders/backprop-neural-test.pdf binary
• renders/extensible-neural.md +552/−0
• renders/extensible-neural.pdf binary
• samples/backprop-mnist.lhs +321/−195
• samples/backprop-monotest.hs +0/−18
• samples/backprop-neural-test.lhs +0/−405
• samples/extensible-neural.lhs +530/−0
• src/Data/Type/Util.hs +50/−146
• src/Numeric/Backprop.hs +394/−1480
• src/Numeric/Backprop/Implicit.hs +0/−393
• src/Numeric/Backprop/Internal.hs +596/−232
• src/Numeric/Backprop/Iso.hs +0/−209
• src/Numeric/Backprop/Mono.hs +0/−825
• src/Numeric/Backprop/Mono/Implicit.hs +0/−155
• src/Numeric/Backprop/Op.hs +283/−434
• src/Numeric/Backprop/Op/Mono.hs +0/−653

@@ -8,7 +8,7 @@ opts = shakeOptions { shakeFiles     = ".shake"                     , shakeVersion   = "1.0"                     , shakeVerbosity = Normal-                    , shakeThreads   = 0+                    , shakeThreads   = 1                     }  data Doc = Lab@@ -21,7 +21,7 @@     want ["all"]      "all" ~>-      need ["pdf", "md", "haddocks", "gentags", "install", "exe"]+      need ["pdf", "md", "gentags", "install", "exe"]      "pdf" ~>       need [ "renders" </> takeFileName f -<.> ".pdf"@@ -56,7 +56,7 @@       cmd "pandoc" "-V geometry:margin=1in"                    "-V fontfamily:palatino,cmtt"                    "-V links-as-notes"-                   "-sS"+                   "-s"                    "--highlight-style tango"                    "--reference-links"                    "--reference-location block"@@ -70,15 +70,20 @@         createDirectoryIfMissing True "samples-exe"         createDirectoryIfMissing True ".build"       removeFilesAfter "samples" ["/*.o"]-      cmd "stack ghc --" ("samples" </> src)-                         "-o" f-                         "-hidir" ".build"-                         "-threaded"-                         "-rtsopts"-                         "-with-rtsopts=-N"-                         "-Wall"-                         "-O2"-                         "-package backprop"+      cmd "stack" "ghc"+                  "--resolver lts-10"+                  "--package backprop"+                  "--package hmatrix"+                  "--package lens"+                  "--package mnist-idx"+                  "--package one-liner-instances"+                  "--package split"+                  "--"+                  ("samples" </> src)+                  "-o" f+                  "-hidir" ".build"+                  "-Wall"+                  "-O2"      ["tags","TAGS"] &%> \_ -> do       need (("src" </>) <$> allSrc)

@@ -1,9 +1,33 @@ Changelog ========= +Version 0.1.0.0+---------------++*Feb 5, 2018*++<https://github.com/mstksg/backprop/releases/tag/v0.1.0.0>++*   First non-alpha release.+*   More or less complete redesign of library.  The entire API is completely+    changed, and there is no backwards compatibility!+    *   Everything is now "implicit" style, and there is no more `BP` monad.+    *   Accessing items in `BVar`s is now lens-, prism-, and traversal- based,+        instead of iso- and generics-based.+    *   `Op` is no longer monadic+    *   *Mono* modules are removed.+    *   *Implicit* modules are removed, since they are the default+    *   *Iso* module is removed, since `Iso`s no longer play major role in the+        implementation of the library.+*   Removed dependency on *ad* and *ad*-based ops, which had been pulling in+    the vast majority of dependencies.+*   Moved from *.cabal* file to *hpack* system.+ Version 0.0.3.0 --------------- +*Alpha*+ <https://github.com/mstksg/backprop/releases/tag/v0.0.3.0>  *   Removed samples as registered executables in the cabal file, to reduce@@ -13,11 +37,13 @@ *   Added experimental (unsafe) combinators for working with GADTs with     existential types, `withGADT`, to *Numeric.Backprop* module. -*   Fixed broken links in Changelog.+*   Fixed broken links in changelog.  Version 0.0.2.0 --------------- +*Alpha*+ <https://github.com/mstksg/backprop/releases/tag/v0.0.2.0>  *   Added optimized numeric `Op`s, and re-write `Num`/`Fractional`/`Floating`@@ -32,6 +58,8 @@  Version 0.0.1.0 ---------------++*Alpha*  <https://github.com/mstksg/backprop/releases/tag/v0.0.1.0> 

@@ -1,4 +1,4 @@-Copyright Justin Le (c) 2017+Copyright Justin Le (c) 2018  All rights reserved. 

@@ -7,14 +7,15 @@ [**Literate Haskell Tutorial/Demo on MNIST data set**][mnist-lhs] (and [PDF rendering][mnist-pdf]) -Automatic *heterogeneous* back-propagation that can be used either *implicitly*-(in the style of the [ad][] library) or using *explicit* graphs built in-monadic style.  Implements reverse-mode automatic differentiation.  Differs-from [ad][] by offering full heterogeneity -- each intermediate step and the-resulting value can have different types.  Mostly intended for usage with-tensor manipulation libraries to implement automatic back-propagation for-gradient descent and other optimization techniques.+Automatic *heterogeneous* back-propagation. +Write your functions to compute your result, and the library will automatically+generate functions to compute your gradient.++Differs from [ad][] by offering full heterogeneity -- each intermediate step+and the resulting value can have different types.  Mostly intended for usage+with gradient descent and other numeric optimization techniques.+ [ad]: http://hackage.haskell.org/package/ad  Currently up on [hackage][] (with 100% documentation coverage), but more@@ -23,15 +24,6 @@ [hackage]: http://hackage.haskell.org/package/backprop [docs]: https://mstksg.github.io/backprop -At the moment this project is in **pre-alpha** (*v0.0.1.0*), and is-published/put up on Hackage as a call for comments and thoughts.  It has 100%-documentation coverage at the moment.  Performance was not yet a priority-before this, but will be from now on.  (Previously, highest priority was-API/usability). See [the todos section][todos] for more information on what's-missing, and how one would be able to contribute!--[todos]: https://github.com/mstksg/backprop#todo- MNIST Digit Classifier Example ------------------------------ @@ -50,10 +42,19 @@  [stack]: http://haskellstack.org/ -~~~bash+```bash $ ./Build.hs exe-~~~+``` +After the MNIST tutorial, there is a follow-up tutorial on using the library+with more advanced types, with extensible neural networks a la [this blog+post][blog], [available as literate haskell][neural-lhs] and also [rendered as+a PDF][neural-pdf].++[blog]: https://blog.jle.im/entries/series/+practical-dependent-types-in-haskell.html+[neural-lhs]: https://github.com/mstksg/backprop/blob/master/samples/extensible-neural.lhs+[neural-pdf]: https://github.com/mstksg/backprop/blob/master/renders/extensible-neural.pdf+ Brief example ------------- @@ -62,163 +63,238 @@ which is parameterized by two weight matrices and two bias vectors. Vector/matrix types are from the *hmatrix* package. -~~~haskell+Let's make a data type to store our parameters, with convenient accessors using+*[lens][]*:++[lens]: http://hackage.haskell.org/package/lens++```haskell+data Network i h o = Net { _weight1 :: L h i+                         , _bias1   :: R h+                         , _weight2 :: L o h+                         , _bias2   :: R o+                         }++makeLenses ''Network+```++Normally, we might write code to "run" a neural network on an input like this:++```haskell+neuralNet+    :: R i+    -> Network i h o+    -> R h+neuralNet x n = z+  where+    y = logistic $ (n ^. weight1) #> x + (n ^. bias1)+    z = logistic $ (n ^. weight2) #> y + (n ^. bias2)+ logistic :: Floating a => a -> a logistic x = 1 / (1 + exp (-x))+``` -matVec-    :: (KnownNat m, KnownNat n)-    => Op '[ L m n, R n ] (R m)+(`R i` is an i-length vector, `L h i` is an h-by-i matrix, etc., `#>` is+matrix-vector multiplication, and `^.` is access to a field via lens.) -neuralNetImplicit-      :: (KnownNat m, KnownNat n, KnownNat o)-      => R m-      -> BPOpI s '[ L n m, R n, L o n, R o ] (R o)-neuralNetImplicit inp = \(w1 :< b1 :< w2 :< b2 :< Ø) ->-    let z = logistic (liftB2 matVec w1 x + b1)-    in  logistic (liftB2 matVec w2 z + b2)-  where-    x = constRef inp+When given an input vector and the network, we compute the result of the neural+network ran on the input vector. -neuralNetExplicit-      :: (KnownNat m, KnownNat n, KnownNat o)-      => R m-      -> BPOp s '[ L n m, R n, L o n, R o ] (R o)-neuralNetExplicit inp = withInps $ \(w1 :< b1 :< w2 :< b2 :< Ø) -> do-    y1  <- matVec ~$ (w1 :< x1 :< Ø)-    let x2 = logistic (y1 + b1)-    y2  <- matVec ~$ (w2 :< x2 :< Ø)-    return $ logistic (y2 + b2)+We can write it, instead, using *backprop*:++```haskell+neuralNet+    :: Reifies s W+    => BVar s (R i)+    -> BVar s (Network i h o)+    -> BVar s (R o)+neuralNet x n = z   where-    x1 = constVar inp-~~~+    y = logistic $ (n ^^. weight1) #>! x + (n ^^. bias1)+    z = logistic $ (n ^^. weight2) #>! y + (n ^^. bias2) -Now `neuralNetExplicit` and `neuralNetImplicit` can be "run" with the input-vectors and parameters (a `L n m`, `R n`, `L o n`, and `R o`) and calculate the-output of the neural net.+logistic :: Floating a => a -> a+logistic x = 1 / (1 + exp (-x))+``` -~~~haskell-runNet-    :: (KnownNat m, KnownNat n, KnownNat o)-    => R m-    -> Tuple '[ L n m, R n, L o n, R o ]-    -> R o-runNet inp = evalBPOp (neuralNetExplicit inp)-~~~+(`#>!` is a backprop-aware version of `#>`, and `^^.` is access to a field via+lens in a `BVar`) -But, in defining `neuralNet`, we also generated a graph that *backprop* can-use to do back-propagation, too!+And that's it!  `neuralNet` is now backpropagatable! -~~~haskell-dot :: KnownNat n-    => Op '[ R n  , R n ] Double+We can "run" it using `evalBP`: -netGrad-    :: forall m n o. (KnownNat m, KnownNat n, KnownNat o)-    => R m+```haskell+evalBP (neuralNet (constVar x)) :: Network i h o -> R o+```++And we can find the gradient using `gradBP`:++```haskell+gradBP (neuralNet (constVar x)) :: Network i h o -> Network i h o+```++If we write a function to compute errors:++```haskell+netError+    :: Reifies s W+    => BVar s (R i)+    -> BVar s (R o)+    -> BVar s (Network i h o)+    -> BVar s Double+netError x targ n = sum' (err <.>! err)+  where+    err = neuralNet x - t+```++(`sum'` is a backprop-aware vector sum, and `<.>!` is a backprop-aware dot+product)++Now, we can perform gradient descent!++```haskell+gradDescent+    :: R i     -> R o-    -> Tuple '[ L n m, R n, L o n, R o ]-    -> Tuple '[ L n m, R n, L o n, R o ]-netGrad inp targ params = gradBPOp opError params+    -> Network i h o+    -> Network i h o+gradDescent x targ n0 = n0 - 0.1 * gradient   where-    -- calculate squared error, in *explicit* style-    opError :: BPOp s '[ L n m, R n, L o n, R o ] Double-    opError = do-        res <- neuralNetExplicit inp-        err <- bindRef (res - t)-        dot ~$ (err :< err :< Ø)-      where-        t = constRef targ-~~~+    gradient = gradBP (netError (constVar x) (constVar targ)) n0+``` -The result is the gradient of the input tuple's components, with respect-to the `Double` result of `opError` (the squared error).  We can then use-this gradient to do gradient descent.+Ta dah!  We were able to compute the gradient of our error function, just by+only saying how to compute *the error itself*.  For a more fleshed out example, see the [MNIST tutorial][mnist-lhs] (also [rendered as a pdf][mnist-pdf]) +Lens Access+-----------++A lot of the friction of dealing with `BVar s a`s instead of `a`s directly is+alleviated with the lens interface.++With a lens, you can "view" and "set" items inside a `BVar`, as if they were+the actual values:++```haskell+(^.)  ::        a -> Lens' a b ->        b+(^^.) :: BVar s a -> Lens' a b -> BVar s b++(.~)  :: Lens' a b ->        b ->        a ->        a+(.~~) :: Lens' a b -> BVar s b -> BVar s a -> BVar s a+```++And you can also extract multiple potential targets, as well, using+`Traversal`s and `Prism`s:++```haskell+-- | Actually takes a Traversal, to be more general.+-- Can be used to implement "pattern matching" on BVars+(^?)  ::        a -> Prism' a b -> Maybe (       b)+(^^?) :: BVar s a -> Prism' a b -> Maybe (BVar s b)++(^..)  ::        a -> Traversal' a b -> [       b]+(^^..) :: BVar s a -> Traversal' a b -> [BVar s b]+```++Note that the library itself has no *lens* dependency, using *[microlens][]*+instead.++[microlens]: http://hackage.haskell.org/package/microlens+ Benchmarks ---------- -The current version isn't optimized, but here are some basic benchmarks-comparing the library's automatic differentiation process to "manual"-differentiation by hand.  When using the [MNIST tutorial][bench] as an-example:+Here are some basic benchmarks comparing the library's automatic+differentiation process to "manual" differentiation by hand.  When using the+[MNIST tutorial][bench] as an example:  [bench]: https://github.com/mstksg/backprop/blob/master/bench/MNISTBench.hs -![benchmarks](http://i.imgur.com/xIZbhHa.png)+![benchmarks](https://i.imgur.com/9DXUaOI.png) -Calculating the gradient using *backprop* and calculating it by hand (by manual-symbolic differentiation) are within an order of magnitude of each-other,-time-wise.  Using the *backprop* library takes about *6.5x* as long-in this case.+*   For computing the gradient, there is about a 2.5ms overhead (or about 3.5x)+    compared to computing the gradients by hand.  Some more profiling and+    investigation can be done, since there are two main sources of potential+    slow-downs: -However, a full *update* step (calculate the gradient and update the neural-net) adds a lot of constant costs, so for a full training step, the *backprop*-library takes only *2.7x* as long as manual symbolic differentation.+    1.  "Inefficient" gradient computations, because of automated+        differentiation not being as efficient as what you might get from doing+        things by hand and simplifying.  This sort of cost is probably not+        avoidable.+    2.  Overhead incurred by the book-keeping and actual automatic+        differentiating system, which involves keeping track of a dependency+        graph and propagating gradients backwards in memory.  This sort of+        overhead is what we would be aiming to reduce. -This means using this library only slows down your program by a factor of-about 2.5x, compared to using only *hmatrix*.+    It is unclear which one dominates the current slowdown. -It's still definitely not ideal that more than half of the computation time is-overhead from the library, but this is just where we stand at the moment.-Optimization is just now starting!+*   However, it may be worth noting that this isn't necessarily a significant+    bottleneck.  *Updating* the networks using *hmatrix* actually dominates the+    runtime of the training.  Manual gradient descent takes 3.2ms, so the extra+    overhead is about 60%-70%. -Note that at the moment, simply running the network is only slightly slower-when using *backprop*.+*   Running the network (and the backprop-aware functions) incurs virtually+    zero overhead (about 4%), meaning that library authors could actually+    export backprop-aware functions by default and not lose any performance.  Todo ---- -1.  Profiling, to gauge where the overhead comes from (compared to "manual"-    back-propagation) and how to bring it down.+1.  Benchmark against competing back-propagation libraries like *ad*, and+    auto-differentiating tensor libraries like *[grenade][]* -2.  Some simple performance and API tweaks that are probably possible now and-    would clearly benefit: (if you want to contribute)+    [grenade]: https://github.com/HuwCampbell/grenade -    a.  ~~Providing optimized `Num`/`Fractional`/`Floating` instances for `BVal`-        by supplying known gradients directly instead of relying on *ad*.~~-        (Now finished, since [b3898ae][optnum])+2.  Write tests! -[optnum]: https://github.com/mstksg/backprop/commit/b3898ae676b8048e03709fb5d3d38a6fedb48e1e+3.  Explore potentially ditching `Num` for another typeclass that only has `+`,+    `0`, and `1`.  Currently, `Num` is required for all backpropagated types,+    but only `+`, `fromInteger 0`, and `fromInteger 1` are ever used. -    b.  Switch from `ST s` to `IO`, and use `unsafePerformIO` to automatically-        bind `BVal`s (like *ad* does) when using `liftB`.  This might remove-        some overhead during graph building, and, from an API standpoint,-        remove the need for explicit binding.+    The main upside to using `Num` is that it integrates well with the rest of+    the Haskell ecosystem, and many things already have useful `Num` instances. -    c.  Switch from `STRef`s/`IORef`s to `Array`.  (This one I'm unclear if it-        would help any)+    There are two downsides -- one minor and one major. -3.  Benchmark against competing back-propagation libraries like *ad*, and-    auto-differentiating tensor libraries like *[grenade][]*+    *   It requires more work to make a type backpropagatable.  Instead of+        writing only `+`, `0` and `1`, users must also define `*`, `-` or+        `negate`, `abs`, `signum`, and all of `fromInteger`.  However, I don't+        see this being a big issue in practice, since most values that will be+        used with *backprop* would presumably also benefit from having a full+        `Num` instance even without the need to backprop. -[grenade]: https://github.com/HuwCampbell/grenade+    *   Automatically generated prisms (used with `^^?`) work with tuples, and+        so cannot work out-of-the-box without a `Num` instance for tuples.  In+        addition, it's often useful to have anonymous products and tuples in+        general. -4.  Explore opportunities for parallelization.  There are some naive ways of-    directly parallelizing right now, but potential overhead should be-    investigated.+        However, this can be resolved by using the orphan instances in the+        *[NumInstances][]* package.  Still, there might be some headache for+        application developers if different libraries using *backprop*+        accidentally pull in their orphan instances from different places. -5.  Some open questions:+        [NumInstances]: https://hackage.haskell.org/package/NumInstances -    a. Is it possible to offer pattern matching on sum types/with different-       constructors for implicit-graph backprop?  It's possible for-       explicit-graph versions already, with `choicesVar`, but not yet with-       the implicit-graph interface.  Could be similar to an "Applicative vs.-       Monad" issue where you can only have pre-determined fixed computation-       paths when using `Applicative`, but I'm not sure.  Still, it would be-       nice, because if this was possible, we could possibly do away with-       explicit-graph mode completely.+    The extra complexity that would come from adding a custom typeclass just+    for `+` / `0` / `1`, though, I feel, might not be worth the benefit.  The+    entire numeric Haskell ecosystem, at the time, revolves around `Num`. -    b. Though we already have safe sum type support with explicit-graph mode,-       we can't support GADTs yet safely.  It'd be nice to see if this is-       possible, because a lot of dependently typed neural network stuff is-       made much simpler with GADTs.+    However, it is worth noting that it wouldn't be too hard to add "Additive+    Typeclass" instances for any custom types -- one would just need to define+    `(<+>) = (+)`, `zero = fromInteger 1`, and `one = fromInteger 1` (a+    three-liner), so it might not be too bad. -       As of v0.0.3.0, we have a way of dealing with GADTs in explicit-graph-       mode (using `withGADT`) that is *unsafe*, and requires some ugly manual-       plumbing by the user that could potentially be confusing.  But it would-       still be nice to have a way that is safe and doesn't require the manual-       plumbing and isn't as easy to mess up.+    But really, a lot of this would all resolve itself if we got `Num`+    instances for tuples in base :)++3.  Explore opportunities for parallelization.  There are some naive ways of+    directly parallelizing right now, but potential overhead should be+    investigated.++4.  Some open questions:++    a. Is it possible to support constructors with existential types?

@@ -1,102 +1,83 @@-name:                backprop-version:             0.0.3.0-synopsis:            Heterogeneous, type-safe automatic backpropagation in Haskell-description:         See <https://github.com/mstksg/backprop#readme README.md>-                     .-                     At the moment, this project is in pre-alpha, and is-                     published and put up on Hackage with 100% documentation-                     coverage as a call for comments and thoughts.  See-                     <https://github.com/mstksg/backprop#todo TODO.md> section-                     in the README for more information on what's missing and-                     potential avenues for contribution.-homepage:            https://github.com/mstksg/backprop-bug-reports:         https://github.com/mstksg/backprop/issues-license:             BSD3-license-file:        LICENSE-author:              Justin Le-maintainer:          justin@jle.im-copyright:           (c) Justin Le 2017-category:            Web-build-type:          Simple-extra-source-files:  README.md-                     CHANGELOG.md-                     Build.hs-                     renders/backprop-mnist.md-                     renders/backprop-mnist.pdf-                     renders/backprop-neural-test.md-                     renders/backprop-neural-test.pdf-                     samples/backprop-mnist.lhs-                     samples/backprop-monotest.hs-                     samples/backprop-neural-test.lhs-cabal-version:       >=1.10--library-  hs-source-dirs:      src-  exposed-modules:     Numeric.Backprop-                       Numeric.Backprop.Implicit-                       Numeric.Backprop.Iso-                       Numeric.Backprop.Mono-                       Numeric.Backprop.Mono.Implicit-                       Numeric.Backprop.Op-                       Numeric.Backprop.Op.Mono-  other-modules:       Numeric.Backprop.Internal-                       Data.Type.Util-  build-depends:       base >= 4.7 && < 5-                     , ad-                     , generics-sop-                     , microlens-                     , microlens-mtl-                     , microlens-th-                     , mtl-                     , profunctors-                     , reflection-                     , tagged-                     , transformers-base-                     , type-combinators-  default-language:    Haskell2010-  ghc-options:         -Wall--benchmark backprop-mnist-bench-  type:                exitcode-stdio-1.0-  hs-source-dirs:      bench-  main-is:             MNISTBench.hs-  ghc-options:         -threaded -rtsopts -with-rtsopts=-N -Wall -O2-  build-depends:       base-                     , backprop-                     , bifunctors-                     , criterion-                     , deepseq-                     , directory-                     , generics-sop-                     , hmatrix    >= 0.18-                     , mnist-idx-                     , mwc-random-                     , time-                     , transformers-                     , type-combinators-                     , vector-  default-language:    Haskell2010+-- This file has been generated from package.yaml by hpack version 0.20.0.+--+-- see: https://github.com/sol/hpack+--+-- hash: 78d9facc552fa43f3b324fa4ffe1c526e33e67c5cf55fdde75b05f60f81e423c --- test-suite backprop-doctest---   type:                exitcode-stdio-1.0---   hs-source-dirs:      doctest---   main-is:             doctest.hs---   build-depends:       base---                      , backprop---                      , doctest---                      , Glob---   ghc-options:         -threaded -rtsopts -with-rtsopts=-N---   default-language:    Haskell2010+name:           backprop+version:        0.1.0.0+synopsis:       Heterogeneous automatic backpropagation in Haskell+description:    Write your functions to compute your result, and the library will+                automatically generate functions to compute your gradient.+                .+                See <https://github.com/mstksg/backprop#readme README.md>+category:       Math+homepage:       https://github.com/mstksg/backprop#readme+bug-reports:    https://github.com/mstksg/backprop/issues+author:         Justin Le+maintainer:     justin@jle.im+copyright:      (c) Justin Le 2018+license:        BSD3+license-file:   LICENSE+build-type:     Simple+cabal-version:  >= 1.10 --- test-suite backprop-test---   type:                exitcode-stdio-1.0---   hs-source-dirs:      test---   main-is:             Spec.hs---   build-depends:       base---                      , backprop---   ghc-options:         -threaded -rtsopts -with-rtsopts=-N---   default-language:    Haskell2010+extra-source-files:+    Build.hs+    CHANGELOG.md+    README.md+    renders/backprop-mnist.md+    renders/backprop-mnist.pdf+    renders/extensible-neural.md+    renders/extensible-neural.pdf+    samples/backprop-mnist.lhs+    samples/extensible-neural.lhs  source-repository head-  type:     git+  type: git   location: https://github.com/mstksg/backprop++library+  hs-source-dirs:+      src+  ghc-options: -Wall -fprint-explicit-kinds -fwarn-redundant-constraints+  build-depends:+      base >=4.7 && <5+    , deepseq+    , microlens+    , primitive+    , reflection+    , transformers+    , type-combinators+    , vector+  exposed-modules:+      Numeric.Backprop+      Numeric.Backprop.Op+  other-modules:+      Numeric.Backprop.Internal+      Data.Type.Util+  default-language: Haskell2010++benchmark backprop-mnist-bench+  type: exitcode-stdio-1.0+  main-is: MNISTBench.hs+  hs-source-dirs:+      bench+  ghc-options: -Wall -fprint-explicit-kinds -fwarn-redundant-constraints -threaded -rtsopts -with-rtsopts=-N -O2+  build-depends:+      backprop+    , base >=4.7 && <5+    , bifunctors+    , criterion+    , deepseq+    , directory+    , hmatrix >=0.18+    , lens+    , mnist-idx+    , mwc-random+    , time+    , transformers+    , vector+  other-modules:+      Paths_backprop+  default-language: Haskell2010

@@ -5,36 +5,32 @@ {-# LANGUAGE GADTs                #-} {-# LANGUAGE LambdaCase           #-} {-# LANGUAGE ScopedTypeVariables  #-}+{-# LANGUAGE TemplateHaskell      #-} {-# LANGUAGE TypeApplications     #-} {-# LANGUAGE ViewPatterns         #-} {-# OPTIONS_GHC -fno-warn-orphans #-}  import           Control.DeepSeq import           Control.Exception+import           Control.Lens hiding          ((:<), (<.>)) import           Control.Monad.IO.Class import           Control.Monad.Trans.Maybe import           Criterion.Main import           Criterion.Types import           Data.Bitraversable import           Data.IDX-import           Data.Maybe-import           Data.Time.Format-import           Data.Time.LocalTime+import           Data.Time import           Data.Traversable import           Data.Tuple-import           Data.Type.Combinator-import           Data.Type.Index-import           Data.Type.Product-import           GHC.Generics                        (Generic)+import           GHC.Generics                 (Generic) import           GHC.TypeLits import           Numeric.Backprop-import           Numeric.LinearAlgebra.Static hiding (dot)+import           Numeric.LinearAlgebra.Static import           System.Directory-import qualified Data.Vector.Generic                 as VG-import qualified Data.Vector.Unboxed                 as VU-import qualified Generics.SOP                        as SOP-import qualified Numeric.LinearAlgebra               as HM-import qualified System.Random.MWC                   as MWC+import qualified Data.Vector.Generic          as VG+import qualified Data.Vector.Unboxed          as VU+import qualified Numeric.LinearAlgebra        as HM+import qualified System.Random.MWC            as MWC  data Layer i o =     Layer { _lWeights :: !(L o i)@@ -42,9 +38,10 @@           }   deriving (Show, Generic) -instance SOP.Generic (Layer i o) instance NFData (Layer i o) +makeLenses ''Layer+ data Network i h1 h2 o =     Net { _nLayer1 :: !(Layer i  h1)         , _nLayer2 :: !(Layer h1 h2)@@ -52,82 +49,89 @@         }   deriving (Show, Generic) -instance SOP.Generic (Network i h1 h2 o) instance NFData (Network i h1 h2 o) -matVec-    :: (KnownNat m, KnownNat n)-    => Op '[ L m n, R n ] (R m)-matVec = op2' $ \m v ->-  ( m #> v, \(fromMaybe 1 -> g) ->-              (g `outer` v, tr m #> g)-  )+makeLenses ''Network -dot :: KnownNat n-    => Op '[ R n, R n ] Double-dot = op2' $ \x y ->-  ( x <.> y, \case Nothing -> (y, x)-                   Just g  -> (konst g * y, x * konst g)-  )+infixr 8 #>!+(#>!)+    :: (KnownNat m, KnownNat n, Reifies s W)+    => BVar s (L m n)+    -> BVar s (R n)+    -> BVar s (R m)+(#>!) = liftOp2 . op2 $ \m v ->+  ( m #> v, \g -> (g `outer` v, tr m #> g) ) -scale-    :: KnownNat n-    => Op '[ Double, R n ] (R n)-scale = op2' $ \a x ->-  ( konst a * x-  , \case Nothing -> (HM.sumElements (extract x      ), konst a    )-          Just g  -> (HM.sumElements (extract (x * g)), konst a * g)++infixr 8 <.>!+(<.>!)+    :: (KnownNat n, Reifies s W)+    => BVar s (R n)+    -> BVar s (R n)+    -> BVar s Double+(<.>!) = liftOp2 . op2 $ \x y ->+  ( x <.> y, \g -> (konst g * y, x * konst g)   ) -vsum-    :: KnownNat n-    => Op '[ R n ] Double-vsum = op1' $ \x -> (HM.sumElements (extract x), maybe 1 konst)+konst'+    :: (KnownNat n, Reifies s W)+    => BVar s Double+    -> BVar s (R n)+konst' = liftOp1 . op1 $ \c -> (konst c, HM.sumElements . extract) +sumElements :: KnownNat n => R n -> Double+sumElements = HM.sumElements . extract++sumElements'+    :: (KnownNat n, Reifies s W)+    => BVar s (R n)+    -> BVar s Double+sumElements' = liftOp1 . op1 $ \x -> (sumElements x, konst)+ logistic :: Floating a => a -> a logistic x = 1 / (1 + exp (-x))+{-# INLINE logistic #-}  runLayer-    :: (KnownNat i, KnownNat o)-    => BPOp s '[ R i, Layer i o ] (R o)-runLayer = withInps $ \(x :< l :< Ø) -> do-    w :< b :< Ø <- gTuple #<~ l-    y <- matVec ~$ (w :< x :< Ø)-    return $ y + b+    :: (KnownNat i, KnownNat o, Reifies s W)+    => BVar s (Layer i o)+    -> BVar s (R i)+    -> BVar s (R o)+runLayer l x = (l ^^. lWeights) #>! x + (l ^^. lBiases)+{-# INLINE runLayer #-} -runNetwork-    :: (KnownNat i, KnownNat h1, KnownNat h2, KnownNat o)-    => BPOp s '[ R i, Network i h1 h2 o ] (R o)-runNetwork = withInps $ \(x :< n :< Ø) -> do-    l1 :< l2 :< l3 :< Ø <- gTuple #<~ n-    y <- runLayer -$ (x          :< l1 :< Ø)-    z <- runLayer -$ (logistic y :< l2 :< Ø)-    r <- runLayer -$ (logistic z :< l3 :< Ø)-    softmax       -$ (r          :< Ø)+softMax :: (KnownNat n, Reifies s W) => BVar s (R n) -> BVar s (R n)+softMax x = konst' (1 / sumElements' expx) * expx+  where+    expx = exp x+{-# INLINE softMax #-} -softmax :: KnownNat n => BPOp s '[ R n ] (R n)-softmax = withInps $ \(x :< Ø) -> do-    expX <- bindVar (exp x)-    totX <- vsum ~$ (expX   :< Ø)-    scale        ~$ (1/totX :< expX :< Ø)+runNetwork+    :: (KnownNat i, KnownNat h1, KnownNat h2, KnownNat o, Reifies s W)+    => BVar s (Network i h1 h2 o)+    -> R i+    -> BVar s (R o)+runNetwork n = softMax+             . runLayer (n ^^. nLayer3)+             . logistic+             . runLayer (n ^^. nLayer2)+             . logistic+             . runLayer (n ^^. nLayer1)+             . constVar+{-# INLINE runNetwork #-} -crossEntropy-    :: KnownNat n-    => R n-    -> BPOpI s '[ R n ] Double-crossEntropy targ (r :< Ø) = negate (dot .$ (log r :< t :< Ø))-  where-    t = constVar targ+crossEntropy :: (KnownNat n, Reifies s W) => R n -> BVar s (R n) -> BVar s Double+crossEntropy t r = negate $ log r <.>! constVar t+{-# INLINE crossEntropy #-} -softMaxCrossEntropy-    :: KnownNat n-    => R n-    -> BPOpI s '[ R n ] Double-softMaxCrossEntropy targ (r :< Ø) =  realToFrac tsum * log (vsum .$ (r :< Ø))-                                       - (dot .$ (r :< t :< Ø))-  where-    tsum = HM.sumElements . extract $ targ-    t    = constVar targ+netErr+    :: (KnownNat i, KnownNat h1, KnownNat h2, KnownNat o, Reifies s W)+    => R i+    -> R o+    -> BVar s (Network i h1 h2 o)+    -> BVar s Double+netErr x t n = crossEntropy t (runNetwork n x)+{-# INLINE netErr #-}  trainStep     :: forall i h1 h2 o. (KnownNat i, KnownNat h1, KnownNat h2, KnownNat o)@@ -136,28 +140,35 @@     -> R o     -> Network i h1 h2 o     -> Network i h1 h2 o-trainStep r !x !t !n = case gradBPOp (netErr t) (x ::< n ::< Ø) of-    _ :< I gN :< Ø ->-        n - (realToFrac r * gN)+trainStep r !x !t !n = n - realToFrac r * gradBP (netErr x t) n+{-# INLINE trainStep #-} +runLayerManual+    :: (KnownNat i, KnownNat o)+    => Layer i o+    -> R i+    -> R o+runLayerManual l x = (l ^. lWeights) #> x + (l ^. lBiases)+{-# INLINE runLayerManual #-}++softMaxManual :: KnownNat n => R n -> R n+softMaxManual x = konst (1 / sumElements expx) * expx+  where+    expx = exp x+{-# INLINE softMaxManual #-}+ runNetManual-    :: forall i h1 h2 o. (KnownNat i, KnownNat h1, KnownNat h2, KnownNat o)+    :: (KnownNat i, KnownNat h1, KnownNat h2, KnownNat o)     => Network i h1 h2 o     -> R i     -> R o-runNetManual (Net (Layer w1 b1) (Layer w2 b2) (Layer w3 b3)) x =-    let y1 = w1 #> x-        z1 = y1 + b1-        x2 = logistic z1-        y2 = w2 #> x2-        z2 = y2 + b2-        x3 = logistic z2-        y3 = w3 #> x3-        z3 = y3 + b3-        o0 = exp z3-        o1 = HM.sumElements (extract o0)-        o2 = o0 / konst o1-    in  o2+runNetManual n = softMaxManual+               . runLayerManual (n ^. nLayer3)+               . logistic+               . runLayerManual (n ^. nLayer2)+               . logistic+               . runLayerManual (n ^. nLayer1)+{-# INLINE runNetManual #-}  gradNetManual     :: forall i h1 h2 o. (KnownNat i, KnownNat h1, KnownNat h2, KnownNat o)@@ -181,7 +192,7 @@         dEdO3 = 1         dEdO2 = dEdO3 * (- t / o2)         dEdO1 = - (dEdO2 <.> o0) / (o1 ** 2)-        dEdO0 = konst dEdO1 * 1 + dEdO2 / konst o1+        dEdO0 = konst dEdO1 + dEdO2 / konst o1         dEdZ3 = dEdO0 * o0         dEdY3 = dEdZ3         dEdX3 = tr w3 #> dEdY3@@ -190,7 +201,6 @@         dEdX2 = tr w2 #> dEdY2         dEdZ1 = dEdX2 * (x2 * (1 - x2))         dEdY1 = dEdZ1-         dEdB3 = dEdZ3         dEdW3 = dEdY3 `outer` x3         dEdB2 = dEdZ2@@ -198,6 +208,7 @@         dEdB1 = dEdZ1         dEdW1 = dEdY1 `outer` x     in  Net (Layer dEdW1 dEdB1) (Layer dEdW2 dEdB2) (Layer dEdW3 dEdB3)+{-# INLINE gradNetManual #-}  trainStepManual     :: forall i h1 h2 o. (KnownNat i, KnownNat h1, KnownNat h2, KnownNat o)@@ -210,14 +221,6 @@     let gN = gradNetManual x t n     in  n - (realToFrac r * gN) -netErr-    :: (KnownNat i, KnownNat h1, KnownNat h2, KnownNat o)-    => R o-    -> BPOp s '[ R i, Network i h1 h2 o ] Double-netErr t = do-    y <- runNetwork-    implicitly (crossEntropy t) -$ (y :< Ø)- main :: IO () main = MWC.withSystemRandom $ \g -> do     Just test  <- loadMNIST "data/t10k-images-idx3-ubyte"  "data/t10k-labels-idx1-ubyte"@@ -234,8 +237,7 @@         bgroup "gradient" [             let testManual x y = gradNetManual x y net0             in  bench "manual" $ nf (uncurry testManual) test0-          , let testBP     x y = getI . index (IS IZ) $-                  gradBPOp (netErr y) (x ::< net0 ::< Ø)+          , let testBP     x y = gradBP (netErr x y) net0             in  bench "bp"     $ nf (uncurry testBP) test0           ]       , bgroup "descent" [@@ -245,9 +247,9 @@             in  bench "bp"     $ nf (uncurry testBP) test0           ]       , bgroup "run" [-            let testManual x   = runNetManual net0 x+            let testManual     = runNetManual net0             in  bench "manual" $ nf testManual (fst test0)-          , let testBP     x   = evalBPOp runNetwork (x ::< net0 ::< Ø)+          , let testBP     x   = evalBP (`runNetwork` x) net0             in  bench "bp"     $ nf testBP (fst test0)           ]       ]@@ -311,4 +313,22 @@ instance (KnownNat i, KnownNat h1, KnownNat h2, KnownNat o) => MWC.Variate (Network i h1 h2 o) where     uniform g = Net <$> MWC.uniform g <*> MWC.uniform g <*> MWC.uniform g     uniformR (l, h) g = (\x -> x * (h - l) + l) <$> MWC.uniform g++instance (Num a, Num b) => Num (a, b) where+    (x1,y1) + (x2,y2) = (x1 + x2, y1 + y2)+    (x1,y1) * (x2,y2) = (x1 * x2, y1 * y2)+    (x1,y1) - (x2,y2) = (x1 - x2, y1 - y2)+    abs (x, y)        = (abs x, abs y)+    signum (x, y)     = (signum x, signum y)+    fromInteger x     = (fromInteger x, fromInteger x)++-- softMaxCrossEntropy+--     :: KnownNat n+--     => R n+--     -> BPOpI s '[ R n ] Double+-- softMaxCrossEntropy targ (r :< Ø) =  realToFrac tsum * log (vsum .$ (r :< Ø))+--                                        - (dot .$ (r :< t :< Ø))+--   where+--     tsum = HM.sumElements . extract $ targ+--     t    = constVar targ 

@@ -1,25 +1,23 @@ --- author: - Justin Le-fontfamily: 'palatino,cmtt'-geometry: margin=1in-links-as-notes: true title: Learning MNIST with Neural Networks with backprop library ---  The *backprop* library performs back-propagation over a *hetereogeneous*-system of relationships. It offers both an implicit (*[ad]*-like) and-explicit graph building usage style. Let’s use it to build neural-networks and learn mnist!+system of relationships. back-propagation is done automatically (as+reverse-mode automatic differentiation), and you work with your values+as if you were writing normal functions with them, with the help of+[lens]. -  [ad]: http://hackage.haskell.org/package/ad+  [lens]: http://hackage.haskell.org/package/lens  Repository source is [on github], and docs are [on hackage].    [on github]: https://github.com/mstksg/backprop   [on hackage]: http://hackage.haskell.org/package/backprop -If you’re reading this as a literate haskell file, you should know that+If you're reading this as a literate haskell file, you should know that a [rendered pdf version is available on github.]. If you are reading this as a pdf file, you should know that a [literate haskell version that you can run] is also available on github!@@ -27,22 +25,36 @@   [rendered pdf version is available on github.]: https://github.com/mstksg/backprop/blob/master/renders/backprop-mnist.pdf   [literate haskell version that you can run]: https://github.com/mstksg/backprop/blob/master/samples/backprop-mnist.lhs +The packages involved are:++-   deepseq+-   hmatrix+-   lens+-   mnist-idx+-   mwc-random+-   one-liner-instances+-   split+-   vector+ ``` {.sourceCode .literate .haskell} {-# LANGUAGE BangPatterns                     #-} {-# LANGUAGE DataKinds                        #-} {-# LANGUAGE DeriveGeneric                    #-}+{-# LANGUAGE FlexibleContexts                 #-} {-# LANGUAGE GADTs                            #-} {-# LANGUAGE LambdaCase                       #-} {-# LANGUAGE ScopedTypeVariables              #-}+{-# LANGUAGE TemplateHaskell                  #-} {-# LANGUAGE TupleSections                    #-} {-# LANGUAGE TypeApplications                 #-} {-# LANGUAGE ViewPatterns                     #-}-{-# OPTIONS_GHC -fno-warn-orphans             #-} {-# OPTIONS_GHC -fno-warn-incomplete-patterns #-}+{-# OPTIONS_GHC -fno-warn-orphans             #-} {-# OPTIONS_GHC -fno-warn-unused-top-binds    #-}  import           Control.DeepSeq import           Control.Exception+import           Control.Lens hiding ((<.>)) import           Control.Monad import           Control.Monad.IO.Class import           Control.Monad.Trans.Maybe@@ -51,34 +63,65 @@ import           Data.Foldable import           Data.IDX import           Data.List.Split-import           Data.Maybe import           Data.Time.Clock import           Data.Traversable import           Data.Tuple import           GHC.Generics                        (Generic) import           GHC.TypeLits import           Numeric.Backprop-import           Numeric.LinearAlgebra.Static hiding (dot)+import           Numeric.LinearAlgebra.Static+import           Numeric.OneLiner import           Text.Printf import qualified Data.Vector                         as V import qualified Data.Vector.Generic                 as VG import qualified Data.Vector.Unboxed                 as VU-import qualified Generics.SOP                        as SOP import qualified Numeric.LinearAlgebra               as HM import qualified System.Random.MWC                   as MWC import qualified System.Random.MWC.Distributions     as MWC ``` +Introduction+============++In this walkthrough, we'll be building a classifier for the *[MNIST]*+data set. This is meant to mirror the [Tensorflow Tutorial] for+beginners.++  [MNIST]: http://yann.lecun.com/exdb/mnist/+  [Tensorflow Tutorial]: https://www.tensorflow.org/versions/r1.2/get_started/mnist/beginners++Essentially, we use a two-layer artificial neural network -- or a series+of matrix multiplications, differentiable function applications, and+vector additions. We feed our input image to the ANN and then try to get+a label from it. Training an ANN is a matter of finding the right+matrices to multiply by, and the right vectors to add.++To do that, we train our network by treating our network's accuracy as a+function `Network -> Error`. If we can find the gradient of the input+network with respect to the error, we can perform [gradient descent],+and slowly make our network better and better.++  [gradient descent]: https://en.wikipedia.org/wiki/Gradient_descent++Finding the gradient is usually complicated, but *backprop* makes it+simpler:++1.  Write a function to compute the error from the network+2.  That's it!++Hooray! Once you do that, the library finds the gradient function+*automatically*, without any further intervention!+ Types ===== -For the most part, we’re going to be using the great *[hmatrix]* library+For the most part, we're going to be using the great *[hmatrix]* library and its vector and matrix types. It offers a type `L m n` for $m \times n$ matrices, and a type `R n` for an $n$ vector.    [hmatrix]: http://hackage.haskell.org/package/hmatrix -First things first: let’s define our neural networks as simple+First things first: let's define our neural networks as simple containers of parameters (weight matrices and bias vectors).  First, a type for layers:@@ -90,8 +133,8 @@           }   deriving (Show, Generic) -instance SOP.Generic (Layer i o) instance NFData (Layer i o)+makeLenses ''Layer ```  And a type for a simple feed-forward network with two hidden layers:@@ -104,59 +147,64 @@         }   deriving (Show, Generic) -instance SOP.Generic (Network i h1 h2 o) instance NFData (Network i h1 h2 o)+makeLenses ''Network ``` -These are pretty straightforward container types…pretty much exactly the-type you’d make to represent these networks! Note that, following true-Haskell form, we separate out logic from data. This should be all we-need.--We derive an instance of `SOP.Generic` from the *[generics-sop]*-package, which *backprop* uses to propagate derivatives on values inside-product types.--  [generics-sop]: http://hackage.haskell.org/package/generics-sop+These are pretty straightforward container types...pretty much exactly+the type you'd make to represent these networks! Note that, following+true Haskell form, we separate out logic from data. This should be all+we need.  Instances ---------  Things are much simplier if we had `Num` and `Fractional` instances for-everything, so let’s just go ahead and define that now, as well. Just a-little bit of boilerplate.+everything, so let's just go ahead and define that now, as well. Just a+little bit of boilerplate, made easier using *[one-liner-instances]* to+auto-derive instances using Generics. +  [one-liner-instances]: http://hackage.haskell.org/package/one-liner-instances+ ``` {.sourceCode .literate .haskell} instance (KnownNat i, KnownNat o) => Num (Layer i o) where-    Layer w1 b1 + Layer w2 b2 = Layer (w1 + w2) (b1 + b2)-    Layer w1 b1 - Layer w2 b2 = Layer (w1 - w2) (b1 - b2)-    Layer w1 b1 * Layer w2 b2 = Layer (w1 * w2) (b1 * b2)-    abs    (Layer w b)        = Layer (abs    w) (abs    b)-    signum (Layer w b)        = Layer (signum w) (signum b)-    negate (Layer w b)        = Layer (negate w) (negate b)-    fromInteger x             = Layer (fromInteger x) (fromInteger x)+    (+)         = gPlus+    (-)         = gMinus+    (*)         = gTimes+    negate      = gNegate+    abs         = gAbs+    signum      = gSignum+    fromInteger = gFromInteger -instance (KnownNat i, KnownNat h1, KnownNat h2, KnownNat o) => Num (Network i h1 h2 o) where-    Net a b c + Net d e f = Net (a + d) (b + e) (c + f)-    Net a b c - Net d e f = Net (a - d) (b - e) (c - f)-    Net a b c * Net d e f = Net (a * d) (b * e) (c * f)-    abs    (Net a b c)    = Net (abs    a) (abs    b) (abs    c)-    signum (Net a b c)    = Net (signum a) (signum b) (signum c)-    negate (Net a b c)    = Net (negate a) (negate b) (negate c)-    fromInteger x         = Net (fromInteger x) (fromInteger x) (fromInteger x)+instance ( KnownNat i+         , KnownNat h1+         , KnownNat h2+         , KnownNat o+         ) => Num (Network i h1 h2 o) where+    (+)         = gPlus+    (-)         = gMinus+    (*)         = gTimes+    negate      = gNegate+    abs         = gAbs+    signum      = gSignum+    fromInteger = gFromInteger  instance (KnownNat i, KnownNat o) => Fractional (Layer i o) where-    Layer w1 b1 / Layer w2 b2 = Layer (w1 / w2) (b1 / b2)-    recip (Layer w b)         = Layer (recip w) (recip b)-    fromRational x            = Layer (fromRational x) (fromRational x)+    (/)          = gDivide+    recip        = gRecip+    fromRational = gFromRational -instance (KnownNat i, KnownNat h1, KnownNat h2, KnownNat o) => Fractional (Network i h1 h2 o) where-    Net a b c / Net d e f = Net (a / d) (b / e) (c / f)-    recip (Net a b c)     = Net (recip a) (recip b) (recip c)-    fromRational x        = Net (fromRational x) (fromRational x) (fromRational x)+instance ( KnownNat i+         , KnownNat h1+         , KnownNat h2+         , KnownNat o+         ) => Fractional (Network i h1 h2 o) where+    (/)          = gDivide+    recip        = gRecip+    fromRational = gFromRational ``` -`KnownNat` comes from *base*; it’s a typeclass that *hmatrix* uses to+`KnownNat` comes from *base*; it's a typeclass that *hmatrix* uses to refer to the numbers in its type and use it to go about its normal hmatrixy business. @@ -164,13 +212,13 @@ ===  Now, *backprop* does require *primitive* differentiable operations on-our relevant types to be defined. *backprop* uses these primitive `Op`s-to tie everything together. Ideally we’d import these from a library-that implements these for you, and the end-user never has to make `Op`-primitives.+our relevant types to be defined. *backprop* uses these primitive+operations to tie everything together. Ideally we'd import these from a+library that implements these for you, and the end-user never has to+make these primitives. -But in this case, I’m going to put the definitions here to show that-there isn’t any magic going on. If you’re curious, refer to+But in this case, I'm going to put the definitions here to show that+there isn't any magic going on. If you're curious, refer to [documentation for `Op`] for more details on how `Op` is implemented and how this works. @@ -180,225 +228,297 @@ gradient function.  ``` {.sourceCode .literate .haskell}-matVec-    :: (KnownNat m, KnownNat n)-    => Op '[ L m n, R n ] (R m)-matVec = op2' $ \m v ->-  ( m #> v, \(fromMaybe 1 -> g) ->-              (g `outer` v, tr m #> g)-  )+infixr 8 #>!+(#>!)+    :: (KnownNat m, KnownNat n, Reifies s W)+    => BVar s (L m n)+    -> BVar s (R n)+    -> BVar s (R m)+(#>!) = liftOp2 . op2 $ \m v ->+  ( m #> v, \g -> (g `outer` v, tr m #> g) ) ```  Dot products would be nice too.  ``` {.sourceCode .literate .haskell}-dot :: KnownNat n-    => Op '[ R n, R n ] Double-dot = op2' $ \x y ->-  ( x <.> y, \case Nothing -> (y, x)-                   Just g  -> (konst g * y, x * konst g)+infixr 8 <.>!+(<.>!)+    :: (KnownNat n, Reifies s W)+    => BVar s (R n)+    -> BVar s (R n)+    -> BVar s Double+(<.>!) = liftOp2 . op2 $ \x y ->+  ( x <.> y, \g -> (konst g * y, x * konst g)   ) ``` -Also a “scaling” function, scales a vector by a given factor.+Also a function to fill a vector with the same element:  ``` {.sourceCode .literate .haskell}-scale-    :: KnownNat n-    => Op '[ Double, R n ] (R n)-scale = op2' $ \a x ->-  ( konst a * x-  , \case Nothing -> (HM.sumElements (extract x      ), konst a    )-          Just g  -> (HM.sumElements (extract (x * g)), konst a * g)-  )+konst'+    :: (KnownNat n, Reifies s W)+    => BVar s Double+    -> BVar s (R n)+konst' = liftOp1 . op1 $ \c -> (konst c, HM.sumElements . extract) ```  Finally, an operation to sum all of the items in the vector.  ``` {.sourceCode .literate .haskell}-vsum-    :: KnownNat n-    => Op '[ R n ] Double-vsum = op1' $ \x -> (HM.sumElements (extract x), maybe 1 konst)+sumElements'+    :: (KnownNat n, Reifies s W)+    => BVar s (R n)+    -> BVar s Double+sumElements' = liftOp1 . op1 $ \x -> (HM.sumElements (extract x), konst) ``` -And why not, here’s the [logistic function], which we’ll use as an-activation function for internal layers. We don’t need to define this as-an `Op` up-front right now, because the library can automatically-promote any numeric polymorphic function (an `a -> a` or `a -> a -> a`,-etc.) to an `Op` anyways.+Again, these are not intended to be used by end-users of *backprop*, but+rather are meant to be provided by libraries as primitive operations for+users of the library to use. -  [logistic function]: https://en.wikipedia.org/wiki/Logistic_function+Running our Network+=================== +Now that we have our primitives in place, let's actually write a+function to run our network! And, once we do this, we automatically also+have functions to back-propagate our network!++Normally, to write this function, we'd write:+ ``` {.sourceCode .literate .haskell}-logistic :: Floating a => a -> a-logistic x = 1 / (1 + exp (-x))+runLayerNormal+    :: (KnownNat i, KnownNat o)+    => Layer i o+    -> R i+    -> R o+runLayerNormal l x = (l ^. lWeights) #> x + (l ^. lBiases)+{-# INLINE runLayerNormal #-} ``` -Running our Network-===================--Now that we have our primitives in place, let’s actually write a-function to run our network!+Using the `lWeights` and `lBiases` lenses to access the weights and+biases of our layer. However, we can translate this to *backprop* by+operating on `BVar`s instead of the type directly, and using our+backprop-aware `#>!`:  ``` {.sourceCode .literate .haskell} runLayer-    :: (KnownNat i, KnownNat o)-    => BPOp s '[ R i, Layer i o ] (R o)-runLayer = withInps $ \(x :< l :< Ø) -> do-    w :< b :< Ø <- gTuple #<~ l-    y <- matVec ~$ (w :< x :< Ø)-    return $ y + b+    :: (KnownNat i, KnownNat o, Reifies s W)+    => BVar s (Layer i o)+    -> BVar s (R i)+    -> BVar s (R o)+runLayer l x = (l ^^. lWeights) #>! x + (l ^^. lBiases)+{-# INLINE runLayer #-} ``` -A `BPOp s '[ R i, Layer i o ] (R o)` is a backpropagatable function that-produces an `R o` (a vector with `o` elements, from the *[hmatrix]*-library) given an input environment of an `R i` (the “input” of the-layer) and a layer.+`^.` lets to access data within a value using a lens, and `^^.` lets you+access data within a `BVar` using a lens: -  [hmatrix]: http://hackage.haskell.org/package/hmatrix+``` {.haskell}+(^.)  ::        a -> Lens' a b ->        b+(^^.) :: BVar s a -> Lens' a b -> BVar s b+``` -We use `withInps` to bring the environment into scope as a bunch of-`BVar`s. `x` is a `BVar` containing the input vector, and `l` is a-`BVar` containing the layer.+(There is also `^^?`, which can use a `Prism` or `Traversal` to extract+a target that might not exist, `^^..`, which uses a `Traversal` to+extract all targets, and `.~~`, which uses a `Lens` to update a value+inside `BVar`) -The first thing we do is split out the parts of the layer so we can work-with the internal matrices. We can use `#<~` to “split out” the-components of a `BVar`, splitting on `gTuple` (which uses `GHC.Generics`-to automatically figure out how to split up a product type).+Now `runLayer` is a function on two inputs that can be backpropagated,+automatically! We can find its gradient given any input, and also run it+to get our expected output as well. -Then we apply `matVec` (our primitive `Op` that does matrix-vector-multiplication) to `w` and `x`, and then the result is that added to the-bias vector `b`.+Before writing our final network runner, we need a function to compute+the "softmax" of our output vector. Writing it normally would look like: -We can write the `runNetwork` function pretty much the same way.+``` {.sourceCode .literate .haskell}+softMaxNormal :: KnownNat n => R n -> R n+softMaxNormal x = konst (1 / HM.sumElements (extract expx)) * expx+  where+    expx = exp x+{-# INLINE softMaxNormal #-}+``` +But we can make the mechanical shift to the backpropagatable version:+ ``` {.sourceCode .literate .haskell}-runNetwork-    :: (KnownNat i, KnownNat h1, KnownNat h2, KnownNat o)-    => BPOp s '[ R i, Network i h1 h2 o ] (R o)-runNetwork = withInps $ \(x :< n :< Ø) -> do-    l1 :< l2 :< l3 :< Ø <- gTuple #<~ n-    y <- runLayer -$ (x          :< l1 :< Ø)-    z <- runLayer -$ (logistic y :< l2 :< Ø)-    r <- runLayer -$ (logistic z :< l3 :< Ø)-    softmax       -$ (r          :< Ø)+softMax :: (KnownNat n, Reifies s W) => BVar s (R n) -> BVar s (R n)+softMax x = konst' (1 / sumElements' expx) * expx   where-    softmax :: KnownNat n => BPOp s '[ R n ] (R n)-    softmax = withInps $ \(x :< Ø) -> do-        expX <- bindVar (exp x)-        totX <- vsum ~$ (expX   :< Ø)-        scale        ~$ (1/totX :< expX :< Ø)+    expx = exp x+{-# INLINE softMax #-} ``` -After splitting out the layers in the input `Network`, we run each layer-successively using our previously defined `runLayer`, giving inputs-using `-$`. We can directly apply `logistic` to `BVar`s. At the end, we-run a [softmax function] because MNIST is a classification challenge.-The softmax is done by applying $e^x$ for every item in the input-vector, and dividing each element by the total.+We also need the [logistic function], which is our activation function+between layer outputs. Because `BVar`s have a `Floating` instance, we+can just write it using typeclass functions. -  [softmax function]: https://en.wikipedia.org/wiki/Softmax_function+  [logistic function]: https://en.wikipedia.org/wiki/Logistic_function -The Magic----------+``` {.sourceCode .literate .haskell}+logistic :: Floating a => a -> a+logistic x = 1 / (1 + exp (-x))+{-# INLINE logistic #-}+``` -What did we just define? Well, with a `BPOp s rs a`, we can *run* it and-get the output:+With those in hand, let's compare how we would normally write a function+to run our network:  ``` {.sourceCode .literate .haskell}-runNetOnInp+runNetNormal     :: (KnownNat i, KnownNat h1, KnownNat h2, KnownNat o)     => Network i h1 h2 o     -> R i     -> R o-runNetOnInp n x = evalBPOp runNetwork (x ::< n ::< Ø)+runNetNormal n = softMaxNormal+               . runLayerNormal (n ^. nLayer3)+               . logistic+               . runLayerNormal (n ^. nLayer2)+               . logistic+               . runLayerNormal (n ^. nLayer1)+{-# INLINE runNetNormal #-} ``` -But, the magic part is that we can also get the gradient!+Basic function composition, neat. We use our lenses `nLayer1`,+`nLayer2`, and `nLayer3` to extract the first, second, and third layers+from our network. +Writing it in a way that backprop can use is also very similar:+ ``` {.sourceCode .literate .haskell}-gradNet-    :: (KnownNat i, KnownNat h1, KnownNat h2, KnownNat o)-    => Network i h1 h2 o+runNetwork+    :: (KnownNat i, KnownNat h1, KnownNat h2, KnownNat o, Reifies s W)+    => BVar s (Network i h1 h2 o)     -> R i-    -> Network i h1 h2 o-gradNet n x = case gradBPOp runNetwork (x ::< n ::< Ø) of-    _gradX ::< gradN ::< Ø -> gradN+    -> BVar s (R o)+runNetwork n = softMax+             . runLayer (n ^^. nLayer3)+             . logistic+             . runLayer (n ^^. nLayer2)+             . logistic+             . runLayer (n ^^. nLayer1)+             . constVar+{-# INLINE runNetwork #-} ``` -This gives the gradient of all of the parameters in the matrices and-vectors inside the `Network`, which we can use to “train”!+We use `constVar` on the input vector, because we don't care about its+gradient and so treat it as a constant. -Training-========+And now here again we use `^^.` (instead of `^.`) to extract a value+from our `BVar` of a `Network`, using a lens. -Now for the real work. To train a network, we can do gradient descent-based on the gradient of some type of *error function* with respect to-the network parameters. Let’s use the [cross entropy], which is popular-for classification problems.+Computing Errors+---------------- +Now, training a neural network is about calculating its gradient with+respect to some error function. The library calculatues the gradient for+us -- we just need to tell it how to compute the error function.++For classification problems, we usually use a [cross entropy] error.+Given a target vector, how does our neural network's output differ from+what is expected? Lower numbers are better!+   [cross entropy]: https://en.wikipedia.org/wiki/Cross_entropy +Again, let's look at a "normal" implementation, regular variables and no+backprop:+ ``` {.sourceCode .literate .haskell}+crossEntropyNormal :: KnownNat n => R n -> R n -> Double+crossEntropyNormal targ res = -(log res <.> targ)+{-# INLINE crossEntropyNormal #-}+```++And we can see that the backpropable version is pretty similar. We see+`constVar t`, to introduce a `BVar` that is a constant value (that we+don't care about the gradient of).++``` {.sourceCode .literate .haskell} crossEntropy-    :: KnownNat n+    :: (KnownNat n, Reifies s W)     => R n-    -> BPOpI s '[ R n ] Double-crossEntropy targ (r :< Ø) = negate (dot .$ (log r :< t :< Ø))-  where-    t = constVar targ+    -> BVar s (R n)+    -> BVar s Double+crossEntropy targ res = -(log res <.>! constVar targ)+{-# INLINE crossEntropy #-} ``` -Given a target vector and a `BVar` referring to the result of the-network, we can directly apply:+Our final "error function", then, is: -$$-H(\mathbf{r}, \mathbf{t}) = - (log(\mathbf{r}) \cdot \mathbf{t})-$$+``` {.sourceCode .literate .haskell}+netErr+    :: (KnownNat i, KnownNat h1, KnownNat h2, KnownNat o, Reifies s W)+    => R i+    -> R o+    -> BVar s (Network i h1 h2 o)+    -> BVar s Double+netErr x targ n = crossEntropy targ (runNetwork n x)+{-# INLINE netErr #-}+``` -Just for fun, I implemented `crossEntropy` in “implicit-graph” mode, so-you don’t see any binds or returns.+The Magic+========= -Now, a function to make one gradient descent step based on an input-vector and a target, using `gradBPOp`:+The actual "magic" of the library happens with the functions to "run"+the functions we defined earlier: +``` {.haskell}+evalBP   :: (forall s. Reifies s W => BVar s a -> BVar s b) -> a -> b+gradBP   :: (forall s. Reifies s W => BVar s a -> BVar s b) -> a -> a+backprop :: (forall s. Reifies s W => BVar s a -> BVar s b) -> a -> (b, a)+```++`evalBP` "runs" the function like normal, `gradBP` computes the gradient+of the function, and `backprop` computes both the result and the+gradient.++So, if we have a network `net0`, an input vector `x`, and a target+vector `t`, we could compute its error using:++``` {.haskell}+evalBP (netErr x targ) net0 :: Double+```++And we can calculate its *gradient* using:++``` {.haskell}+gradBP (netErr x targ) net0 :: (Network i h1 h2 o, R i)+```++Pulling it all together+=======================++Let's write a simple function to step our network in the direction+opposite of the gradient to train our model:+ ``` {.sourceCode .literate .haskell} trainStep     :: forall i h1 h2 o. (KnownNat i, KnownNat h1, KnownNat h2, KnownNat o)-    => Double-    -> R i-    -> R o-    -> Network i h1 h2 o+    => Double             -- ^ learning rate+    -> R i                -- ^ input+    -> R o                -- ^ target+    -> Network i h1 h2 o  -- ^ initial network     -> Network i h1 h2 o-trainStep r !x !t !n = case gradBPOp o (x ::< n ::< Ø) of-    _ ::< gN ::< Ø ->-        n - (realToFrac r * gN)-  where-    o :: BPOp s '[ R i, Network i h1 h2 o ] Double-    o = do-      y <- runNetwork-      implicitly (crossEntropy t) -$ (y :< Ø)+trainStep r !x !targ !n = n - realToFrac r * gradBP (netErr x targ) n+{-# INLINE trainStep #-} ``` -A convenient wrapper for training over all of the observations in a-list:+Here's a convenient wrapper for training over all of the observations in+a list:  ``` {.sourceCode .literate .haskell} trainList     :: (KnownNat i, KnownNat h1, KnownNat h2, KnownNat o)-    => Double-    -> [(R i, R o)]-    -> Network i h1 h2 o+    => Double             -- ^ learning rate+    -> [(R i, R o)]       -- ^ input and target pairs+    -> Network i h1 h2 o  -- ^ initial network     -> Network i h1 h2 o trainList r = flip $ foldl' (\n (x,y) -> trainStep r x y n)+{-# INLINE trainList #-} ``` -Pulling it all together-=======================- `testNet` will be a quick way to test our net by computing the-percentage of correct guesses: (mostly using *hmatrix* stuff)+percentage of correct guesses: (mostly using *hmatrix* stuff, so don't+mind too much)  ``` {.sourceCode .literate .haskell} testNet@@ -408,13 +528,13 @@     -> Double testNet xs n = sum (map (uncurry test) xs) / fromIntegral (length xs)   where-    test :: R i -> R o -> Double+    test :: R i -> R o -> Double          -- test if the max index is correct     test x (extract->t)         | HM.maxIndex t == HM.maxIndex (extract r) = 1         | otherwise                                = 0       where         r :: R o-        r = evalBPOp runNetwork (x ::< n ::< Ø)+        r = evalBP (`runNetwork` x) n ```  And now, a main loop!@@ -431,7 +551,7 @@     Just train <- loadMNIST "data/train-images-idx3-ubyte" "data/train-labels-idx1-ubyte"     Just test  <- loadMNIST "data/t10k-images-idx3-ubyte"  "data/t10k-labels-idx1-ubyte"     putStrLn "Loaded data."-    net0 <- MWC.uniformR @(Network 784 300 100 9) (-0.5, 0.5) g+    net0 <- MWC.uniformR @(Network 784 300 100 10) (-0.5, 0.5) g     flip evalStateT net0 . forM_ [1..] $ \e -> do       train' <- liftIO . fmap V.toList $ MWC.uniformShuffle (V.fromList train) g       liftIO $ printf "[Epoch %d]\n" (e :: Int)@@ -464,19 +584,24 @@     the test set 5.  Prints out the results -And, that’s really it!+And, that's really it! -Result-------+Performance+----------- -I haven’t put much into optimizing the library yet, but the network-(with hidden layer sizes 300 and 100) seems to take 25s on my computer-to finish a batch of 5000 training points. It’s slow (five minutes per-60000 point epooch), but it’s a first unoptimized run and a proof of-concept! It’s my goal to get this down to a point where the result has-the same performance characteristics as the actual backend (*hmatrix*),-and so overhead is 0.+Currently, benchmarks show that *running* the network has virtually zero+overhead (\~ 4%) over writing the running function directly. The actual+gradient descent process (compute gradient, then descend) carries about+60% overhead over writing the gradients manually, but it is unclear how+much of this is because of the library, and how much of it is just+because of automatic differentation giving slightly less efficient+matrix/vector multiplication operations. +The [README] has some more detailed benchmarks and statistics, if you+want to get more detailed information.++  [README]: https://github.com/mstksg/backprop+ Main takeaways ============== @@ -485,7 +610,7 @@ matrices.  Basically, all that *backprop* did was give you an API to define *how to-run* a neural net — how to *run* a net based on a `Network` and `R i`+run* a neural net --- how to *run* a net based on a `Network` and `R i` input you were given. The goal of the library is to let you write down how to run things in as natural way as possible. @@ -494,18 +619,27 @@  Because the heavy lifting is done by the data types themselves, we can presumably plug in *any* type and any tensor/numerical backend, and reap-the benefits of those libraries’ optimizations and parallelizations.+the benefits of those libraries' optimizations and parallelizations. *Any* type can be backpropagated! :D  What now? ---------  Check out the docs for the [Numeric.Backprop] module for a more detailed-picture of what’s going on, or find more examples at the [github repo]!+picture of what's going on, or find more examples at the [github repo]!    [Numeric.Backprop]: http://hackage.haskell.org/package/backprop/docs/Numeric-Backprop.html   [github repo]: https://github.com/mstksg/backprop +Also, check out follow-up writeup to this tutorial, expanding on using+the library with more advanced extensible neural network types, like the+ones described in [this blog post]. Check out the [literate haskell+here], and the [rendered PDF here].++  [this blog post]: https://blog.jle.im/entries/series/+practical-dependent-types-in-haskell.html+  [literate haskell here]: https://github.com/mstksg/backprop/blob/master/samples/extensible-neural.lhs+  [rendered PDF here]: https://github.com/mstksg/backprop/blob/master/renders/extensible-neural.pdf+ Boring stuff ============ @@ -518,7 +652,7 @@ loadMNIST     :: FilePath     -> FilePath-    -> IO (Maybe [(R 784, R 9)])+    -> IO (Maybe [(R 784, R 10)]) loadMNIST fpI fpL = runMaybeT $ do     i <- MaybeT          $ decodeIDXFile       fpI     l <- MaybeT          $ decodeIDXLabelsFile fpL@@ -528,8 +662,8 @@   where     mkImage :: VU.Vector Int -> Maybe (R 784)     mkImage = create . VG.convert . VG.map (\i -> fromIntegral i / 255)-    mkLabel :: Int -> Maybe (R 9)-    mkLabel n = create $ HM.build 9 (\i -> if round i == n then 1 else 0)+    mkLabel :: Int -> Maybe (R 10)+    mkLabel n = create $ HM.build 10 (\i -> if round i == n then 1 else 0) ```  And here are instances to generating random@@ -548,7 +682,12 @@     uniform g = Layer <$> MWC.uniform g <*> MWC.uniform g     uniformR (l, h) g = (\x -> x * (h - l) + l) <$> MWC.uniform g -instance (KnownNat i, KnownNat h1, KnownNat h2, KnownNat o) => MWC.Variate (Network i h1 h2 o) where+instance ( KnownNat i+         , KnownNat h1+         , KnownNat h2+         , KnownNat o+         )+      => MWC.Variate (Network i h1 h2 o) where     uniform g = Net <$> MWC.uniform g <*> MWC.uniform g <*> MWC.uniform g     uniformR (l, h) g = (\x -> x * (h - l) + l) <$> MWC.uniform g ```

binary file changed (143090 → 133909 bytes)

@@ -1,447 +0,0 @@-----author:-- Justin Le-fontfamily: 'palatino,cmtt'-geometry: margin=1in-links-as-notes: true-title: Neural networks with backprop library------The *backprop* library performs back-propagation over a *hetereogeneous*-system of relationships. It offers both an implicit ([ad]-like) and-explicit graph building usage style. Let’s use it to build neural-networks!--  [ad]: http://hackage.haskell.org/package/ad--Repository source is [on github], and so are the [rendered unstable-docs].--  [on github]: https://github.com/mstksg/backprop-  [rendered unstable docs]: https://mstksg.github.io/backprop--``` {.sourceCode .literate .haskell}-{-# LANGUAGE DeriveGeneric                 #-}-{-# LANGUAGE GADTs                         #-}-{-# LANGUAGE LambdaCase                    #-}-{-# LANGUAGE RankNTypes                    #-}-{-# LANGUAGE ScopedTypeVariables           #-}-{-# LANGUAGE StandaloneDeriving            #-}-{-# LANGUAGE TypeApplications              #-}-{-# LANGUAGE TypeInType                    #-}-{-# LANGUAGE TypeOperators                 #-}-{-# LANGUAGE ViewPatterns                  #-}-{-# OPTIONS_GHC -fno-warn-orphans          #-}-{-# OPTIONS_GHC -fno-warn-unused-top-binds #-}--import           Data.Functor-import           Data.Kind-import           Data.Maybe-import           Data.Singletons-import           Data.Singletons.Prelude-import           Data.Singletons.TypeLits-import           Data.Type.Combinator-import           Data.Type.Product-import           GHC.Generics                        (Generic)-import           Numeric.Backprop-import           Numeric.Backprop.Iso-import           Numeric.LinearAlgebra.Static hiding (dot)-import           System.Random.MWC-import qualified Generics.SOP                        as SOP-```--Ops-===--First, we define values of `Op` for the operations we want to do. `Op`s-are bundles of functions packaged with their hetereogeneous gradients.-For simple numeric functions, *backprop* can derive `Op`s automatically.-But for matrix operations, we have to derive them ourselves.--The types help us with matching up the dimensions, but we still need to-be careful that our gradients are calculated correctly.--`L` and `R` are matrix and vector types from the great *hmatrix*-library.--First, matrix-vector multiplication:--``` {.sourceCode .literate .haskell}-matVec-    :: (KnownNat m, KnownNat n)-    => Op '[ L m n, R n ] (R m)-matVec = op2' $ \m v -> ( m #> v-                        , \(fromMaybe 1 -> g) ->-                             (g `outer` v, tr m #> g)-                        )-```--Now, dot products:--``` {.sourceCode .literate .haskell}-dot :: KnownNat n-    => Op '[ R n, R n ] Double-dot = op2' $ \x y -> ( x <.> y-                     , \case Nothing -> (y, x)-                             Just g  -> (konst g * y, x * konst g)-                     )-```--Polymorphic functions can be easily turned into `Op`s with `op1`/`op2`-etc., but they can also be run directly on graph nodes.--``` {.sourceCode .literate .haskell}-logistic :: Floating a => a -> a-logistic x = 1 / (1 + exp (-x))-```--A Simple Complete Example-=========================--At this point, we already have enough to train a simple-single-hidden-layer neural network:--``` {.sourceCode .literate .haskell}-simpleOp-      :: (KnownNat m, KnownNat n, KnownNat o)-      => R m-      -> BPOpI s '[ L n m, R n, L o n, R o ] (R o)-simpleOp inp = \(w1 :< b1 :< w2 :< b2 :< Ø) ->-    let z = logistic $ liftB2 matVec w1 x + b1-    in  logistic $ liftB2 matVec w2 z + b2-  where-    x = constVar inp-```--Here, `simpleOp` is defined in implicit (non-monadic) style, given a-tuple of inputs and returning outputs. Now `simpleOp` can be “run” with-the input vectors and parameters (a `L n m`, `R n`, `L o n`, and `R o`)-and calculate the output of the neural net.--``` {.sourceCode .literate .haskell}-runSimple-    :: (KnownNat m, KnownNat n, KnownNat o)-    => R m-    -> Tuple '[ L n m, R n, L o n, R o ]-    -> R o-runSimple inp = evalBPOp (implicitly $ simpleOp inp)-```--Alternatively, we can define `simpleOp` in explicit monadic style, were-we specify our graph nodes explicitly. The results should be the same.--``` {.sourceCode .literate .haskell}-simpleOpExplicit-      :: (KnownNat m, KnownNat n, KnownNat o)-      => R m-      -> BPOp s '[ L n m, R n, L o n, R o ] (R o)-simpleOpExplicit inp = withInps $ \(w1 :< b1 :< w2 :< b2 :< Ø) -> do-    -- First layer-    y1  <- matVec ~$ (w1 :< x1 :< Ø)-    let x2 = logistic (y1 + b1)-    -- Second layer-    y2  <- matVec ~$ (w2 :< x2 :< Ø)-    return $ logistic (y2 + b2)-  where-    x1 = constVar inp-```--Now, for the magic of *backprop*: the library can now take advantage of-the implicit (or explicit) graph and use it to do back-propagation, too!--``` {.sourceCode .literate .haskell}-simpleGrad-    :: forall m n o. (KnownNat m, KnownNat n, KnownNat o)-    => R m-    -> R o-    -> Tuple '[ L n m, R n, L o n, R o ]-    -> Tuple '[ L n m, R n, L o n, R o ]-simpleGrad inp targ params = gradBPOp opError params-  where-    opError :: BPOp s '[ L n m, R n, L o n, R o ] Double-    opError = do-        res <- implicitly $ simpleOp inp-        -- we explicitly bind err to prevent recomputation-        err <- bindVar $ res - t-        dot ~$ (err :< err :< Ø)-      where-        t = constVar targ-```--The result is the gradient of the input tuple’s components, with respect-to the `Double` result of `opError` (the squared error). We can then use-this gradient to do gradient descent.--With Parameter Containers-=========================--This method doesn’t quite scale, because we might want to make networks-with multiple layers and parameterize networks by layers. Let’s make-some basic container data types to help us organize our types, including-a recursive `Network` type that lets us chain multiple layers.--``` {.sourceCode .literate .haskell}-data Layer :: Nat -> Nat -> Type where-    Layer :: { _lWeights :: L m n-             , _lBiases  :: R m-             }-          -> Layer n m-      deriving (Show, Generic)---data Network :: Nat -> [Nat] -> Nat -> Type where-    NØ   :: !(Layer a b) -> Network a '[] b-    (:&) :: !(Layer a b) -> Network b bs c -> Network a (b ': bs) c-```--A `Layer n m` is a layer taking an n-vector and returning an m-vector. A-`Network a '[b, c, d] e` would be a Network that takes in an a-vector-and outputs an e-vector, with hidden layers of sizes b, c, and d.--Isomorphisms---------------The *backprop* library lets you apply operations on “parts” of data-types (like on the weights and biases of a `Layer`) by using `Iso`’s-(isomorphisms), like the ones from the *lens* library. The library-doesn’t depend on lens, but it can use the `Iso`s from the library and-also custom-defined ones.--First, we can auto-generate isomorphisms using the *generics-sop*-library:--``` {.sourceCode .literate .haskell}-instance SOP.Generic (Layer n m)-```--And then can create isomorphisms by hand for the two `Network`-constructors:--``` {.sourceCode .literate .haskell}-netExternal :: Iso' (Network a '[] b) (Tuple '[Layer a b])-netExternal = iso (\case NØ x     -> x ::< Ø)-                  (\case I x :< Ø -> NØ x   )--netInternal :: Iso' (Network a (b ': bs) c) (Tuple '[Layer a b, Network b bs c])-netInternal = iso (\case x :& xs          -> x ::< xs ::< Ø)-                  (\case I x :< I xs :< Ø -> x :& xs       )-```--An `Iso' a (Tuple as)` means that an `a` can really just be seen as a-tuple of `as`.--Running a network-=================--Now, we can write the `BPOp` that reprenents running the network and-getting a result. We pass in a `Sing bs` (a singleton list of the hidden-layer sizes) so that we can “pattern match” on the list and handle the-different network constructors differently.--``` {.sourceCode .literate .haskell}-netOp-    :: forall s a bs c. (KnownNat a, KnownNat c)-    => Sing bs-    -> BPOp s '[ R a, Network a bs c ] (R c)-netOp sbs = go sbs-  where-    go :: forall d es. KnownNat d-        => Sing es-        -> BPOp s '[ R d, Network d es c ] (R c)-    go = \case-      SNil -> withInps $ \(x :< n :< Ø) -> do-        -- peek into the NØ using netExternal iso-        l :< Ø <- netExternal #<~ n-        -- run the 'layerOp' BP, with x and l as inputs-        bpOp layerOp ~$ (x :< l :< Ø)-      SNat `SCons` ses -> withInps $ \(x :< n :< Ø) -> withSingI ses $ do-        -- peek into the (:&) using the netInternal iso-        l :< n' :< Ø <- netInternal #<~ n-        -- run the 'layerOp' BP, with x and l as inputs-        z <- bpOp layerOp  ~$ (x :< l :< Ø)-        -- run the 'go ses' BP, with z and n as inputs-        bpOp (go ses)      ~$ (z :< n' :< Ø)-    layerOp-        :: forall d e. (KnownNat d, KnownNat e)-        => BPOp s '[ R d, Layer d e ] (R e)-    layerOp = withInps $ \(x :< l :< Ø) -> do-        -- peek into the layer using the gTuple iso, auto-generated with SOP.Generic-        w :< b :< Ø <- gTuple #<~ l-        y           <- matVec  ~$ (w :< x :< Ø)-        return $ logistic (y + b)-```--There’s some singletons work going on here, but it’s fairly standard-singletons stuff. Most of the complexity here is from the static typing-in our neural network type, and *not* from *backprop*.--From *backprop* specifically, the only elements are `#<~` lets you-“split” an input ref with the given iso, and `bpOp`, which converts a-`BPOp` into an `Op` that you can bind with `~$`.--Note that this library doesn’t support truly pattern matching on GADTs,-and that we had to pass in `Sing bs` as a reference to the structure of-our networks.--Gradient Descent-------------------Now we can do simple gradient descent. Defining an error function:--``` {.sourceCode .literate .haskell}-errOp-    :: KnownNat m-    => R m-    -> BVar s rs (R m)-    -> BPOp s rs Double-errOp targ r = do-    err <- bindVar $ r - t-    dot ~$ (err :< err :< Ø)-  where-    t = constVar targ-```--And now, we can use `backprop` to generate the gradient, and shift the-`Network`! Things are made a bit cleaner from the fact that-`Network a bs c` has a `Num` instance, so we can use `(-)` and `(*)`-etc.--``` {.sourceCode .literate .haskell}-train-    :: (KnownNat a, SingI bs, KnownNat c)-    => Double-    -> R a-    -> R c-    -> Network a bs c-    -> Network a bs c-train r x t n = case backprop (errOp t =<< netOp sing) (x ::< n ::< Ø) of-    (_, _ :< I g :< Ø) -> n - (realToFrac r * g)-```--(`(::<)` is cons and `Ø` is nil for tuples.)--Main-====--`main`, which will train on sample data sets, is still in progress!-Right now it just generates a random network using the *mwc-random*-library and prints each internal layer.--``` {.sourceCode .literate .haskell}-main :: IO ()-main = withSystemRandom $ \g -> do-    n <- uniform @(Network 4 '[3,2] 1) g-    void $ traverseNetwork sing (\l -> l <$ print l) n-```--Appendix: Boilerplate-=====================--And now for some typeclass instances and boilerplates unrelated to the-*backprop* library that makes our custom types easier to use.--``` {.sourceCode .literate .haskell}-instance KnownNat n => Variate (R n) where-    uniform g = randomVector <$> uniform g <*> pure Uniform-    uniformR (l, h) g = (\x -> x * (h - l) + l) <$> uniform g--instance (KnownNat m, KnownNat n) => Variate (L m n) where-    uniform g = uniformSample <$> uniform g <*> pure 0 <*> pure 1-    uniformR (l, h) g = (\x -> x * (h - l) + l) <$> uniform g--instance (KnownNat n, KnownNat m) => Variate (Layer n m) where-    uniform g = subtract 1 . (* 2) <$> (Layer <$> uniform g <*> uniform g)-    uniformR (l, h) g = (\x -> x * (h - l) + l) <$> uniform g--instance (KnownNat m, KnownNat n) => Num (Layer n m) where-    Layer w1 b1 + Layer w2 b2 = Layer (w1 + w2) (b1 + b2)-    Layer w1 b1 - Layer w2 b2 = Layer (w1 - w2) (b1 - b2)-    Layer w1 b1 * Layer w2 b2 = Layer (w1 * w2) (b1 * b2)-    abs    (Layer w b) = Layer (abs w) (abs b)-    signum (Layer w b) = Layer (signum w) (signum b)-    negate (Layer w b) = Layer (negate w) (negate b)-    fromInteger x = Layer (fromInteger x) (fromInteger x)--instance (KnownNat m, KnownNat n) => Fractional (Layer n m) where-    Layer w1 b1 / Layer w2 b2 = Layer (w1 / w2) (b1 / b2)-    recip (Layer w b) = Layer (recip w) (recip b)-    fromRational x = Layer (fromRational x) (fromRational x)--instance (KnownNat a, SingI bs, KnownNat c) => Variate (Network a bs c) where-    uniform g = genNet sing (uniform g)-    uniformR (l, h) g = (\x -> x * (h - l) + l) <$> uniform g--genNet-    :: forall f a bs c. (Applicative f, KnownNat a, KnownNat c)-    => Sing bs-    -> (forall d e. (KnownNat d, KnownNat e) => f (Layer d e))-    -> f (Network a bs c)-genNet sbs f = go sbs-  where-    go :: forall d es. KnownNat d => Sing es -> f (Network d es c)-    go = \case-      SNil             -> NØ <$> f-      SNat `SCons` ses -> (:&) <$> f <*> go ses--mapNetwork0-    :: forall a bs c. (KnownNat a, KnownNat c)-    => Sing bs-    -> (forall d e. (KnownNat d, KnownNat e) => Layer d e)-    -> Network a bs c-mapNetwork0 sbs f = getI $ genNet sbs (I f)--traverseNetwork-    :: forall a bs c f. (KnownNat a, KnownNat c, Applicative f)-    => Sing bs-    -> (forall d e. (KnownNat d, KnownNat e) => Layer d e -> f (Layer d e))-    -> Network a bs c-    -> f (Network a bs c)-traverseNetwork sbs f = go sbs-  where-    go :: forall d es. KnownNat d => Sing es -> Network d es c -> f (Network d es c)-    go = \case-      SNil -> \case-        NØ x -> NØ <$> f x-      SNat `SCons` ses -> \case-        x :& xs -> (:&) <$> f x <*> go ses xs--mapNetwork1-    :: forall a bs c. (KnownNat a, KnownNat c)-    => Sing bs-    -> (forall d e. (KnownNat d, KnownNat e) => Layer d e -> Layer d e)-    -> Network a bs c-    -> Network a bs c-mapNetwork1 sbs f = getI . traverseNetwork sbs (I . f)--mapNetwork2-    :: forall a bs c. (KnownNat a, KnownNat c)-    => Sing bs-    -> (forall d e. (KnownNat d, KnownNat e) => Layer d e -> Layer d e -> Layer d e)-    -> Network a bs c-    -> Network a bs c-    -> Network a bs c-mapNetwork2 sbs f = go sbs-  where-    go :: forall d es. KnownNat d => Sing es -> Network d es c -> Network d es c -> Network d es c-    go = \case-      SNil -> \case-        NØ x -> \case-          NØ y -> NØ (f x y)-      SNat `SCons` ses -> \case-        x :& xs -> \case-          y :& ys -> f x y :& go ses xs ys--instance (KnownNat a, SingI bs, KnownNat c) => Num (Network a bs c) where-    (+)           = mapNetwork2 sing (+)-    (-)           = mapNetwork2 sing (-)-    (*)           = mapNetwork2 sing (*)-    negate        = mapNetwork1 sing negate-    abs           = mapNetwork1 sing abs-    signum        = mapNetwork1 sing signum-    fromInteger x = mapNetwork0 sing (fromInteger x)--instance (KnownNat a, SingI bs, KnownNat c) => Fractional (Network a bs c) where-    (/)            = mapNetwork2 sing (/)-    recip          = mapNetwork1 sing recip-    fromRational x = mapNetwork0 sing (fromRational x)-```

binary file changed (104794 → absent bytes)

@@ -0,0 +1,552 @@+---+author:+- Justin Le+title: Extensible Neural Networks with Backprop+---++This write-up is a follow-up to the *MNIST* tutorial ([rendered] here,+and [literate haskell] here). This write-up itself is available as a+[literate haskell file], and also [rendered as a pdf].++  [rendered]: https://github.com/mstksg/backprop/blob/master/renders/backprop-mnist.pdf+  [literate haskell]: https://github.com/mstksg/backprop/blob/master/samples/backprop-mnist.lhs+  [literate haskell file]: https://github.com/mstksg/backprop/blob/master/samples/extensible-neural.lhs+  [rendered as a pdf]: https://github.com/mstksg/backprop/blob/master/renders/extensible-neural.pdf++The packages involved are:++-   deepseq+-   hmatrix+-   lens+-   mnist-idx+-   mwc-random+-   one-liner-instances+-   singletons+-   split+-   vector++``` {.sourceCode .literate .haskell}+{-# LANGUAGE BangPatterns         #-}+{-# LANGUAGE DataKinds            #-}+{-# LANGUAGE DeriveGeneric        #-}+{-# LANGUAGE FlexibleContexts     #-}+{-# LANGUAGE GADTs                #-}+{-# LANGUAGE InstanceSigs         #-}+{-# LANGUAGE LambdaCase           #-}+{-# LANGUAGE LambdaCase           #-}+{-# LANGUAGE RankNTypes           #-}+{-# LANGUAGE ScopedTypeVariables  #-}+{-# LANGUAGE TemplateHaskell      #-}+{-# LANGUAGE TypeApplications     #-}+{-# LANGUAGE TypeInType           #-}+{-# LANGUAGE TypeOperators        #-}+{-# LANGUAGE ViewPatterns         #-}+{-# OPTIONS_GHC -fno-warn-orphans #-}++import           Control.DeepSeq+import           Control.Exception+import           Control.Lens hiding             ((<.>))+import           Control.Monad+import           Control.Monad.IO.Class+import           Control.Monad.Primitive+import           Control.Monad.Trans.Maybe+import           Control.Monad.Trans.State+import           Data.Bitraversable+import           Data.Foldable+import           Data.IDX+import           Data.Kind+import           Data.List.Split+import           Data.Singletons+import           Data.Singletons.Prelude+import           Data.Singletons.TypeLits+import           Data.Time.Clock+import           Data.Traversable+import           Data.Tuple+import           GHC.Generics                    (Generic)+import           Numeric.Backprop+import           Numeric.LinearAlgebra.Static+import           Numeric.OneLiner+import           Text.Printf+import qualified Data.Vector                     as V+import qualified Data.Vector.Generic             as VG+import qualified Data.Vector.Unboxed             as VU+import qualified Numeric.LinearAlgebra           as HM+import qualified System.Random.MWC               as MWC+import qualified System.Random.MWC.Distributions as MWC+```++Introduction+============++The *[backprop]* library lets us manipulate our values in a natural way.+We write the function to compute our result, and the library then+automatically finds the *gradient* of that function, which we can use+for gradient descent.++  [backprop]: http://hackage.haskell.org/package/backprop++In the last post, we looked at using a fixed-structure neural network.+However, in [this blog series], I discuss a system of extensible neural+networks that can be chained and composed.++  [this blog series]: https://blog.jle.im/entries/series/+practical-dependent-types-in-haskell.html++One issue, however, in naively translating the implementations, is that+we normally run the network by pattern matching on each layer. However,+we cannot directly pattern match on `BVar`s.++We *could* get around it by being smart with prisms and `^^?`, to+extract a "Maybe BVar". However, we can do better! This is because the+*shape* of a `Net i hs o` is known already at compile-time, so there is+no need for runtime checks like prisms and `^^?`.++Instead, we can just directly use lenses, since we know *exactly* what+constructor will be present! We can use singletons to determine which+constructor is present, and so always just directly use lenses without+any runtime nondeterminism.++Types+=====++First, our types:++``` {.sourceCode .literate .haskell}+data Layer i o =+    Layer { _lWeights :: !(L o i)+          , _lBiases  :: !(R o)+          }+  deriving (Show, Generic)++instance NFData (Layer i o)+makeLenses ''Layer++data Net :: Nat -> [Nat] -> Nat -> Type where+    NO   :: !(Layer i o) -> Net i '[] o+    (:~) :: !(Layer i h) -> !(Net h hs o) -> Net i (h ': hs) o+```++Unfortunately, we can't automatically generate lenses for GADTs, so we+have to make them by hand.\[\^poly\]++with type safety via paraemtric polymorphism.++``` {.sourceCode .literate .haskell}+_NO :: Lens (Net i '[] o) (Net i' '[] o')+            (Layer i o  ) (Layer i' o'  )+_NO f (NO l) = NO <$> f l++_NIL :: Lens (Net i (h ': hs) o) (Net i' (h ': hs) o)+             (Layer i h        ) (Layer i' h        )+_NIL f (l :~ n) = (:~ n) <$> f l++_NIN :: Lens (Net i (h ': hs) o) (Net i (h ': hs') o')+             (Net h hs        o) (Net h hs'        o')+_NIN f (l :~ n) = (l :~) <$> f n+```++You can read `_NO` as:++``` {.haskell}+_NO :: Lens' (Net i '[] o) (Layer i o)+```++A lens into a single-layer network, and++``` {.haskell}+_NIL :: Lens' (Net i (h ': hs) o) (Layer i h )+_NIN :: Lens' (Net i (h ': hs) o) (Net h hs o)+```++Lenses into a multiple-layer network, getting the first layer and the+tail of the network.++If we pattern match on `Sing hs`, we can always determine exactly which+lenses we can use, and so never fumble around with prisms or+nondeterminism.++Running the network+===================++Here's the meat of process, then: specifying how to run the network. We+re-use our `BVar`-based combinators defined in the last write-up:++``` {.sourceCode .literate .haskell}+runLayer+    :: (KnownNat i, KnownNat o, Reifies s W)+    => BVar s (Layer i o)+    -> BVar s (R i)+    -> BVar s (R o)+runLayer l x = (l ^^. lWeights) #>! x + (l ^^. lBiases)+{-# INLINE runLayer #-}+```++For `runNetwork`, we pattern match on `hs` using singletons, so we+always know exactly what type of network we have:++``` {.sourceCode .literate .haskell}+runNetwork+    :: (KnownNat i, KnownNat o, Reifies s W)+    => BVar s (Net i hs o)+    -> Sing hs+    -> BVar s (R i)+    -> BVar s (R o)+runNetwork n = \case+    SNil          -> softMax . runLayer (n ^^. _NO)+    SCons SNat hs -> withSingI hs (runNetwork (n ^^. _NIN) hs)+                   . logistic+                   . runLayer (n ^^. _NIL)+{-# INLINE runNetwork #-}+```++The rest of it is the same as before.++``` {.sourceCode .literate .haskell}+netErr+    :: (KnownNat i, KnownNat o, SingI hs, Reifies s W)+    => R i+    -> R o+    -> BVar s (Net i hs o)+    -> BVar s Double+netErr x targ n = crossEntropy targ (runNetwork n sing (constVar x))+{-# INLINE netErr #-}++trainStep+    :: forall i hs o. (KnownNat i, KnownNat o, SingI hs)+    => Double             -- ^ learning rate+    -> R i                -- ^ input+    -> R o                -- ^ target+    -> Net i hs o         -- ^ initial network+    -> Net i hs o+trainStep r !x !targ !n = n - realToFrac r * gradBP (netErr x targ) n+{-# INLINE trainStep #-}++trainList+    :: (KnownNat i, SingI hs, KnownNat o)+    => Double             -- ^ learning rate+    -> [(R i, R o)]       -- ^ input and target pairs+    -> Net i hs o         -- ^ initial network+    -> Net i hs o+trainList r = flip $ foldl' (\n (x,y) -> trainStep r x y n)+{-# INLINE trainList #-}++testNet+    :: forall i hs o. (KnownNat i, KnownNat o, SingI hs)+    => [(R i, R o)]+    -> Net i hs o+    -> Double+testNet xs n = sum (map (uncurry test) xs) / fromIntegral (length xs)+  where+    test :: R i -> R o -> Double          -- test if the max index is correct+    test x (extract->t)+        | HM.maxIndex t == HM.maxIndex (extract r) = 1+        | otherwise                                = 0+      where+        r :: R o+        r = evalBP (\n' -> runNetwork n' sing (constVar x)) n+```++And that's it!++Running+=======++Everything here is the same as before, except now we can dynamically+pick the network size. Here we pick `'[300,100]` for the hidden layer+sizes.++``` {.sourceCode .literate .haskell}+main :: IO ()+main = MWC.withSystemRandom $ \g -> do+    Just train <- loadMNIST "data/train-images-idx3-ubyte" "data/train-labels-idx1-ubyte"+    Just test  <- loadMNIST "data/t10k-images-idx3-ubyte"  "data/t10k-labels-idx1-ubyte"+    putStrLn "Loaded data."+    net0 <- MWC.uniformR @(Net 784 '[300,100] 10) (-0.5, 0.5) g+    flip evalStateT net0 . forM_ [1..] $ \e -> do+      train' <- liftIO . fmap V.toList $ MWC.uniformShuffle (V.fromList train) g+      liftIO $ printf "[Epoch %d]\n" (e :: Int)++      forM_ ([1..] `zip` chunksOf batch train') $ \(b, chnk) -> StateT $ \n0 -> do+        printf "(Batch %d)\n" (b :: Int)++        t0 <- getCurrentTime+        n' <- evaluate . force $ trainList rate chnk n0+        t1 <- getCurrentTime+        printf "Trained on %d points in %s.\n" batch (show (t1 `diffUTCTime` t0))++        let trainScore = testNet chnk n'+            testScore  = testNet test n'+        printf "Training error:   %.2f%%\n" ((1 - trainScore) * 100)+        printf "Validation error: %.2f%%\n" ((1 - testScore ) * 100)++        return ((), n')+  where+    rate  = 0.02+    batch = 5000+```++Looking Forward+===============++One common thing people might do is want to be able to mix different+types of layers. This could also be easily encoded as different+constructors in `Layer`, and so `runLayer` will now be different+depending on what constructor is present.++In this case, we can either:++1.  Have a different indexed type for layers, so that we can always know+    exactly what layer is involved, so we don't have to runtime pattern+    match:++    ``` {.haskell}+    data LayerType = FullyConnected | Convolutional++    data Layer :: LayerType -> Nat -> Nat -> Type where+        LayerFC :: .... -> Layer 'FullyConnected i o+        LayerC  :: .... -> Layer 'Convolutional  i o+    ```++    We would then have `runLayer` take `Sing (t :: LayerType)`, so we+    can again use `^^.` and directly pattern match.++2.  Use a typeclass-based approach, so users can add their own layer+    types. In this situation, layer types would all be different types,+    and running them would be a typeclass method that would give our+    `BVar s (Layer i o) -> BVar s (R i) -> BVar s (R o)` operation as a+    typeclass method.++    ``` {.haskell}+    class Layer (l :: Nat -> Nat -> Type) where+        runLayer+            :: forall s. Reifies s W+            => BVar s (l i o)+            -> BVar s (R i)+            -> BVar s (R o)+    ```++In all cases, it shouldn't be much more cognitive overhead to use+*backprop* to build your neural network framework!++And, remember that `evalBP` (directly running the function) introduces+virtually zero overhead, so if you only provided `BVar` functions, you+could easily get the original non-`BVar` functions with `evalBP` without+any loss.++What now?+---------++Ready to start? Check out the docs for the [Numeric.Backprop] module for+the full technical specs, and find more examples and updates at the+[github repo]!++  [Numeric.Backprop]: http://hackage.haskell.org/package/backprop/docs/Numeric-Backprop.html+  [github repo]: https://github.com/mstksg/backprop++Internals+=========++That's it for the post! Now for the internal plumbing :)++``` {.sourceCode .literate .haskell}+loadMNIST+    :: FilePath+    -> FilePath+    -> IO (Maybe [(R 784, R 10)])+loadMNIST fpI fpL = runMaybeT $ do+    i <- MaybeT          $ decodeIDXFile       fpI+    l <- MaybeT          $ decodeIDXLabelsFile fpL+    d <- MaybeT . return $ labeledIntData l i+    r <- MaybeT . return $ for d (bitraverse mkImage mkLabel . swap)+    liftIO . evaluate $ force r+  where+    mkImage :: VU.Vector Int -> Maybe (R 784)+    mkImage = create . VG.convert . VG.map (\i -> fromIntegral i / 255)+    mkLabel :: Int -> Maybe (R 10)+    mkLabel n = create $ HM.build 10 (\i -> if round i == n then 1 else 0)+```++HMatrix Operations+------------------++``` {.sourceCode .literate .haskell}+infixr 8 #>!+(#>!)+    :: (KnownNat m, KnownNat n, Reifies s W)+    => BVar s (L m n)+    -> BVar s (R n)+    -> BVar s (R m)+(#>!) = liftOp2 . op2 $ \m v ->+  ( m #> v, \g -> (g `outer` v, tr m #> g) )++infixr 8 <.>!+(<.>!)+    :: (KnownNat n, Reifies s W)+    => BVar s (R n)+    -> BVar s (R n)+    -> BVar s Double+(<.>!) = liftOp2 . op2 $ \x y ->+  ( x <.> y, \g -> (konst g * y, x * konst g)+  )++konst'+    :: (KnownNat n, Reifies s W)+    => BVar s Double+    -> BVar s (R n)+konst' = liftOp1 . op1 $ \c -> (konst c, HM.sumElements . extract)++sumElements'+    :: (KnownNat n, Reifies s W)+    => BVar s (R n)+    -> BVar s Double+sumElements' = liftOp1 . op1 $ \x -> (HM.sumElements (extract x), konst)++softMax :: (KnownNat n, Reifies s W) => BVar s (R n) -> BVar s (R n)+softMax x = konst' (1 / sumElements' expx) * expx+  where+    expx = exp x+{-# INLINE softMax #-}++crossEntropy+    :: (KnownNat n, Reifies s W)+    => R n+    -> BVar s (R n)+    -> BVar s Double+crossEntropy targ res = -(log res <.>! constVar targ)+{-# INLINE crossEntropy #-}++logistic :: Floating a => a -> a+logistic x = 1 / (1 + exp (-x))+{-# INLINE logistic #-}+```++Instances+---------++``` {.sourceCode .literate .haskell}+instance (KnownNat i, KnownNat o) => Num (Layer i o) where+    (+)         = gPlus+    (-)         = gMinus+    (*)         = gTimes+    negate      = gNegate+    abs         = gAbs+    signum      = gSignum+    fromInteger = gFromInteger++instance (KnownNat i, KnownNat o) => Fractional (Layer i o) where+    (/)          = gDivide+    recip        = gRecip+    fromRational = gFromRational+++liftNet0+    :: forall i hs o. (KnownNat i, KnownNat o)+    => (forall m n. (KnownNat m, KnownNat n) => Layer m n)+    -> Sing hs+    -> Net i hs o+liftNet0 x = go+  where+    go :: forall w ws. KnownNat w => Sing ws -> Net w ws o+    go = \case+      SNil          -> NO x+      SCons SNat hs -> x :~ go hs++liftNet1+    :: forall i hs o. (KnownNat i, KnownNat o)+    => (forall m n. (KnownNat m, KnownNat n)+          => Layer m n+          -> Layer m n+       )+    -> Sing hs+    -> Net i hs o+    -> Net i hs o+liftNet1 f = go+  where+    go  :: forall w ws. KnownNat w+        => Sing ws+        -> Net w ws o+        -> Net w ws o+    go = \case+      SNil          -> \case+        NO x -> NO (f x)+      SCons SNat hs -> \case+        x :~ xs -> f x :~ go hs xs++liftNet2+    :: forall i hs o. (KnownNat i, KnownNat o)+    => (forall m n. (KnownNat m, KnownNat n)+          => Layer m n+          -> Layer m n+          -> Layer m n+       )+    -> Sing hs+    -> Net i hs o+    -> Net i hs o+    -> Net i hs o+liftNet2 f = go+  where+    go  :: forall w ws. KnownNat w+        => Sing ws+        -> Net w ws o+        -> Net w ws o+        -> Net w ws o+    go = \case+      SNil          -> \case+        NO x -> \case+          NO y -> NO (f x y)+      SCons SNat hs -> \case+        x :~ xs -> \case+          y :~ ys -> f x y :~ go hs xs ys++instance ( KnownNat i+         , KnownNat o+         , SingI hs+         )+      => Num (Net i hs o) where+    (+)           = liftNet2 (+) sing+    (-)           = liftNet2 (-) sing+    (*)           = liftNet2 (*) sing+    negate        = liftNet1 negate sing+    abs           = liftNet1 abs sing+    signum        = liftNet1 signum sing+    fromInteger x = liftNet0 (fromInteger x) sing++instance ( KnownNat i+         , KnownNat o+         , SingI hs+         )+      => Fractional (Net i hs o) where+    (/)            = liftNet2 (/) sing+    recip          = liftNet1 negate sing+    fromRational x = liftNet0 (fromRational x) sing++instance KnownNat n => MWC.Variate (R n) where+    uniform g = randomVector <$> MWC.uniform g <*> pure Uniform+    uniformR (l, h) g = (\x -> x * (h - l) + l) <$> MWC.uniform g++instance (KnownNat m, KnownNat n) => MWC.Variate (L m n) where+    uniform g = uniformSample <$> MWC.uniform g <*> pure 0 <*> pure 1+    uniformR (l, h) g = (\x -> x * (h - l) + l) <$> MWC.uniform g++instance (KnownNat i, KnownNat o) => MWC.Variate (Layer i o) where+    uniform g = Layer <$> MWC.uniform g <*> MWC.uniform g+    uniformR (l, h) g = (\x -> x * (h - l) + l) <$> MWC.uniform g++instance ( KnownNat i+         , KnownNat o+         , SingI hs+         )+      => MWC.Variate (Net i hs o) where+    uniform :: forall m. PrimMonad m => MWC.Gen (PrimState m) -> m (Net i hs o)+    uniform g = go sing+      where+        go :: forall w ws. KnownNat w => Sing ws -> m (Net w ws o)+        go = \case+          SNil          -> NO <$> MWC.uniform g+          SCons SNat hs -> (:~) <$> MWC.uniform g <*> go hs+    uniformR (l, h) g = (\x -> x * (h - l) + l) <$> MWC.uniform g++instance NFData (Net i hs o) where+    rnf = \case+      NO l    -> rnf l+      x :~ xs -> rnf x `seq` rnf xs+```

binary file changed (absent → 108272 bytes)

@@ -2,11 +2,12 @@ % Justin Le  The *backprop* library performs back-propagation over a *hetereogeneous*-system of relationships.  It offers both an implicit (*[ad][]*-like) and explicit graph-building usage style.  Let's use it to build neural networks and learn-mnist!+system of relationships.  back-propagation is done automatically (as+reverse-mode automatic differentiation), and you work with your values as if+you were writing normal functions with them, with the help of [lens][].  [ad]: http://hackage.haskell.org/package/ad+[lens]: http://hackage.haskell.org/package/lens  Repository source is [on github][repo], and docs are [on hackage][hackage]. @@ -21,22 +22,35 @@ [rendered]: https://github.com/mstksg/backprop/blob/master/renders/backprop-mnist.pdf [lhs]: https://github.com/mstksg/backprop/blob/master/samples/backprop-mnist.lhs +The packages involved are: +*   deepseq+*   hmatrix+*   lens+*   mnist-idx+*   mwc-random+*   one-liner-instances+*   split+*   vector+ > {-# LANGUAGE BangPatterns                     #-} > {-# LANGUAGE DataKinds                        #-} > {-# LANGUAGE DeriveGeneric                    #-}+> {-# LANGUAGE FlexibleContexts                 #-} > {-# LANGUAGE GADTs                            #-} > {-# LANGUAGE LambdaCase                       #-} > {-# LANGUAGE ScopedTypeVariables              #-}+> {-# LANGUAGE TemplateHaskell                  #-} > {-# LANGUAGE TupleSections                    #-} > {-# LANGUAGE TypeApplications                 #-} > {-# LANGUAGE ViewPatterns                     #-}-> {-# OPTIONS_GHC -fno-warn-orphans             #-} > {-# OPTIONS_GHC -fno-warn-incomplete-patterns #-}+> {-# OPTIONS_GHC -fno-warn-orphans             #-} > {-# OPTIONS_GHC -fno-warn-unused-top-binds    #-} > > import           Control.DeepSeq > import           Control.Exception+> import           Control.Lens hiding ((<.>)) > import           Control.Monad > import           Control.Monad.IO.Class > import           Control.Monad.Trans.Maybe@@ -45,23 +59,52 @@ > import           Data.Foldable > import           Data.IDX > import           Data.List.Split-> import           Data.Maybe > import           Data.Time.Clock > import           Data.Traversable > import           Data.Tuple > import           GHC.Generics                        (Generic) > import           GHC.TypeLits > import           Numeric.Backprop-> import           Numeric.LinearAlgebra.Static hiding (dot)+> import           Numeric.LinearAlgebra.Static+> import           Numeric.OneLiner > import           Text.Printf > import qualified Data.Vector                         as V > import qualified Data.Vector.Generic                 as VG > import qualified Data.Vector.Unboxed                 as VU-> import qualified Generics.SOP                        as SOP > import qualified Numeric.LinearAlgebra               as HM > import qualified System.Random.MWC                   as MWC > import qualified System.Random.MWC.Distributions     as MWC +Introduction+============++In this walkthrough, we'll be building a classifier for the *[MNIST][]* data+set.  This is meant to mirror the [Tensorflow Tutorial][tf-intro] for+beginners.++[tf-intro]: https://www.tensorflow.org/versions/r1.2/get_started/mnist/beginners++Essentially, we use a two-layer artificial neural network -- or a series of+matrix multiplications, differentiable function applications, and vector+additions.  We feed our input image to the ANN and then try to get a label+from it.  Training an ANN is a matter of finding the right matrices to+multiply by, and the right vectors to add.++To do that, we train our network by treating our network's accuracy as a+function `Network -> Error`.  If we can find the gradient of the input network+with respect to the error, we can perform [gradient descent][], and slowly+make our network better and better.++[gradient descent]: https://en.wikipedia.org/wiki/Gradient_descent++Finding the gradient is usually complicated, but *backprop* makes it simpler:++1.  Write a function to compute the error from the network+2.  That's it!++Hooray!  Once you do that, the library finds the gradient function+*automatically*, without any further intervention!+ Types ===== @@ -82,8 +125,8 @@ >           } >   deriving (Show, Generic) >-> instance SOP.Generic (Layer i o) > instance NFData (Layer i o)+> makeLenses ''Layer  And a type for a simple feed-forward network with two hidden layers: @@ -94,53 +137,58 @@ >         } >   deriving (Show, Generic) >-> instance SOP.Generic (Network i h1 h2 o) > instance NFData (Network i h1 h2 o)+> makeLenses ''Network  These are pretty straightforward container types...pretty much exactly the type you'd make to represent these networks!  Note that, following true Haskell form, we separate out logic from data.  This should be all we need. -We derive an instance of `SOP.Generic` from the *[generics-sop][]* package,-which *backprop* uses to propagate derivatives on values inside product-types.--[generics-sop]: http://hackage.haskell.org/package/generics-sop- Instances ---------  Things are much simplier if we had `Num` and `Fractional` instances for everything, so let's just go ahead and define that now, as well.  Just a-little bit of boilerplate.+little bit of boilerplate, made easier using *[one-liner-instances][]* to+auto-derive instances using Generics. +[one-liner-instances]: http://hackage.haskell.org/package/one-liner-instances+ > instance (KnownNat i, KnownNat o) => Num (Layer i o) where->     Layer w1 b1 + Layer w2 b2 = Layer (w1 + w2) (b1 + b2)->     Layer w1 b1 - Layer w2 b2 = Layer (w1 - w2) (b1 - b2)->     Layer w1 b1 * Layer w2 b2 = Layer (w1 * w2) (b1 * b2)->     abs    (Layer w b)        = Layer (abs    w) (abs    b)->     signum (Layer w b)        = Layer (signum w) (signum b)->     negate (Layer w b)        = Layer (negate w) (negate b)->     fromInteger x             = Layer (fromInteger x) (fromInteger x)+>     (+)         = gPlus+>     (-)         = gMinus+>     (*)         = gTimes+>     negate      = gNegate+>     abs         = gAbs+>     signum      = gSignum+>     fromInteger = gFromInteger >-> instance (KnownNat i, KnownNat h1, KnownNat h2, KnownNat o) => Num (Network i h1 h2 o) where->     Net a b c + Net d e f = Net (a + d) (b + e) (c + f)->     Net a b c - Net d e f = Net (a - d) (b - e) (c - f)->     Net a b c * Net d e f = Net (a * d) (b * e) (c * f)->     abs    (Net a b c)    = Net (abs    a) (abs    b) (abs    c)->     signum (Net a b c)    = Net (signum a) (signum b) (signum c)->     negate (Net a b c)    = Net (negate a) (negate b) (negate c)->     fromInteger x         = Net (fromInteger x) (fromInteger x) (fromInteger x)+> instance ( KnownNat i+>          , KnownNat h1+>          , KnownNat h2+>          , KnownNat o+>          ) => Num (Network i h1 h2 o) where+>     (+)         = gPlus+>     (-)         = gMinus+>     (*)         = gTimes+>     negate      = gNegate+>     abs         = gAbs+>     signum      = gSignum+>     fromInteger = gFromInteger > > instance (KnownNat i, KnownNat o) => Fractional (Layer i o) where->     Layer w1 b1 / Layer w2 b2 = Layer (w1 / w2) (b1 / b2)->     recip (Layer w b)         = Layer (recip w) (recip b)->     fromRational x            = Layer (fromRational x) (fromRational x)+>     (/)          = gDivide+>     recip        = gRecip+>     fromRational = gFromRational >-> instance (KnownNat i, KnownNat h1, KnownNat h2, KnownNat o) => Fractional (Network i h1 h2 o) where->     Net a b c / Net d e f = Net (a / d) (b / e) (c / f)->     recip (Net a b c)     = Net (recip a) (recip b) (recip c)->     fromRational x        = Net (fromRational x) (fromRational x) (fromRational x)+> instance ( KnownNat i+>          , KnownNat h1+>          , KnownNat h2+>          , KnownNat o+>          ) => Fractional (Network i h1 h2 o) where+>     (/)          = gDivide+>     recip        = gRecip+>     fromRational = gFromRational  `KnownNat` comes from *base*; it's a typeclass that *hmatrix* uses to refer to the numbers in its type and use it to go about its normal hmatrixy@@ -150,216 +198,275 @@ ===  Now, *backprop* does require *primitive* differentiable operations on our-relevant types to be defined.  *backprop* uses these primitive `Op`s to tie-everything together.  Ideally we'd import these from a library that-implements these for you, and the end-user never has to make `Op`-primitives.+relevant types to be defined.  *backprop* uses these primitive operations to+tie everything together.  Ideally we'd import these from a library that+implements these for you, and the end-user never has to make these primitives.  But in this case, I'm going to put the definitions here to show that there isn't any magic going on.  If you're curious, refer to [documentation for-`Op`][opdoc] for more details on how `Op` is implemented and how this-works.+`Op`][opdoc] for more details on how `Op` is implemented and how this works.  [opdoc]: http://hackage.haskell.org/package/backprop/docs/Numeric-Backprop-Op.html  First, matrix-vector multiplication primitive, giving an explicit gradient function. -> matVec->     :: (KnownNat m, KnownNat n)->     => Op '[ L m n, R n ] (R m)-> matVec = op2' $ \m v ->->   ( m #> v, \(fromMaybe 1 -> g) ->->               (g `outer` v, tr m #> g)->   )+> infixr 8 #>!+> (#>!)+>     :: (KnownNat m, KnownNat n, Reifies s W)+>     => BVar s (L m n)+>     -> BVar s (R n)+>     -> BVar s (R m)+> (#>!) = liftOp2 . op2 $ \m v ->+>   ( m #> v, \g -> (g `outer` v, tr m #> g) )  Dot products would be nice too. -> dot :: KnownNat n->     => Op '[ R n, R n ] Double-> dot = op2' $ \x y ->->   ( x <.> y, \case Nothing -> (y, x)->                    Just g  -> (konst g * y, x * konst g)+> infixr 8 <.>!+> (<.>!)+>     :: (KnownNat n, Reifies s W)+>     => BVar s (R n)+>     -> BVar s (R n)+>     -> BVar s Double+> (<.>!) = liftOp2 . op2 $ \x y ->+>   ( x <.> y, \g -> (konst g * y, x * konst g) >   ) -Also a "scaling" function, scales a vector by a given factor.+Also a function to fill a vector with the same element: -> scale->     :: KnownNat n->     => Op '[ Double, R n ] (R n)-> scale = op2' $ \a x ->->   ( konst a * x->   , \case Nothing -> (HM.sumElements (extract x      ), konst a    )->           Just g  -> (HM.sumElements (extract (x * g)), konst a * g)->   )+> konst'+>     :: (KnownNat n, Reifies s W)+>     => BVar s Double+>     -> BVar s (R n)+> konst' = liftOp1 . op1 $ \c -> (konst c, HM.sumElements . extract)  Finally, an operation to sum all of the items in the vector. -> vsum->     :: KnownNat n->     => Op '[ R n ] Double-> vsum = op1' $ \x -> (HM.sumElements (extract x), maybe 1 konst)--And why not, here's the [logistic function][], which we'll use as an-activation function for internal layers.  We don't need to define this as-an `Op` up-front right now, because the library can automatically promote-any numeric polymorphic function (an `a -> a` or `a -> a -> a`, etc.) to an-`Op` anyways.--[logistic function]: https://en.wikipedia.org/wiki/Logistic_function+> sumElements'+>     :: (KnownNat n, Reifies s W)+>     => BVar s (R n)+>     -> BVar s Double+> sumElements' = liftOp1 . op1 $ \x -> (HM.sumElements (extract x), konst) -> logistic :: Floating a => a -> a-> logistic x = 1 / (1 + exp (-x))+Again, these are not intended to be used by end-users of *backprop*, but+rather are meant to be provided by libraries as primitive operations for users+of the library to use.  Running our Network ===================  Now that we have our primitives in place, let's actually write a function-to run our network!+to run our network!  And, once we do this, we automatically also have+functions to back-propagate our network! -> runLayer+Normally, to write this function, we'd write:++> runLayerNormal >     :: (KnownNat i, KnownNat o)->     => BPOp s '[ R i, Layer i o ] (R o)-> runLayer = withInps $ \(x :< l :< Ø) -> do->     w :< b :< Ø <- gTuple #<~ l->     y <- matVec ~$ (w :< x :< Ø)->     return $ y + b+>     => Layer i o+>     -> R i+>     -> R o+> runLayerNormal l x = (l ^. lWeights) #> x + (l ^. lBiases)+> {-# INLINE runLayerNormal #-} -A `BPOp s '[ R i, Layer i o ] (R o)` is a backpropagatable function that-produces an `R o` (a vector with `o` elements, from the *[hmatrix][]*-library) given an input environment of an `R i` (the "input" of the layer)-and a layer.+Using the `lWeights` and `lBiases` lenses to access the weights and biases of+our layer.  However, we can translate this to *backprop* by operating on+`BVar`s instead of the type directly, and using our backprop-aware `#>!`: -We use `withInps` to bring the environment into scope as a bunch of-`BVar`s.  `x` is a `BVar` containing the input vector, and `l` is a `BVar`-containing the layer.+> runLayer+>     :: (KnownNat i, KnownNat o, Reifies s W)+>     => BVar s (Layer i o)+>     -> BVar s (R i)+>     -> BVar s (R o)+> runLayer l x = (l ^^. lWeights) #>! x + (l ^^. lBiases)+> {-# INLINE runLayer #-} -The first thing we do is split out the parts of the layer so we can work-with the internal matrices.  We can use `#<~` to "split out" the components-of a `BVar`, splitting on `gTuple` (which uses `GHC.Generics` to-automatically figure out how to split up a product type).+`^.` lets to access data within a value using a lens, and `^^.` lets you+access data within a `BVar` using a lens: -Then we apply `matVec` (our primitive `Op` that does matrix-vector-multiplication) to `w` and `x`, and then the result is that added to the-bias vector `b`.+```haskell+(^.)  ::        a -> Lens' a b ->        b+(^^.) :: BVar s a -> Lens' a b -> BVar s b+``` -We can write the `runNetwork` function pretty much the same way.+(There is also `^^?`, which can use a `Prism` or `Traversal` to extract a+target that might not exist, `^^..`, which uses a `Traversal` to extract all+targets, and `.~~`, which uses a `Lens` to update a value inside `BVar`) -> runNetwork->     :: (KnownNat i, KnownNat h1, KnownNat h2, KnownNat o)->     => BPOp s '[ R i, Network i h1 h2 o ] (R o)-> runNetwork = withInps $ \(x :< n :< Ø) -> do->     l1 :< l2 :< l3 :< Ø <- gTuple #<~ n->     y <- runLayer -$ (x          :< l1 :< Ø)->     z <- runLayer -$ (logistic y :< l2 :< Ø)->     r <- runLayer -$ (logistic z :< l3 :< Ø)->     softmax       -$ (r          :< Ø)+Now `runLayer` is a function on two inputs that can be backpropagated,+automatically!  We can find its gradient given any input, and also run it to+get our expected output as well.++Before writing our final network runner, we need a function to compute the+"softmax" of our output vector.  Writing it normally would look like:++> softMaxNormal :: KnownNat n => R n -> R n+> softMaxNormal x = konst (1 / HM.sumElements (extract expx)) * expx >   where->     softmax :: KnownNat n => BPOp s '[ R n ] (R n)->     softmax = withInps $ \(x :< Ø) -> do->         expX <- bindVar (exp x)->         totX <- vsum ~$ (expX   :< Ø)->         scale        ~$ (1/totX :< expX :< Ø)+>     expx = exp x+> {-# INLINE softMaxNormal #-} +But we can make the mechanical shift to the backpropagatable version: -After splitting out the layers in the input `Network`, we run each layer-successively using our previously defined `runLayer`, giving inputs using-`-$`.  We can directly apply `logistic` to `BVar`s.  At the end, we run a-[softmax function][] because MNIST is a classification challenge.  The softmax-is done by applying $e^x$ for every item in the input vector, and dividing-each element by the total.+> softMax :: (KnownNat n, Reifies s W) => BVar s (R n) -> BVar s (R n)+> softMax x = konst' (1 / sumElements' expx) * expx+>   where+>     expx = exp x+> {-# INLINE softMax #-} -[softmax function]: https://en.wikipedia.org/wiki/Softmax_function+We also need the [logistic function][], which is our activation function+between layer outputs. Because `BVar`s have a `Floating` instance, we can just+write it using typeclass functions. +[logistic function]: https://en.wikipedia.org/wiki/Logistic_function -The Magic----------+> logistic :: Floating a => a -> a+> logistic x = 1 / (1 + exp (-x))+> {-# INLINE logistic #-} -What did we just define?  Well, with a `BPOp s rs a`, we can *run* it and-get the output:+With those in hand, let's compare how we would normally write a function to run+our network: -> runNetOnInp+> runNetNormal >     :: (KnownNat i, KnownNat h1, KnownNat h2, KnownNat o) >     => Network i h1 h2 o >     -> R i >     -> R o-> runNetOnInp n x = evalBPOp runNetwork (x ::< n ::< Ø)+> runNetNormal n = softMaxNormal+>                . runLayerNormal (n ^. nLayer3)+>                . logistic+>                . runLayerNormal (n ^. nLayer2)+>                . logistic+>                . runLayerNormal (n ^. nLayer1)+> {-# INLINE runNetNormal #-} -But, the magic part is that we can also get the gradient!+Basic function composition, neat.  We use our lenses `nLayer1`, `nLayer2`, and+`nLayer3` to extract the first, second, and third layers from our network. -> gradNet->     :: (KnownNat i, KnownNat h1, KnownNat h2, KnownNat o)->     => Network i h1 h2 o+Writing it in a way that backprop can use is also very similar:++> runNetwork+>     :: (KnownNat i, KnownNat h1, KnownNat h2, KnownNat o, Reifies s W)+>     => BVar s (Network i h1 h2 o) >     -> R i->     -> Network i h1 h2 o-> gradNet n x = case gradBPOp runNetwork (x ::< n ::< Ø) of->     _gradX ::< gradN ::< Ø -> gradN+>     -> BVar s (R o)+> runNetwork n = softMax+>              . runLayer (n ^^. nLayer3)+>              . logistic+>              . runLayer (n ^^. nLayer2)+>              . logistic+>              . runLayer (n ^^. nLayer1)+>              . constVar+> {-# INLINE runNetwork #-} -This gives the gradient of all of the parameters in the matrices and-vectors inside the `Network`, which we can use to "train"!+We use `constVar` on the input vector, because we don't care about its+gradient and so treat it as a constant. -Training-========+And now here again we use `^^.` (instead of `^.`) to extract a value from our+`BVar` of a `Network`, using a lens. -Now for the real work.  To train a network, we can do gradient descent-based on the gradient of some type of *error function* with respect to the-network parameters.  Let's use the [cross entropy][], which is popular for-classification problems.+Computing Errors+---------------- +Now, training a neural network is about calculating its gradient with respect+to some error function.  The library calculatues the gradient for us -- we+just need to tell it how to compute the error function.++For classification problems, we usually use a [cross entropy][] error.  Given+a target vector, how does our neural network's output differ from what is+expected?  Lower numbers are better!+ [cross entropy]: https://en.wikipedia.org/wiki/Cross_entropy +Again, let's look at a "normal" implementation, regular variables and no+backprop:++> crossEntropyNormal :: KnownNat n => R n -> R n -> Double+> crossEntropyNormal targ res = -(log res <.> targ)+> {-# INLINE crossEntropyNormal #-}++And we can see that the backpropable version is pretty similar.  We see+`constVar t`, to introduce a `BVar` that is a constant value (that we don't+care about the gradient of).+ > crossEntropy->     :: KnownNat n+>     :: (KnownNat n, Reifies s W) >     => R n->     -> BPOpI s '[ R n ] Double-> crossEntropy targ (r :< Ø) = negate (dot .$ (log r :< t :< Ø))->   where->     t = constVar targ+>     -> BVar s (R n)+>     -> BVar s Double+> crossEntropy targ res = -(log res <.>! constVar targ)+> {-# INLINE crossEntropy #-} -Given a target vector and a `BVar` referring to the result of the network,-we can directly apply:+Our final "error function", then, is: -$$-H(\mathbf{r}, \mathbf{t}) = - (log(\mathbf{r}) \cdot \mathbf{t})-$$+> netErr+>     :: (KnownNat i, KnownNat h1, KnownNat h2, KnownNat o, Reifies s W)+>     => R i+>     -> R o+>     -> BVar s (Network i h1 h2 o)+>     -> BVar s Double+> netErr x targ n = crossEntropy targ (runNetwork n x)+> {-# INLINE netErr #-} -Just for fun, I implemented `crossEntropy` in "implicit-graph" mode, so you-don't see any binds or returns.+The Magic+========= -Now, a function to make one gradient descent step based on an input vector-and a target, using `gradBPOp`:+The actual "magic" of the library happens with the functions to "run" the+functions we defined earlier: +```haskell+evalBP   :: (forall s. Reifies s W => BVar s a -> BVar s b) -> a -> b+gradBP   :: (forall s. Reifies s W => BVar s a -> BVar s b) -> a -> a+backprop :: (forall s. Reifies s W => BVar s a -> BVar s b) -> a -> (b, a)+```++`evalBP` "runs" the function like normal, `gradBP` computes the gradient of+the function, and `backprop` computes both the result and the gradient.++So, if we have a network `net0`, an input vector `x`, and a target vector `t`,+we could compute its error using:++```haskell+evalBP (netErr x targ) net0 :: Double+```++And we can calculate its *gradient* using:++```haskell+gradBP (netErr x targ) net0 :: (Network i h1 h2 o, R i)+```++Pulling it all together+=======================++Let's write a simple function to step our network in the direction opposite of+the gradient to train our model:+ > trainStep >     :: forall i h1 h2 o. (KnownNat i, KnownNat h1, KnownNat h2, KnownNat o)->     => Double->     -> R i->     -> R o->     -> Network i h1 h2 o+>     => Double             -- ^ learning rate+>     -> R i                -- ^ input+>     -> R o                -- ^ target+>     -> Network i h1 h2 o  -- ^ initial network >     -> Network i h1 h2 o-> trainStep r !x !t !n = case gradBPOp o (x ::< n ::< Ø) of->     _ ::< gN ::< Ø ->->         n - (realToFrac r * gN)->   where->     o :: BPOp s '[ R i, Network i h1 h2 o ] Double->     o = do->       y <- runNetwork->       implicitly (crossEntropy t) -$ (y :< Ø)+> trainStep r !x !targ !n = n - realToFrac r * gradBP (netErr x targ) n+> {-# INLINE trainStep #-} -A convenient wrapper for training over all of the observations in a list:+Here's a convenient wrapper for training over all of the observations in a+list:  > trainList >     :: (KnownNat i, KnownNat h1, KnownNat h2, KnownNat o)->     => Double->     -> [(R i, R o)]->     -> Network i h1 h2 o+>     => Double             -- ^ learning rate+>     -> [(R i, R o)]       -- ^ input and target pairs+>     -> Network i h1 h2 o  -- ^ initial network >     -> Network i h1 h2 o > trainList r = flip $ foldl' (\n (x,y) -> trainStep r x y n)+> {-# INLINE trainList #-} -Pulling it all together-=======================  `testNet` will be a quick way to test our net by computing the percentage-of correct guesses: (mostly using *hmatrix* stuff)+of correct guesses: (mostly using *hmatrix* stuff, so don't mind too much)  > testNet >     :: forall i h1 h2 o. (KnownNat i, KnownNat h1, KnownNat h2, KnownNat o)@@ -368,28 +475,27 @@ >     -> Double > testNet xs n = sum (map (uncurry test) xs) / fromIntegral (length xs) >   where->     test :: R i -> R o -> Double+>     test :: R i -> R o -> Double          -- test if the max index is correct >     test x (extract->t) >         | HM.maxIndex t == HM.maxIndex (extract r) = 1 >         | otherwise                                = 0 >       where >         r :: R o->         r = evalBPOp runNetwork (x ::< n ::< Ø)+>         r = evalBP (`runNetwork` x) n  And now, a main loop! -If you are following along at home, download the [mnist data set-files][mnist] and uncompress them into the folder `data`, and everything-should work fine.+If you are following along at home, download the [mnist data set files][MNIST]+and uncompress them into the folder `data`, and everything should work fine. -[mnist]: http://yann.lecun.com/exdb/mnist/+[MNIST]: http://yann.lecun.com/exdb/mnist/  > main :: IO () > main = MWC.withSystemRandom $ \g -> do >     Just train <- loadMNIST "data/train-images-idx3-ubyte" "data/train-labels-idx1-ubyte" >     Just test  <- loadMNIST "data/t10k-images-idx3-ubyte"  "data/t10k-labels-idx1-ubyte" >     putStrLn "Loaded data."->     net0 <- MWC.uniformR @(Network 784 300 100 9) (-0.5, 0.5) g+>     net0 <- MWC.uniformR @(Network 784 300 100 10) (-0.5, 0.5) g >     flip evalStateT net0 . forM_ [1..] $ \e -> do >       train' <- liftIO . fmap V.toList $ MWC.uniformShuffle (V.fromList train) g >       liftIO $ printf "[Epoch %d]\n" (e :: Int)@@ -423,16 +529,20 @@  And, that's really it! -Result-------+Performance+----------- -I haven't put much into optimizing the library yet, but the network (with-hidden layer sizes 300 and 100) seems to take 25s on my computer to finish-a batch of 5000 training points.  It's slow (five minutes per 60000 point-epooch), but it's a first unoptimized run and a proof of concept!  It's my-goal to get this down to a point where the result has the same performance-characteristics as the actual backend (*hmatrix*), and so overhead is 0.+Currently, benchmarks show that *running* the network has virtually zero+overhead (~ 4%) over writing the running function directly.  The actual+gradient descent process (compute gradient, then descend) carries about 60%+overhead over writing the gradients manually, but it is unclear how much of+this is because of the library, and how much of it is just because of+automatic differentation giving slightly less efficient matrix/vector+multiplication operations. +The [README][repo] has some more detailed benchmarks and statistics, if you+want to get more detailed information.+ Main takeaways ============== @@ -460,6 +570,16 @@  [Numeric.Backprop]: http://hackage.haskell.org/package/backprop/docs/Numeric-Backprop.html +Also, check out follow-up writeup to this tutorial, expanding on using the+library with more advanced extensible neural network types, like the ones+described in [this blog post][blog].  Check out the [literate haskell+here][neural-lhs], and the [rendered PDF here][neural-pdf].++[blog]: https://blog.jle.im/entries/series/+practical-dependent-types-in-haskell.html++[neural-lhs]: https://github.com/mstksg/backprop/blob/master/samples/extensible-neural.lhs+[neural-pdf]: https://github.com/mstksg/backprop/blob/master/renders/extensible-neural.pdf+ Boring stuff ============ @@ -471,7 +591,7 @@ > loadMNIST >     :: FilePath >     -> FilePath->     -> IO (Maybe [(R 784, R 9)])+>     -> IO (Maybe [(R 784, R 10)]) > loadMNIST fpI fpL = runMaybeT $ do >     i <- MaybeT          $ decodeIDXFile       fpI >     l <- MaybeT          $ decodeIDXLabelsFile fpL@@ -481,8 +601,8 @@ >   where >     mkImage :: VU.Vector Int -> Maybe (R 784) >     mkImage = create . VG.convert . VG.map (\i -> fromIntegral i / 255)->     mkLabel :: Int -> Maybe (R 9)->     mkLabel n = create $ HM.build 9 (\i -> if round i == n then 1 else 0)+>     mkLabel :: Int -> Maybe (R 10)+>     mkLabel n = create $ HM.build 10 (\i -> if round i == n then 1 else 0)  And here are instances to generating random vectors/matrices/layers/networks, used for the initialization step.@@ -499,6 +619,12 @@ >     uniform g = Layer <$> MWC.uniform g <*> MWC.uniform g >     uniformR (l, h) g = (\x -> x * (h - l) + l) <$> MWC.uniform g >-> instance (KnownNat i, KnownNat h1, KnownNat h2, KnownNat o) => MWC.Variate (Network i h1 h2 o) where+> instance ( KnownNat i+>          , KnownNat h1+>          , KnownNat h2+>          , KnownNat o+>          )+>       => MWC.Variate (Network i h1 h2 o) where >     uniform g = Net <$> MWC.uniform g <*> MWC.uniform g <*> MWC.uniform g >     uniformR (l, h) g = (\x -> x * (h - l) + l) <$> MWC.uniform g+

@@ -1,18 +0,0 @@-{-# LANGUAGE GADTs #-}--import           Numeric.Backprop.Mono--testImplicit :: BPOp s N3 Double Double-testImplicit = implicitly $ \(x :* y :* z :* ØV) ->-    ((x * y) + y) * z--testExplicit :: BPOp s N3 Double Double-testExplicit = withInps $ \(x :* y :* z :* ØV) -> do-    xy  <- op2 (*) ~$ (x   :* y :* ØV)-    xyy <- op2 (+) ~$ (xy  :* y :* ØV)-    op2 (*)        ~$ (xyy :* z :* ØV)--main :: IO ()-main = do-    print $ backprop testImplicit (2 :+ 3 :+ 4 :+ ØV)-    print $ backprop testExplicit (2 :+ 3 :+ 4 :+ ØV)

@@ -1,405 +0,0 @@-% Neural networks with backprop library-% Justin Le--The *backprop* library performs back-propagation over a *hetereogeneous*-system of relationships.  It offers both an implicit ([ad][]-like) and explicit graph-building usage style.  Let's use it to build neural networks!--[ad]: http://hackage.haskell.org/package/ad--Repository source is [on github][repo], and so are the [rendered unstable-docs][docs].--[repo]: https://github.com/mstksg/backprop-[docs]: https://mstksg.github.io/backprop--> {-# LANGUAGE DeriveGeneric                 #-}-> {-# LANGUAGE GADTs                         #-}-> {-# LANGUAGE LambdaCase                    #-}-> {-# LANGUAGE RankNTypes                    #-}-> {-# LANGUAGE ScopedTypeVariables           #-}-> {-# LANGUAGE StandaloneDeriving            #-}-> {-# LANGUAGE TypeApplications              #-}-> {-# LANGUAGE TypeInType                    #-}-> {-# LANGUAGE TypeOperators                 #-}-> {-# LANGUAGE ViewPatterns                  #-}-> {-# OPTIONS_GHC -fno-warn-orphans          #-}-> {-# OPTIONS_GHC -fno-warn-unused-top-binds #-}-> -> import           Data.Functor-> import           Data.Kind-> import           Data.Maybe-> import           Data.Singletons-> import           Data.Singletons.Prelude-> import           Data.Singletons.TypeLits-> import           Data.Type.Combinator-> import           Data.Type.Product-> import           GHC.Generics                        (Generic)-> import           Numeric.Backprop-> import           Numeric.Backprop.Iso-> import           Numeric.LinearAlgebra.Static hiding (dot)-> import           System.Random.MWC-> import qualified Generics.SOP                        as SOP--Ops-===--First, we define values of `Op` for the operations we want to do.  `Op`s-are bundles of functions packaged with their hetereogeneous gradients.  For-simple numeric functions, *backprop* can derive `Op`s automatically.  But-for matrix operations, we have to derive them ourselves.--The types help us with matching up the dimensions, but we still need to be-careful that our gradients are calculated correctly.--`L` and `R` are matrix and vector types from the great *hmatrix* library.--First, matrix-vector multiplication:--> matVec->     :: (KnownNat m, KnownNat n)->     => Op '[ L m n, R n ] (R m)-> matVec = op2' $ \m v -> ( m #> v->                         , \(fromMaybe 1 -> g) ->->                              (g `outer` v, tr m #> g)->                         )--Now, dot products:--> dot :: KnownNat n->     => Op '[ R n, R n ] Double-> dot = op2' $ \x y -> ( x <.> y->                      , \case Nothing -> (y, x)->                              Just g  -> (konst g * y, x * konst g)->                      )--Polymorphic functions can be easily turned into `Op`s with `op1`/`op2`-etc., but they can also be run directly on graph nodes.--> logistic :: Floating a => a -> a-> logistic x = 1 / (1 + exp (-x))--A Simple Complete Example-=========================--At this point, we already have enough to train a simple single-hidden-layer-neural network:--> simpleOp->       :: (KnownNat m, KnownNat n, KnownNat o)->       => R m->       -> BPOpI s '[ L n m, R n, L o n, R o ] (R o)-> simpleOp inp = \(w1 :< b1 :< w2 :< b2 :< Ø) ->->     let z = logistic $ liftB2 matVec w1 x + b1->     in  logistic $ liftB2 matVec w2 z + b2->   where->     x = constVar inp--Here, `simpleOp` is defined in implicit (non-monadic) style, given a tuple-of inputs and returning outputs.  Now `simpleOp` can be "run" with the-input vectors and parameters (a `L n m`, `R n`, `L o n`, and `R o`) and-calculate the output of the neural net.--> runSimple->     :: (KnownNat m, KnownNat n, KnownNat o)->     => R m->     -> Tuple '[ L n m, R n, L o n, R o ]->     -> R o-> runSimple inp = evalBPOp (implicitly $ simpleOp inp)--Alternatively, we can define `simpleOp` in explicit monadic style, were we-specify our graph nodes explicitly.  The results should be the same.--> simpleOpExplicit->       :: (KnownNat m, KnownNat n, KnownNat o)->       => R m->       -> BPOp s '[ L n m, R n, L o n, R o ] (R o)-> simpleOpExplicit inp = withInps $ \(w1 :< b1 :< w2 :< b2 :< Ø) -> do->     -- First layer->     y1  <- matVec ~$ (w1 :< x1 :< Ø)->     let x2 = logistic (y1 + b1)->     -- Second layer->     y2  <- matVec ~$ (w2 :< x2 :< Ø)->     return $ logistic (y2 + b2)->   where->     x1 = constVar inp--Now, for the magic of *backprop*:  the library can now take advantage of-the implicit (or explicit) graph and use it to do back-propagation, too!--> simpleGrad->     :: forall m n o. (KnownNat m, KnownNat n, KnownNat o)->     => R m->     -> R o->     -> Tuple '[ L n m, R n, L o n, R o ]->     -> Tuple '[ L n m, R n, L o n, R o ]-> simpleGrad inp targ params = gradBPOp opError params->   where->     opError :: BPOp s '[ L n m, R n, L o n, R o ] Double->     opError = do->         res <- implicitly $ simpleOp inp->         -- we explicitly bind err to prevent recomputation->         err <- bindVar $ res - t->         dot ~$ (err :< err :< Ø)->       where->         t = constVar targ--The result is the gradient of the input tuple's components, with respect-to the `Double` result of `opError` (the squared error).  We can then use-this gradient to do gradient descent.--With Parameter Containers-=========================--This method doesn't quite scale, because we might want to make networks-with multiple layers and parameterize networks by layers.  Let's make some-basic container data types to help us organize our types, including a-recursive `Network` type that lets us chain multiple layers.--> data Layer :: Nat -> Nat -> Type where->     Layer :: { _lWeights :: L m n->              , _lBiases  :: R m->              }->           -> Layer n m->       deriving (Show, Generic)-> ->-> data Network :: Nat -> [Nat] -> Nat -> Type where->     NØ   :: !(Layer a b) -> Network a '[] b->     (:&) :: !(Layer a b) -> Network b bs c -> Network a (b ': bs) c--A `Layer n m` is a layer taking an n-vector and returning an m-vector.  A-`Network a '[b, c, d] e` would be a Network that takes in an a-vector and-outputs an e-vector, with hidden layers of sizes b, c, and d.--Isomorphisms---------------The *backprop* library lets you apply operations on "parts" of data types-(like on the weights and biases of a `Layer`) by using `Iso`'s-(isomorphisms), like the ones from the *lens* library.  The library doesn't-depend on lens, but it can use the `Iso`s from the library and also-custom-defined ones.--First, we can auto-generate isomorphisms using the *generics-sop* library:--> instance SOP.Generic (Layer n m)--And then can create isomorphisms by hand for the two `Network`-constructors:--> netExternal :: Iso' (Network a '[] b) (Tuple '[Layer a b])-> netExternal = iso (\case NØ x     -> x ::< Ø)->                   (\case I x :< Ø -> NØ x   )-> -> netInternal :: Iso' (Network a (b ': bs) c) (Tuple '[Layer a b, Network b bs c])-> netInternal = iso (\case x :& xs          -> x ::< xs ::< Ø)->                   (\case I x :< I xs :< Ø -> x :& xs       )--An `Iso' a (Tuple as)` means that an `a` can really just be seen as a tuple-of `as`.--Running a network-=================--Now, we can write the `BPOp` that reprenents running the network and-getting a result.  We pass in a `Sing bs` (a singleton list of the hidden-layer sizes) so that we can "pattern match" on the list and handle the-different network constructors differently.--> netOp->     :: forall s a bs c. (KnownNat a, KnownNat c)->     => Sing bs->     -> BPOp s '[ R a, Network a bs c ] (R c)-> netOp sbs = go sbs->   where->     go :: forall d es. KnownNat d->         => Sing es->         -> BPOp s '[ R d, Network d es c ] (R c)->     go = \case->       SNil -> withInps $ \(x :< n :< Ø) -> do->         -- peek into the NØ using netExternal iso->         l :< Ø <- netExternal #<~ n->         -- run the 'layerOp' BP, with x and l as inputs->         bpOp layerOp ~$ (x :< l :< Ø)->       SNat `SCons` ses -> withInps $ \(x :< n :< Ø) -> withSingI ses $ do->         -- peek into the (:&) using the netInternal iso->         l :< n' :< Ø <- netInternal #<~ n->         -- run the 'layerOp' BP, with x and l as inputs->         z <- bpOp layerOp  ~$ (x :< l :< Ø)->         -- run the 'go ses' BP, with z and n as inputs->         bpOp (go ses)      ~$ (z :< n' :< Ø)->     layerOp->         :: forall d e. (KnownNat d, KnownNat e)->         => BPOp s '[ R d, Layer d e ] (R e)->     layerOp = withInps $ \(x :< l :< Ø) -> do->         -- peek into the layer using the gTuple iso, auto-generated with SOP.Generic->         w :< b :< Ø <- gTuple #<~ l->         y           <- matVec  ~$ (w :< x :< Ø)->         return $ logistic (y + b)--There's some singletons work going on here, but it's fairly standard-singletons stuff.  Most of the complexity here is from the static typing in-our neural network type, and *not* from *backprop*.--From *backprop* specifically, the only elements are `#<~` lets you "split" an-input ref with the given iso, and `bpOp`, which converts a `BPOp` into an `Op`-that you can bind with `~$`.--Note that this library doesn't support truly pattern matching on GADTs, and-that we had to pass in `Sing bs` as a reference to the structure of our-networks.--Gradient Descent-------------------Now we can do simple gradient descent.  Defining an error function:--> errOp->     :: KnownNat m->     => R m->     -> BVar s rs (R m)->     -> BPOp s rs Double-> errOp targ r = do->     err <- bindVar $ r - t->     dot ~$ (err :< err :< Ø)->   where->     t = constVar targ--And now, we can use `backprop` to generate the gradient, and shift the-`Network`!  Things are made a bit cleaner from the fact that `Network a bs c`-has a `Num` instance, so we can use `(-)` and `(*)` etc.--> train->     :: (KnownNat a, SingI bs, KnownNat c)->     => Double->     -> R a->     -> R c->     -> Network a bs c->     -> Network a bs c-> train r x t n = case backprop (errOp t =<< netOp sing) (x ::< n ::< Ø) of->     (_, _ :< I g :< Ø) -> n - (realToFrac r * g)--(`(::<)` is cons and `Ø` is nil for tuples.)--Main-====--`main`, which will train on sample data sets, is still in progress!  Right-now it just generates a random network using the *mwc-random* library and-prints each internal layer.--> main :: IO ()-> main = withSystemRandom $ \g -> do->     n <- uniform @(Network 4 '[3,2] 1) g->     void $ traverseNetwork sing (\l -> l <$ print l) n--Appendix: Boilerplate-=====================--And now for some typeclass instances and boilerplates unrelated to the-*backprop* library that makes our custom types easier to use.--> instance KnownNat n => Variate (R n) where->     uniform g = randomVector <$> uniform g <*> pure Uniform->     uniformR (l, h) g = (\x -> x * (h - l) + l) <$> uniform g-> -> instance (KnownNat m, KnownNat n) => Variate (L m n) where->     uniform g = uniformSample <$> uniform g <*> pure 0 <*> pure 1->     uniformR (l, h) g = (\x -> x * (h - l) + l) <$> uniform g-> -> instance (KnownNat n, KnownNat m) => Variate (Layer n m) where->     uniform g = subtract 1 . (* 2) <$> (Layer <$> uniform g <*> uniform g)->     uniformR (l, h) g = (\x -> x * (h - l) + l) <$> uniform g-> -> instance (KnownNat m, KnownNat n) => Num (Layer n m) where->     Layer w1 b1 + Layer w2 b2 = Layer (w1 + w2) (b1 + b2)->     Layer w1 b1 - Layer w2 b2 = Layer (w1 - w2) (b1 - b2)->     Layer w1 b1 * Layer w2 b2 = Layer (w1 * w2) (b1 * b2)->     abs    (Layer w b) = Layer (abs w) (abs b)->     signum (Layer w b) = Layer (signum w) (signum b)->     negate (Layer w b) = Layer (negate w) (negate b)->     fromInteger x = Layer (fromInteger x) (fromInteger x)-> -> instance (KnownNat m, KnownNat n) => Fractional (Layer n m) where->     Layer w1 b1 / Layer w2 b2 = Layer (w1 / w2) (b1 / b2)->     recip (Layer w b) = Layer (recip w) (recip b)->     fromRational x = Layer (fromRational x) (fromRational x)-> -> instance (KnownNat a, SingI bs, KnownNat c) => Variate (Network a bs c) where->     uniform g = genNet sing (uniform g)->     uniformR (l, h) g = (\x -> x * (h - l) + l) <$> uniform g-> -> genNet->     :: forall f a bs c. (Applicative f, KnownNat a, KnownNat c)->     => Sing bs->     -> (forall d e. (KnownNat d, KnownNat e) => f (Layer d e))->     -> f (Network a bs c)-> genNet sbs f = go sbs->   where->     go :: forall d es. KnownNat d => Sing es -> f (Network d es c)->     go = \case->       SNil             -> NØ <$> f->       SNat `SCons` ses -> (:&) <$> f <*> go ses-> -> mapNetwork0->     :: forall a bs c. (KnownNat a, KnownNat c)->     => Sing bs->     -> (forall d e. (KnownNat d, KnownNat e) => Layer d e)->     -> Network a bs c-> mapNetwork0 sbs f = getI $ genNet sbs (I f)-> -> traverseNetwork->     :: forall a bs c f. (KnownNat a, KnownNat c, Applicative f)->     => Sing bs->     -> (forall d e. (KnownNat d, KnownNat e) => Layer d e -> f (Layer d e))->     -> Network a bs c->     -> f (Network a bs c)-> traverseNetwork sbs f = go sbs->   where->     go :: forall d es. KnownNat d => Sing es -> Network d es c -> f (Network d es c)->     go = \case->       SNil -> \case->         NØ x -> NØ <$> f x->       SNat `SCons` ses -> \case->         x :& xs -> (:&) <$> f x <*> go ses xs-> -> mapNetwork1->     :: forall a bs c. (KnownNat a, KnownNat c)->     => Sing bs->     -> (forall d e. (KnownNat d, KnownNat e) => Layer d e -> Layer d e)->     -> Network a bs c->     -> Network a bs c-> mapNetwork1 sbs f = getI . traverseNetwork sbs (I . f)-> -> mapNetwork2->     :: forall a bs c. (KnownNat a, KnownNat c)->     => Sing bs->     -> (forall d e. (KnownNat d, KnownNat e) => Layer d e -> Layer d e -> Layer d e)->     -> Network a bs c->     -> Network a bs c->     -> Network a bs c-> mapNetwork2 sbs f = go sbs->   where->     go :: forall d es. KnownNat d => Sing es -> Network d es c -> Network d es c -> Network d es c->     go = \case->       SNil -> \case->         NØ x -> \case->           NØ y -> NØ (f x y)->       SNat `SCons` ses -> \case->         x :& xs -> \case->           y :& ys -> f x y :& go ses xs ys-> -> instance (KnownNat a, SingI bs, KnownNat c) => Num (Network a bs c) where->     (+)           = mapNetwork2 sing (+)->     (-)           = mapNetwork2 sing (-)->     (*)           = mapNetwork2 sing (*)->     negate        = mapNetwork1 sing negate->     abs           = mapNetwork1 sing abs->     signum        = mapNetwork1 sing signum->     fromInteger x = mapNetwork0 sing (fromInteger x)-> -> instance (KnownNat a, SingI bs, KnownNat c) => Fractional (Network a bs c) where->     (/)            = mapNetwork2 sing (/)->     recip          = mapNetwork1 sing recip->     fromRational x = mapNetwork0 sing (fromRational x)

@@ -0,0 +1,530 @@+% Extensible Neural Networks with Backprop+% Justin Le++This write-up is a follow-up to the *MNIST* tutorial+([rendered][mnist-rendered] here, and [literate haskell][mnist-lhs] here).+This write-up itself is available as a [literate haskell file][lhs], and+also [rendered as a pdf][rendered].++[mnist-rendered]: https://github.com/mstksg/backprop/blob/master/renders/backprop-mnist.pdf+[mnist-lhs]: https://github.com/mstksg/backprop/blob/master/samples/backprop-mnist.lhs+[rendered]: https://github.com/mstksg/backprop/blob/master/renders/extensible-neural.pdf+[lhs]: https://github.com/mstksg/backprop/blob/master/samples/extensible-neural.lhs++The packages involved are:++*   deepseq+*   hmatrix+*   lens+*   mnist-idx+*   mwc-random+*   one-liner-instances+*   singletons+*   split+*   vector++> {-# LANGUAGE BangPatterns         #-}+> {-# LANGUAGE DataKinds            #-}+> {-# LANGUAGE DeriveGeneric        #-}+> {-# LANGUAGE FlexibleContexts     #-}+> {-# LANGUAGE GADTs                #-}+> {-# LANGUAGE InstanceSigs         #-}+> {-# LANGUAGE LambdaCase           #-}+> {-# LANGUAGE LambdaCase           #-}+> {-# LANGUAGE RankNTypes           #-}+> {-# LANGUAGE ScopedTypeVariables  #-}+> {-# LANGUAGE TemplateHaskell      #-}+> {-# LANGUAGE TypeApplications     #-}+> {-# LANGUAGE TypeInType           #-}+> {-# LANGUAGE TypeOperators        #-}+> {-# LANGUAGE ViewPatterns         #-}+> {-# OPTIONS_GHC -fno-warn-orphans #-}+>+> import           Control.DeepSeq+> import           Control.Exception+> import           Control.Lens hiding             ((<.>))+> import           Control.Monad+> import           Control.Monad.IO.Class+> import           Control.Monad.Primitive+> import           Control.Monad.Trans.Maybe+> import           Control.Monad.Trans.State+> import           Data.Bitraversable+> import           Data.Foldable+> import           Data.IDX+> import           Data.Kind+> import           Data.List.Split+> import           Data.Singletons+> import           Data.Singletons.Prelude+> import           Data.Singletons.TypeLits+> import           Data.Time.Clock+> import           Data.Traversable+> import           Data.Tuple+> import           GHC.Generics                    (Generic)+> import           Numeric.Backprop+> import           Numeric.LinearAlgebra.Static+> import           Numeric.OneLiner+> import           Text.Printf+> import qualified Data.Vector                     as V+> import qualified Data.Vector.Generic             as VG+> import qualified Data.Vector.Unboxed             as VU+> import qualified Numeric.LinearAlgebra           as HM+> import qualified System.Random.MWC               as MWC+> import qualified System.Random.MWC.Distributions as MWC++Introduction+============++The *[backprop][hackage]* library lets us manipulate our values in a+natural way.  We write the function to compute our result, and the library+then automatically finds the *gradient* of that function, which we can use+for gradient descent.++[hackage]: http://hackage.haskell.org/package/backprop++In the last post, we looked at using a fixed-structure neural network.+However, in [this blog series][blog], I discuss a system of extensible+neural networks that can be chained and composed.++[blog]: https://blog.jle.im/entries/series/+practical-dependent-types-in-haskell.html++One issue, however, in naively translating the implementations, is that we+normally run the network by pattern matching on each layer.  However, we+cannot directly pattern match on `BVar`s.++We *could* get around it by being smart with prisms and `^^?`, to extract a+"Maybe BVar".  However, we can do better!  This is because the *shape* of a+`Net i hs o` is known already at compile-time, so there is no need for+runtime checks like prisms and `^^?`.++Instead, we can just directly use lenses, since we know *exactly* what+constructor will be present!  We can use singletons to determine which+constructor is present, and so always just directly use lenses without any+runtime nondeterminism.++Types+=====++First, our types:++> data Layer i o =+>     Layer { _lWeights :: !(L o i)+>           , _lBiases  :: !(R o)+>           }+>   deriving (Show, Generic)+>+> instance NFData (Layer i o)+> makeLenses ''Layer+>+> data Net :: Nat -> [Nat] -> Nat -> Type where+>     NO   :: !(Layer i o) -> Net i '[] o+>     (:~) :: !(Layer i h) -> !(Net h hs o) -> Net i (h ': hs) o++Unfortunately, we can't automatically generate lenses for GADTs, so we have+to make them by hand.[^poly]++[poly]: We write them originally as a polymorphic lens family to help us+with type safety via paraemtric polymorphism.++> _NO :: Lens (Net i '[] o) (Net i' '[] o')+>             (Layer i o  ) (Layer i' o'  )+> _NO f (NO l) = NO <$> f l+>+> _NIL :: Lens (Net i (h ': hs) o) (Net i' (h ': hs) o)+>              (Layer i h        ) (Layer i' h        )+> _NIL f (l :~ n) = (:~ n) <$> f l+>+> _NIN :: Lens (Net i (h ': hs) o) (Net i (h ': hs') o')+>              (Net h hs        o) (Net h hs'        o')+> _NIN f (l :~ n) = (l :~) <$> f n++You can read `_NO` as:++```haskell+_NO :: Lens' (Net i '[] o) (Layer i o)+```++A lens into a single-layer network, and++```haskell+_NIL :: Lens' (Net i (h ': hs) o) (Layer i h )+_NIN :: Lens' (Net i (h ': hs) o) (Net h hs o)+```++Lenses into a multiple-layer network, getting the first layer and the tail+of the network.++If we pattern match on `Sing hs`, we can always determine exactly which+lenses we can use, and so never fumble around with prisms or+nondeterminism.++Running the network+===================++Here's the meat of process, then: specifying how to run the network.  We+re-use our `BVar`-based combinators defined in the last write-up:++> runLayer+>     :: (KnownNat i, KnownNat o, Reifies s W)+>     => BVar s (Layer i o)+>     -> BVar s (R i)+>     -> BVar s (R o)+> runLayer l x = (l ^^. lWeights) #>! x + (l ^^. lBiases)+> {-# INLINE runLayer #-}++For `runNetwork`, we pattern match on `hs` using singletons, so we always+know exactly what type of network we have:++> runNetwork+>     :: (KnownNat i, KnownNat o, Reifies s W)+>     => BVar s (Net i hs o)+>     -> Sing hs+>     -> BVar s (R i)+>     -> BVar s (R o)+> runNetwork n = \case+>     SNil          -> softMax . runLayer (n ^^. _NO)+>     SCons SNat hs -> withSingI hs (runNetwork (n ^^. _NIN) hs)+>                    . logistic+>                    . runLayer (n ^^. _NIL)+> {-# INLINE runNetwork #-}++The rest of it is the same as before.++> netErr+>     :: (KnownNat i, KnownNat o, SingI hs, Reifies s W)+>     => R i+>     -> R o+>     -> BVar s (Net i hs o)+>     -> BVar s Double+> netErr x targ n = crossEntropy targ (runNetwork n sing (constVar x))+> {-# INLINE netErr #-}+>+> trainStep+>     :: forall i hs o. (KnownNat i, KnownNat o, SingI hs)+>     => Double             -- ^ learning rate+>     -> R i                -- ^ input+>     -> R o                -- ^ target+>     -> Net i hs o         -- ^ initial network+>     -> Net i hs o+> trainStep r !x !targ !n = n - realToFrac r * gradBP (netErr x targ) n+> {-# INLINE trainStep #-}+>+> trainList+>     :: (KnownNat i, SingI hs, KnownNat o)+>     => Double             -- ^ learning rate+>     -> [(R i, R o)]       -- ^ input and target pairs+>     -> Net i hs o         -- ^ initial network+>     -> Net i hs o+> trainList r = flip $ foldl' (\n (x,y) -> trainStep r x y n)+> {-# INLINE trainList #-}+>+> testNet+>     :: forall i hs o. (KnownNat i, KnownNat o, SingI hs)+>     => [(R i, R o)]+>     -> Net i hs o+>     -> Double+> testNet xs n = sum (map (uncurry test) xs) / fromIntegral (length xs)+>   where+>     test :: R i -> R o -> Double          -- test if the max index is correct+>     test x (extract->t)+>         | HM.maxIndex t == HM.maxIndex (extract r) = 1+>         | otherwise                                = 0+>       where+>         r :: R o+>         r = evalBP (\n' -> runNetwork n' sing (constVar x)) n++And that's it!++Running+=======++Everything here is the same as before, except now we can dynamically pick+the network size.  Here we pick `'[300,100]` for the hidden layer sizes.++> main :: IO ()+> main = MWC.withSystemRandom $ \g -> do+>     Just train <- loadMNIST "data/train-images-idx3-ubyte" "data/train-labels-idx1-ubyte"+>     Just test  <- loadMNIST "data/t10k-images-idx3-ubyte"  "data/t10k-labels-idx1-ubyte"+>     putStrLn "Loaded data."+>     net0 <- MWC.uniformR @(Net 784 '[300,100] 10) (-0.5, 0.5) g+>     flip evalStateT net0 . forM_ [1..] $ \e -> do+>       train' <- liftIO . fmap V.toList $ MWC.uniformShuffle (V.fromList train) g+>       liftIO $ printf "[Epoch %d]\n" (e :: Int)+>+>       forM_ ([1..] `zip` chunksOf batch train') $ \(b, chnk) -> StateT $ \n0 -> do+>         printf "(Batch %d)\n" (b :: Int)+>+>         t0 <- getCurrentTime+>         n' <- evaluate . force $ trainList rate chnk n0+>         t1 <- getCurrentTime+>         printf "Trained on %d points in %s.\n" batch (show (t1 `diffUTCTime` t0))+>+>         let trainScore = testNet chnk n'+>             testScore  = testNet test n'+>         printf "Training error:   %.2f%%\n" ((1 - trainScore) * 100)+>         printf "Validation error: %.2f%%\n" ((1 - testScore ) * 100)+>+>         return ((), n')+>   where+>     rate  = 0.02+>     batch = 5000++Looking Forward+===============++One common thing people might do is want to be able to mix different types+of layers.  This could also be easily encoded as different constructors in+`Layer`, and so `runLayer` will now be different depending on what+constructor is present.++In this case, we can either:++1.  Have a different indexed type for layers, so that we can always know+    exactly what layer is involved, so we don't have to runtime pattern+    match:++    ```haskell+    data LayerType = FullyConnected | Convolutional++    data Layer :: LayerType -> Nat -> Nat -> Type where+        LayerFC :: .... -> Layer 'FullyConnected i o+        LayerC  :: .... -> Layer 'Convolutional  i o+    ```++    We would then have `runLayer` take `Sing (t :: LayerType)`, so we can+    again use `^^.` and directly pattern match.++2.  Use a typeclass-based approach, so users can add their own layer types.+    In this situation, layer types would all be different types, and+    running them would be a typeclass method that would give our+    `BVar s (Layer i o) -> BVar s (R i) -> BVar s (R o)` operation as a+    typeclass method.++    ```haskell+    class Layer (l :: Nat -> Nat -> Type) where+        runLayer+            :: forall s. Reifies s W+            => BVar s (l i o)+            -> BVar s (R i)+            -> BVar s (R o)+    ```++In all cases, it shouldn't be much more cognitive overhead to use+*backprop* to build your neural network framework!++And, remember that `evalBP` (directly running the function) introduces+virtually zero overhead, so if you only provided `BVar` functions, you+could easily get the original non-`BVar` functions with `evalBP` without+any loss.++What now?+---------++Ready to start?  Check out the docs for the [Numeric.Backprop][] module for+the full technical specs, and find more examples and updates at the [github+repo][repo]!++[Numeric.Backprop]: http://hackage.haskell.org/package/backprop/docs/Numeric-Backprop.html+[repo]: https://github.com/mstksg/backprop++Internals+=========++That's it for the post!  Now for the internal plumbing :)++> loadMNIST+>     :: FilePath+>     -> FilePath+>     -> IO (Maybe [(R 784, R 10)])+> loadMNIST fpI fpL = runMaybeT $ do+>     i <- MaybeT          $ decodeIDXFile       fpI+>     l <- MaybeT          $ decodeIDXLabelsFile fpL+>     d <- MaybeT . return $ labeledIntData l i+>     r <- MaybeT . return $ for d (bitraverse mkImage mkLabel . swap)+>     liftIO . evaluate $ force r+>   where+>     mkImage :: VU.Vector Int -> Maybe (R 784)+>     mkImage = create . VG.convert . VG.map (\i -> fromIntegral i / 255)+>     mkLabel :: Int -> Maybe (R 10)+>     mkLabel n = create $ HM.build 10 (\i -> if round i == n then 1 else 0)++HMatrix Operations+------------------++> infixr 8 #>!+> (#>!)+>     :: (KnownNat m, KnownNat n, Reifies s W)+>     => BVar s (L m n)+>     -> BVar s (R n)+>     -> BVar s (R m)+> (#>!) = liftOp2 . op2 $ \m v ->+>   ( m #> v, \g -> (g `outer` v, tr m #> g) )+>+> infixr 8 <.>!+> (<.>!)+>     :: (KnownNat n, Reifies s W)+>     => BVar s (R n)+>     -> BVar s (R n)+>     -> BVar s Double+> (<.>!) = liftOp2 . op2 $ \x y ->+>   ( x <.> y, \g -> (konst g * y, x * konst g)+>   )+>+> konst'+>     :: (KnownNat n, Reifies s W)+>     => BVar s Double+>     -> BVar s (R n)+> konst' = liftOp1 . op1 $ \c -> (konst c, HM.sumElements . extract)+>+> sumElements'+>     :: (KnownNat n, Reifies s W)+>     => BVar s (R n)+>     -> BVar s Double+> sumElements' = liftOp1 . op1 $ \x -> (HM.sumElements (extract x), konst)+>+> softMax :: (KnownNat n, Reifies s W) => BVar s (R n) -> BVar s (R n)+> softMax x = konst' (1 / sumElements' expx) * expx+>   where+>     expx = exp x+> {-# INLINE softMax #-}+>+> crossEntropy+>     :: (KnownNat n, Reifies s W)+>     => R n+>     -> BVar s (R n)+>     -> BVar s Double+> crossEntropy targ res = -(log res <.>! constVar targ)+> {-# INLINE crossEntropy #-}+>+> logistic :: Floating a => a -> a+> logistic x = 1 / (1 + exp (-x))+> {-# INLINE logistic #-}++Instances+---------++> instance (KnownNat i, KnownNat o) => Num (Layer i o) where+>     (+)         = gPlus+>     (-)         = gMinus+>     (*)         = gTimes+>     negate      = gNegate+>     abs         = gAbs+>     signum      = gSignum+>     fromInteger = gFromInteger+>+> instance (KnownNat i, KnownNat o) => Fractional (Layer i o) where+>     (/)          = gDivide+>     recip        = gRecip+>     fromRational = gFromRational+>+>+> liftNet0+>     :: forall i hs o. (KnownNat i, KnownNat o)+>     => (forall m n. (KnownNat m, KnownNat n) => Layer m n)+>     -> Sing hs+>     -> Net i hs o+> liftNet0 x = go+>   where+>     go :: forall w ws. KnownNat w => Sing ws -> Net w ws o+>     go = \case+>       SNil          -> NO x+>       SCons SNat hs -> x :~ go hs+>+> liftNet1+>     :: forall i hs o. (KnownNat i, KnownNat o)+>     => (forall m n. (KnownNat m, KnownNat n)+>           => Layer m n+>           -> Layer m n+>        )+>     -> Sing hs+>     -> Net i hs o+>     -> Net i hs o+> liftNet1 f = go+>   where+>     go  :: forall w ws. KnownNat w+>         => Sing ws+>         -> Net w ws o+>         -> Net w ws o+>     go = \case+>       SNil          -> \case+>         NO x -> NO (f x)+>       SCons SNat hs -> \case+>         x :~ xs -> f x :~ go hs xs+>+> liftNet2+>     :: forall i hs o. (KnownNat i, KnownNat o)+>     => (forall m n. (KnownNat m, KnownNat n)+>           => Layer m n+>           -> Layer m n+>           -> Layer m n+>        )+>     -> Sing hs+>     -> Net i hs o+>     -> Net i hs o+>     -> Net i hs o+> liftNet2 f = go+>   where+>     go  :: forall w ws. KnownNat w+>         => Sing ws+>         -> Net w ws o+>         -> Net w ws o+>         -> Net w ws o+>     go = \case+>       SNil          -> \case+>         NO x -> \case+>           NO y -> NO (f x y)+>       SCons SNat hs -> \case+>         x :~ xs -> \case+>           y :~ ys -> f x y :~ go hs xs ys+>+> instance ( KnownNat i+>          , KnownNat o+>          , SingI hs+>          )+>       => Num (Net i hs o) where+>     (+)           = liftNet2 (+) sing+>     (-)           = liftNet2 (-) sing+>     (*)           = liftNet2 (*) sing+>     negate        = liftNet1 negate sing+>     abs           = liftNet1 abs sing+>     signum        = liftNet1 signum sing+>     fromInteger x = liftNet0 (fromInteger x) sing+>+> instance ( KnownNat i+>          , KnownNat o+>          , SingI hs+>          )+>       => Fractional (Net i hs o) where+>     (/)            = liftNet2 (/) sing+>     recip          = liftNet1 negate sing+>     fromRational x = liftNet0 (fromRational x) sing+>+> instance KnownNat n => MWC.Variate (R n) where+>     uniform g = randomVector <$> MWC.uniform g <*> pure Uniform+>     uniformR (l, h) g = (\x -> x * (h - l) + l) <$> MWC.uniform g+>+> instance (KnownNat m, KnownNat n) => MWC.Variate (L m n) where+>     uniform g = uniformSample <$> MWC.uniform g <*> pure 0 <*> pure 1+>     uniformR (l, h) g = (\x -> x * (h - l) + l) <$> MWC.uniform g+>+> instance (KnownNat i, KnownNat o) => MWC.Variate (Layer i o) where+>     uniform g = Layer <$> MWC.uniform g <*> MWC.uniform g+>     uniformR (l, h) g = (\x -> x * (h - l) + l) <$> MWC.uniform g+>+> instance ( KnownNat i+>          , KnownNat o+>          , SingI hs+>          )+>       => MWC.Variate (Net i hs o) where+>     uniform :: forall m. PrimMonad m => MWC.Gen (PrimState m) -> m (Net i hs o)+>     uniform g = go sing+>       where+>         go :: forall w ws. KnownNat w => Sing ws -> m (Net w ws o)+>         go = \case+>           SNil          -> NO <$> MWC.uniform g+>           SCons SNat hs -> (:~) <$> MWC.uniform g <*> go hs+>     uniformR (l, h) g = (\x -> x * (h - l) + l) <$> MWC.uniform g+>+> instance NFData (Net i hs o) where+>     rnf = \case+>       NO l    -> rnf l+>       x :~ xs -> rnf x `seq` rnf xs

@@ -1,58 +1,30 @@-{-# LANGUAGE AllowAmbiguousTypes    #-}-{-# LANGUAGE ConstraintKinds        #-} {-# LANGUAGE DataKinds              #-}-{-# LANGUAGE EmptyCase              #-}-{-# LANGUAGE FlexibleContexts       #-}-{-# LANGUAGE KindSignatures         #-} {-# LANGUAGE LambdaCase             #-} {-# LANGUAGE PolyKinds              #-} {-# LANGUAGE RankNTypes             #-} {-# LANGUAGE ScopedTypeVariables    #-}-{-# LANGUAGE TypeApplications       #-}-{-# LANGUAGE TypeFamilies           #-} {-# LANGUAGE TypeFamilyDependencies #-}-{-# LANGUAGE TypeInType             #-} {-# LANGUAGE TypeOperators          #-}  module Data.Type.Util (     Replicate   , unzipP   , zipP-  , tagSum-  , indexP+  , zipWithPM_   , vecToProd+  , vecLen   , prodToVec'-  , prodAlong   , lengthProd-  , prodLength-  , vecLength-  , finIndex-  , replLen-  , replWit-  , itraverse1_-  , ifor1-  , ifor1_-  , for1-  , for1_+  , listToVec+  , fillProd   ) where -import           Control.Applicative import           Data.Bifunctor--- import           Data.Kind-import           Data.Monoid hiding    (Sum) import           Data.Type.Conjunction-import           Data.Type.Fin-import           Data.Type.Index import           Data.Type.Length import           Data.Type.Nat import           Data.Type.Product-import           Data.Type.Sum import           Data.Type.Vector-import           Lens.Micro-import           Type.Class.Higher--- import           Type.Class.Known-import           Type.Class.Witness--- import           Type.Family.List import           Type.Family.Nat  -- | @'Replicate' n a@ is a list of @a@s repeated @n@ times.@@ -72,6 +44,13 @@     ØV      -> Ø     x :* xs -> x :< vecToProd xs +vecLen+    :: VecT n f a+    -> Nat n+vecLen = \case+    ØV      -> Z_+    _ :* xs -> S_ (vecLen xs)+ prodToVec'     :: Nat n     -> Prod f (Replicate n a)@@ -82,71 +61,22 @@     S_ n -> \case       x :< xs -> x :* prodToVec' n xs -prodAlong-    :: VecT n f b-    -> Prod f (Replicate n a)-    -> VecT n f a-prodAlong = \case-    ØV -> \case-      Ø       -> ØV-    _ :* v -> \case-      x :< xs -> x :* prodAlong v xs--finIndex-    :: Fin n-    -> Index (Replicate n a) a-finIndex = \case-    FZ   -> IZ-    FS f -> IS (finIndex f)--traverse1_-    :: (Applicative h, Traversable1 t)-    => (forall a. f a -> h ())-    -> t f b+zipWithPM_+    :: forall h f g as. Applicative h+    => (forall a. f a -> g a -> h ())+    -> Prod f as+    -> Prod g as     -> h ()-traverse1_ f = ($ pure ())-             . appEndo-             . getConst-             . foldMap1 (\y -> Const (Endo (f y *>)))+zipWithPM_ f = go+  where+    go :: forall bs. Prod f bs -> Prod g bs -> h ()+    go = \case+      Ø -> \case+        Ø -> pure ()+      x :< xs -> \case+        y :< ys -> f x y *> go xs ys -itraverse1_-    :: (Applicative h, IxFoldable1 i t)-    => (forall a. i b a -> f a -> h ())-    -> t f b-    -> h ()-itraverse1_ f = ($ pure ())-              . appEndo-              . getConst-              . ifoldMap1 (\i y -> Const (Endo (f i y *>))) -for1-    :: (Applicative h, Traversable1 t)-    => t f b-    -> (forall a. f a -> h (g a))-    -> h (t g b)-for1 x f = traverse1 f x--for1_-    :: (Applicative h, Traversable1 t)-    => t f b-    -> (forall a. f a -> h ())-    -> h ()-for1_ x f = traverse1_ f x--ifor1-    :: (Applicative h, IxTraversable1 i t)-    => t f b-    -> (forall a. i b a -> f a -> h (g a))-    -> h (t g b)-ifor1 x f = itraverse1 f x--ifor1_-    :: (Applicative h, IxFoldable1 i t)-    => t f b-    -> (forall a. i b a -> f a -> h ())-    -> h ()-ifor1_ x f = itraverse1_ f x- zipP     :: Prod f as     -> Prod g as@@ -156,6 +86,7 @@       Ø       -> Ø     x :< xs -> \case       y :< ys -> x :&: y :< zipP xs ys+{-# INLINE zipP #-}  unzipP     :: Prod (f :&: g) as@@ -164,58 +95,6 @@     Ø               -> (Ø, Ø)     (x :&: y) :< zs -> bimap (x :<) (y :<) (unzipP zs) -indexP :: Index as a -> Lens' (Prod g as) (g a)-indexP = \case-    IZ   -> \f -> \case-      x :< xs -> (:< xs) <$> f x-    IS i -> \f -> \case-      x :< xs -> (x :<) <$> indexP i f xs--prodLength-    :: Prod f as-    -> Length as-prodLength = \case-    Ø       -> LZ-    _ :< xs -> LS (prodLength xs)--vecLength-    :: forall n f a. ()-    => VecT n f a-    -> Nat n-vecLength = \case-    ØV      -> Z_-    _ :* xs -> S_ (vecLength xs)---- | Currently not used-tagSum-    :: Prod f as-    -> Sum g as-    -> Sum (f :&: g) as-tagSum = \case-    Ø       -> \case-    x :< xs -> \case-      InL y  -> InL (x :&: y)-      InR ys -> InR (tagSum xs ys)--replWit-    :: Nat n-    -> Wit (c a)-    -> Wit (Every c (Replicate n a))-replWit = \case-    Z_   -> \case-      Wit -> Wit-    S_ n -> \case-      c@Wit -> case replWit n c of-        Wit -> Wit--replLen-    :: forall n a. ()-    => Nat n-    -> Length (Replicate n a)-replLen = \case-    Z_   -> LZ-    S_ n -> LS (replLen @_ @a n)- lengthProd     :: (forall a. f a)     -> Length as@@ -223,3 +102,28 @@ lengthProd x = \case     LZ   -> Ø     LS l -> x :< lengthProd x l++listToVec+    :: Nat n+    -> [f a]+    -> Maybe (VecT n f a)+listToVec = \case+    Z_ -> \_ -> Just ØV+    S_ n -> \case+      []   -> Nothing+      x:xs -> (x :*) <$> listToVec n xs++fillProd+    :: forall f g as c. ()+    => (forall a. f a -> c -> g a)+    -> Prod f as+    -> [c]+    -> Maybe (Prod g as)+fillProd f = go+  where+    go :: Prod f bs -> [c] -> Maybe (Prod g bs)+    go = \case+      Ø -> \_ -> Just Ø+      x :< xs -> \case+        []   -> Nothing+        y:ys -> (f x y :<) <$> go xs ys

@@ -1,1480 +1,394 @@-{-# LANGUAGE FlexibleContexts    #-}-{-# LANGUAGE GADTs               #-}-{-# LANGUAGE LambdaCase          #-}-{-# LANGUAGE PatternSynonyms     #-}-{-# LANGUAGE RankNTypes          #-}-{-# LANGUAGE RecordWildCards     #-}-{-# LANGUAGE ScopedTypeVariables #-}-{-# LANGUAGE TypeApplications    #-}-{-# LANGUAGE TypeInType          #-}-{-# LANGUAGE TypeOperators       #-}---- |--- Module      : Numeric.Backprop--- Copyright   : (c) Justin Le 2017--- License     : BSD3------ Maintainer  : justin@jle.im--- Stability   : experimental--- Portability : non-portable--------- Provides the 'BP' monad and the 'BVar' type; after manipulating 'BVar's--- (inputs to your function) to produce a result, the library tracks internal data--- dependencies, which are used to perform back-propagation (reverse-mode--- automatic differentiation) to calculate the gradient of the output with--- respect to the inputs.------ Similar to automatic differentiation from the /ad/ library and--- "Numeric.AD.Mode.Reverse", except for a few key differences:------ 1. Most importantly, this library implements /heterogeneous/--- back-propagation, so you can manipulate values of different types (like--- different matrix and vector types, and product and sum types).  This is--- essential for things like back-propagation for neural networks.------ 2. This module allows you to /explicitly/ build your data dependency--- graph if you wish, which allows the library to perform optimizations and--- reduce extra allocation, which may or may not provide advantages over--- "Numeric.AD.Mode.Reverse"'s 'System.IO.Unsafe.unsafePerformIO'-based--- implicit graph building.------ See the <https://github.com/mstksg/backprop README> for more information--- and links to demonstrations and tutorials.  If you want to plunge right--- in, you can also look directly at the main types, 'BP', 'BPOp', 'BVar',--- 'Op', and the main functions, 'backprop' and 'opVar'.------ Note that every type involved has to be an instance of 'Num'.  This is--- because gradients all need to be "summable" (which is implemented using--- 'sum' and '+'), and we also need to able to generate gradients of '1'--- and '0'.-----module Numeric.Backprop (-  -- * Types-  -- ** Backprop types-    BP, BPOp, BPOpI, BVar, Op, OpB-  -- ** Tuple types#prod#-  -- $prod-  , Prod(..), Tuple, I(..)-  -- * BP-  -- ** Backprop-  , backprop, evalBPOp, gradBPOp-  -- ** Utility combinators-  , withInps, implicitly-  , withInps', implicitly'-  -- * Vars-  , constVar-  , inpVar, inpVars-  , bpOp-  , bindVar-  , inpVars'-  -- ** From Ops-  , opVar, (~$)-  , opVar1, opVar2, opVar3-  , (-$)-  -- ** Var manipulation-  -- *** As parts-  , partsVar, (#<~), withParts-  , splitVars, gSplit, gTuple-  -- *** As sums-  , choicesVar, (?<~), withChoices-  -- $sum-  , Sum(..)-  -- *** As sums of products-  , sopVar, gSplits, gSOP-  -- *** As GADTs-  , withGADT, BPCont(..)-  -- ** Combining-  , liftB, (.$), liftB1, liftB2, liftB3-  -- * Op-  , op1, op2, op3, opN, composeOp, composeOp1, (~.)-  , op1', op2', op3'-  -- * Utility-  , pattern (:>), only, head'-  , pattern (::<), only_-  -- ** Numeric Ops-  -- | Optimized ops for numeric functions.  See-  -- "Numeric.Backprop.Op#numops" for more information.-  , (+.), (-.), (*.), negateOp, absOp, signumOp-  , (/.), recipOp-  , expOp, logOp, sqrtOp, (**.), logBaseOp-  , sinOp, cosOp, tanOp, asinOp, acosOp, atanOp-  , sinhOp, coshOp, tanhOp, asinhOp, acoshOp, atanhOp-  ) where--import           Control.Monad.Base-import           Control.Monad.Reader-import           Control.Monad.ST-import           Control.Monad.State-import           Data.Kind-import           Data.Maybe-import           Data.Monoid               ((<>))-import           Data.STRef-import           Data.Type.Combinator-import           Data.Type.Conjunction-import           Data.Type.Index-import           Data.Type.Length-import           Data.Type.Product-import           Data.Type.Sum hiding      (index)-import           Data.Type.Util-import           Lens.Micro hiding         (ix)-import           Lens.Micro.Mtl hiding     (view)-import           Numeric.Backprop.Internal-import           Numeric.Backprop.Iso-import           Numeric.Backprop.Op-import           Type.Class.Higher-import           Type.Class.Known-import           Type.Class.Witness-import qualified Generics.SOP              as SOP---- $prod------ 'Prod', from the <http://hackage.haskell.org/package/type-combinators--- type-combinators> library (in "Data.Type.Prod") is a heterogeneous--- list/tuple type, which allows you to tuple together multiple values of--- different types and operate on them generically.------ A @'Prod' f '[a, b, c]@ contains an @f a@, an @f b@, and an @f c@, and--- is constructed by consing them together with ':<' (using 'Ø' as nil):------ @--- 'I' "hello" ':<' I True :< I 7.8 :< Ø    :: 'Prod' 'I' '[String, Bool, Double]--- 'C' "hello" :< C "world" :< C "ok" :< Ø  :: 'Prod' ('C' String) '[a, b, c]--- 'Proxy' :< Proxy :< Proxy :< Ø           :: 'Prod' 'Proxy' '[a, b, c]--- @------ ('I' is the identity functor, and 'C' is the constant functor)------ So, in general:------ @--- x :: f a--- y :: f b--- z :: f c--- x :< y :< z :< Ø :: Prod f '[a, b, c]--- @------ If you're having problems typing 'Ø', you can use 'only':------ @--- only z           :: Prod f '[c]--- x :< y :< only z :: Prod f '[a, b, c]--- @------ 'Tuple' is provided as a convenient type synonym for 'Prod' 'I', and has--- a convenient pattern synonym '::<' (and 'only_'), which can also be used--- for pattern matching:------ @--- x :: a--- y :: b--- z :: c------ 'only_' z             :: 'Tuple' '[c]--- x '::<' y ::< z ::< Ø :: 'Tuple' '[a, b, c]--- x ::< y ::< only_ z :: 'Tuple' '[a, b, c]--- @----- $sum------ #sum#------ Like the 'Prod' type (see mini-tutorial at "Numeric.Backprop#prod"), the--- 'Sum' type (from the--- <http://hackage.haskell.org/package/type-combinators type-combinators>--- library, in "Data.Type.Sum") lets you make arbitrary sum types over--- different types and work with them generically.------ A @'Sum' f '[a, b, c]@ contains /either/ an @f a@, an @f b@, /or/ an @f--- c@, and is constructed with the constructors 'InL' and 'InR', which are--- analogous to 'Left' and 'Right'.------ For a value of type @'Sum' f '[Int, Bool, String]@, there are three--- constructors:------ @--- 'InL'             :: f Int    -> 'Sum' f '[Int, Bool, String]--- InL . InR       :: f Bool   -> Sum f '[Int, Bool, String]--- InL . InR . InR :: f String -> Sum f '[Int, Bool, String]--- @------ Each 'InR' "pushes deeper" into the 'Sum'.------ Likewise, if you have a value of type @'Sum' f '[Int, Bool, String]@,--- you can see which constructor it was made (and what type it contains)--- with by pattern matching:------ @--- foo :: 'Sum' f '[Int, Bool, String]------ case foo of---   'InL' i         -> -- foo contains an "f Int"---   'InR' (InL b)   -> -- foo contains an "f Bool"---   InR (InR (InL s)) -> -- foo contains an "f String"--- @------ | A handy type synonym representing a 'BP' action that returns a 'BVar'.--- This is handy because this is the form of 'BP' actions that--- 'backprop' and 'gradBPOp' (etc.) expects.------ A value of type:------ @--- 'BPOp' s rs a--- @------ is an action that takes an input environment of @rs@ and produces--- a 'BVar' containing a value of type @a@.  Because it returns a 'BVar',--- the library can track the data dependencies between the 'BVar' and the--- input environment and perform back-propagation.------ See documentation for 'BP' for an explanation of the phantom type--- parameter @s@.-type BPOp s rs a  = BP s rs (BVar s rs a)---- | An "implicit" operation on 'BVar's that can be backpropagated.--- A value of type:------ @--- 'BPOpI' s rs a--- @------ takes a bunch of 'BVar's containg @rs@ and uses them to (purely) produce--- a 'BVar' containing an @a@.------ @--- foo :: BPOpI s '[ Double, Double ] Double--- foo (x :< y :< Ø) = x + sqrt y--- @------ If you are exclusively doing implicit back-propagation by combining--- 'BVar's and using 'BPOpI's, you are probably better off just importing--- "Numeric.Backprop.Implicit", which provides better tools.  This type--- synonym exists in "Numeric.Backprop" just for the 'implicitly' function,--- which can convert "implicit" backprop functions like a @'BPOpI' s rs a@--- into an "explicit" graph backprop function, a @'BPOp' s rs a@.-type BPOpI s rs a = Prod (BVar s rs) rs -> BVar s rs a----- | Apply an 'OpB' to a 'Prod' (tupling) of 'BVar's.------ If you had an @'OpB' s '[a, b, c] d@, this function will expect a 3-Prod--- of a @'BVar' s rs a@, a @'BVar' s rs b@, and a @'BVar' s rs c@, and the--- result will be a @'BVar' s rs d@:------ @--- myOp :: 'OpB' s '[a, b, c] d--- x    :: 'BVar' s rs a--- y    :: 'BVar' s rs b--- z    :: 'BVar' s rs c------ x :< y :< z :< Ø              :: 'Prod' ('BVar' s rs) '[a, b, c]--- 'opVar' myOp (x :< y :< z :< Ø) :: 'BP' s rs ('BVar' s rs d)--- @------ Note that 'OpB' is a superclass of 'Op', so you can provide any 'Op'--- here, as well (like those created by 'op1', 'op2', 'constOp', 'op0'--- etc.)------ 'opVar' has an infix alias, '~$', so the above example can also be--- written as:------ @--- myOp '~$' (x :< y :< z :< Ø) :: 'BP' s rs ('BVar' s rs d)--- @------ to let you pretend that you're applying the 'myOp' function to three--- inputs.------ Also note the relation between 'opVar' and 'liftB' and 'bindVar':------ @--- 'opVar' o xs = 'bindVar' ('liftB' o xs)--- @------ 'opVar' can be thought of as a "binding" version of 'liftB'.-opVar-    :: forall s rs as a. Num a-    => OpB s as a-    -> Prod (BVar s rs) as-    -> BP s rs (BVar s rs a)-opVar o i = do-    xs <- traverse1 (fmap I . BP . resolveVar) i-    (res, gf) <- BP . liftBase $ runOpM' o xs-    let bp = BPN { _bpnOut       = only $ FRInternal []-                 , _bpnRes       = only_ res-                 , _bpnGradFunc  = gf . head'-                 , _bpnGradCache = Nothing-                 }-    r <- BP . liftBase $ newSTRef bp-    itraverse1_ (registerVar . flip IRNode r) i-    return (BVNode IZ r)---- | Split out a 'BVar' of a tuple into a tuple ('Prod') of 'BVar's.------ @--- -- the environment is a single Int-Bool tuple, tup--- stuff :: 'BP' s '[ Tuple '[Int, Bool] ] a--- stuff = 'withInps' $ \\(tup :< Ø) -\> do---     i :< b :< Ø <- 'splitVars' tup---     -- now, i is a 'BVar' pointing to the 'Int' inside tup---     -- and b is a 'BVar' pointing to the 'Bool' inside tup---     -- you can do stuff with the i and b here--- @------ Note that------ @--- 'splitVars' = 'partsVar' 'id'--- @-splitVars-    :: forall s rs as. Every Num as-    => BVar s rs (Tuple as)-    -> BP s rs (Prod (BVar s rs) as)-splitVars = partsVar id---- | Use an 'Iso' (or compatible 'Control.Lens.Iso.Iso' from the lens--- library) to "pull out" the parts of a data type and work with each part--- as a 'BVar'.------ If there is an isomorphism between a @b@ and a @'Tuple' as@ (that is, if--- an @a@ is just a container for a bunch of @as@), then it lets you break--- out the @as@ inside and work with those.------ @--- data Foo = F Int Bool------ fooIso :: 'Iso'' Foo (Tuple '[Int, Bool])--- fooIso = 'iso' (\\(F i b)         -\> i ::\< b ::\< Ø)---              (\\(i ::\< b ::\< Ø) -\> F i b        )------ 'partsVar' fooIso :: 'BVar' rs Foo -> 'BP' s rs ('Prod' ('BVar' s rs) '[Int, Bool])------ stuff :: 'BP' s '[Foo] a--- stuff = 'withInps' $ \\(foo :< Ø) -\> do---     i :< b :< Ø <- partsVar fooIso foo---     -- now, i is a 'BVar' pointing to the 'Int' inside foo---     -- and b is a 'BVar' pointing to the 'Bool' inside foo---     -- you can do stuff with the i and b here--- @------ You can use this to pass in product types as the environment to a 'BP',--- and then break out the type into its constituent products.------ Note that for a type like @Foo@, @fooIso@ can be generated automatically--- with 'GHC.Generics.Generic' from "GHC.Generics" and--- 'Generics.SOP.Generic' from "Generics.SOP" and /generics-sop/, using the--- 'gTuple' iso.  See 'gSplit' for more information.------ Also, if you are literally passing a tuple (like--- @'BP' s '[Tuple '[Int, Bool]@) then you can give in the identity--- isomorphism ('id') or use 'splitVars'.-partsVar-    :: forall s rs bs b. Every Num bs-    => Iso' b (Tuple bs)-    -> BVar s rs b-    -> BP s rs (Prod (BVar s rs) bs)-partsVar i = fmap (view sum1) . sopVar (i . resum1)---- | A useful infix alias for 'partsVar'.------ Building on the example from 'partsVar':------ @--- data Foo = F Int Bool------ fooIso :: 'Iso'' Foo (Tuple '[Int, Bool])--- fooIso = 'iso' (\\(F i b)         -\> i ::\< b ::\< Ø)---              (\\(i ::\< b ::\< Ø) -\> F i b        )------ stuff :: 'BP' s '[Foo] a--- stuff = 'withInps' $ \\(foo :< Ø) -\> do---     i :< b :< Ø <- fooIso '#<~' foo---     -- now, i is a 'BVar' pointing to the 'Int' inside foo---     -- and b is a 'BVar' pointing to the 'Bool' inside foo---     -- you can do stuff with the i and b here--- @------ See 'gSplit' for an example usage of splitting up an arbitrary product--- type (like @Foo@) using "GHC.Geneics" and "Generics.SOP".-infixr 1 #<~-(#<~)-    :: (Every Num bs, Known Length bs)-    => Iso' b (Tuple bs)-    -> BVar s rs b-    -> BP s rs (Prod (BVar s rs) bs)-(#<~) = partsVar---- | A continuation-based version of 'partsVar'.  Instead of binding the--- parts and using it in the rest of the block, provide a continuation to--- handle do stuff with the parts inside.------ Building on the example from 'partsVar':------ @--- data Foo = F Int Bool------ fooIso :: 'Iso'' Foo (Tuple '[Int, Bool])--- fooIso = 'iso' (\\(F i b)         -\> i ::\< b ::\< Ø)---              (\\(i ::\< b ::\< Ø) -\> F i b        )------ stuff :: 'BP' s '[Foo] a--- stuff = 'withInps' $ \\(foo :< Ø) -\> do---     'withParts' fooIso foo $ \\(i :< b :< Ø) -\> do---       -- now, i is a 'BVar' pointing to the 'Int' inside foo---       -- and b is a 'BVar' pointing to the 'Bool' inside foo---       -- you can do stuff with the i and b here--- @------ Useful so that you can work with the internal parts of the data type--- in a closure, so the parts don't leak out to the rest of your 'BP'.--- But, mostly just a stylistic choice.-withParts-    :: Every Num bs-    => Iso' b (Tuple bs)-    -> BVar s rs b-    -> (Prod (BVar s rs) bs -> BP s rs a)-    -> BP s rs a-withParts i r f = do-    p <- partsVar i r-    f p---- | Using 'GHC.Generics.Generic' from "GHC.Generics" and--- 'Generics.SOP.Generic' from "Generics.SOP", /split/ a 'BVar' containing--- a product type into a tuple ('Prod') of 'BVar's pointing to each value--- inside.------ Building on the example from 'partsVar':------ @--- import qualified Generics.SOP as SOP------ data Foo = F Int Bool---   deriving Generic------ instance SOP.Generic Foo------ 'gSplit' :: 'BVar' rs Foo -> 'BP' s rs ('Prod' ('BVar' s rs) '[Int, Bool])------ stuff :: 'BP' s '[Foo] a--- stuff = 'withInps' $ \\(foo :< Ø) -\> do---     i :< b :< Ø <- 'gSplit' foo---     -- now, i is a 'BVar' pointing to the 'Int' inside foo---     -- and b is a 'BVar' pointing to the 'Bool' inside foo---     -- you can do stuff with the i and b here--- @------ Because @Foo@ is a straight up product type, 'gSplit' can use--- "GHC.Generics" and take out the items inside.------ Note that because------ @--- 'gSplit' = 'splitVars' 'gTuple'--- @------ Then, you can also use 'gTuple' with '#<~':------ @--- stuff :: 'BP' s '[Foo] a--- stuff = 'withInps' $ \\(foo :< Ø) -\> do---     i :< b :< Ø <- 'gTuple' '#<~' foo---     -- now, i is a 'BVar' pointing to the 'Int' inside foo---     -- and b is a 'BVar' pointing to the 'Bool' inside foo---     -- you can do stuff with the i and b here--- @----gSplit-    :: (Every Num bs, SOP.Generic b, SOP.Code b ~ '[bs])-    => BVar s rs b-    -> BP s rs (Prod (BVar s rs) bs)-gSplit = partsVar gTuple---- | Use an 'Iso' (or compatible 'Control.Lens.Iso.Iso' from the lens--- library) to "pull out" the different constructors of a sum type and--- return a (choice) sum of 'BVar's that you can pattern match on.------ If there is an isomorphism between a @b@ and a @'Sum' 'I' as@ (that is,--- if an @a@ is just a sum type for every type in @as@), then it lets you--- /branch/ on which constructor is used inside the @b@.------ Essentially implements pattern matching on 'BVar' values.------ @--- data Bar = A Int | B Bool | C String------ barIso :: 'Iso'' Bar ('Sum' I '[Int, Bool, String])--- barIso = 'iso' (\\case A i -> 'InL' (I i)---                       B b -> 'InR' ('InL' (I b))---                       C s -> 'InR' ('InR' ('InL' (I s))---                )---                (\\case 'InL' (I i)           -> A i---                       'InR' ('InL' (I b))       -> B b---                       'InR' ('InR' ('InL' (I s))) -> C s---                )------ choicesVar barIso :: BVar rs Bar -> BP s rs (Sum I (BVar s rs) '[Int, Bool, String])------ stuff :: 'BP' s '[Bar] a--- stuff = 'withInps' $ \\(bar :< Ø) -\> do---     c <- 'choicesVar' barIso bar---     case c of---       'InL' i -> do---          -- in this branch, bar was made with the A constructor---          -- i is the Int inside it---       'InR' ('InL' b) -> do---          -- in this branch, bar was made with the B constructor---          -- b is the Bool inside it---       'InR' ('InR' ('InL' s)) -> do---          -- in this branch, bar was made with the B constructor---          -- s is the String inside it--- @------ You can use this to pass in sum types as the environment to a 'BP', and--- then branch on which constructor the value was made with.------ See "Numeric.Backprop#sum" for a mini-tutorial on 'Sum'.-choicesVar-    :: forall s rs bs b. Every Num bs-    => Iso' b (Sum I bs)-    -> BVar s rs b-    -> BP s rs (Sum (BVar s rs) bs)-choicesVar i r = do-    x <- BP $ resolveVar r-    let xs :: Sum I bs-        xs = view i x-    ifor1 xs $ \ix (I (y :: c)) -> every @_ @Num ix // do-      let bp :: BPNode s rs '[b] '[c]-          bp = BPN { _bpnOut       = only $ FRInternal []-                   , _bpnRes       = only_ y-                   , _bpnGradFunc  = return . only_ . review i-                                   . injectSum ix-                                   -- . maybe (I (1 \\ every @_ @Num ix)) I-                                   . maybe (I 1) I-                                   . head'-                   , _bpnGradCache = Nothing-                   }-      r' <- BP . liftBase $ newSTRef bp-      registerVar (IRNode IZ r') r-      return $ BVNode IZ r'--- TODO: cannot implement via sopVar?  oh well.---- | A continuation-based version of 'choicesVar'.  Instead of binding the--- parts and using it in the rest of the block, provide a continuation that--- will handle every possible constructor/case of the type of the value the--- 'BVar' points to.------ Building on the example from 'choicesVar':------ @--- data Bar = A Int | B Bool | C String------ barIso :: 'Iso'' Bar ('Sum' I '[Int, Bool, String])--- barIso = 'iso' (\\case A i -> 'InL' (I i)---                       B b -> 'InR' ('InL' (I b))---                       C s -> 'InR' ('InR' ('InL' (I s))---                )---                (\\case 'InL' (I i)           -> A i---                       'InR' ('InL' (I b))       -> B b---                       'InR' ('InR' ('InL' (I s))) -> C s---                )------ 'choicesVar' barIso :: BVar rs Bar -> BP s rs (Sum I (BVar s rs) '[Int, Bool, String])------ stuff :: 'BP' s '[Bar] a--- stuff = 'withInps' $ \\(bar :< Ø) -\> do---     'withChoices' barIso bar $ \case---       'InL' i -> do---          -- in this branch, bar was made with the A constructor---          -- i is the Int inside it---       'InR' ('InL' b) -> do---          -- in this branch, bar was made with the B constructor---          -- b is the Bool inside it---       'InR' ('InR' ('InL' s)) -> do---          -- in this branch, bar was made with the B constructor---          -- s is the String inside it--- @------ Nicer than 'choicesVar' directly, because you don't have to give the--- result a superfluous name before pattern matching on it.  You can just--- directly pattern match in the lambda, so there's a lot less syntactical--- noise.-withChoices-    :: forall s rs bs b a. Every Num bs-    => Iso' b (Sum I bs)-    -> BVar s rs b-    -> (Sum (BVar s rs) bs -> BP s rs a)-    -> BP s rs a-withChoices i r f = do-    c <- choicesVar i r-    f c---- | A useful infix alias for 'choicesVar'.------ Building on the example from 'choicesVar':------ @--- data Bar = A Int | B Bool | C String------ barIso :: 'Iso'' Bar ('Sum' I '[Int, Bool, String])--- barIso = 'iso' (\\case A i -> 'InL' (I i)---                       B b -> 'InR' ('InL' (I b))---                       C s -> 'InR' ('InR' ('InL' (I s))---                )---                (\\case 'InL' (I i)           -> A i---                       'InR' ('InL' (I b))       -> B b---                       'InR' ('InR' ('InL' (I s))) -> C s---                )------ stuff :: 'BP' s '[Bar] a--- stuff = 'withInps' $ \\(bar :< Ø) -\> do---     c <- barIso '?<~' bar---     case c of---       'InL' i -> do---          -- in this branch, bar was made with the A constructor---          -- i is the Int inside it---       'InR' ('InL' b) -> do---          -- in this branch, bar was made with the B constructor---          -- b is the Bool inside it---       'InR' ('InR' ('InL' s)) -> do---          -- in this branch, bar was made with the B constructor---          -- s is the String inside it--- @-infixr 1 ?<~-(?<~)-    :: (Every Num bs, Known Length bs)-    => Iso' b (Sum I bs)-    -> BVar s rs b-    -> BP s rs (Sum (BVar s rs) bs)-(?<~) = choicesVar---- | A combination of 'partsVar' and 'choicesVar', that lets you split--- a type into a sum of products.  Using an 'Iso' (or compatible--- 'Control.Lens.Iso.Iso' from the lens library), you can pull out a type--- that is a sum of products into a sum of products of 'BVar's.------ Implements branching on the constructors of a value that a 'BVar'--- contains, and also splitting out the different items inside each--- constructor.------ @--- data Baz = A Int    Bool---          | B String Double--------- bazIso :: 'Iso'' Baz ('Sum' 'Tuple' '[ '[Int, Bool], '[String, Double] ])--- bazIso = 'iso' (\\case A i b -> 'InL' (I (i ::< b ::< Ø))---                       B s d -> 'InR' ('InL' (I (s ::< d ::< Ø)))---                )---                (\\case 'InL' (I (i ::< b ::< Ø))     -> A i b---                       'InR' ('InL' (I (s ::< d ::< Ø))) -> B s d---                )------ 'sopVar' bazIso :: 'BVar' rs Baz -> 'BP' s rs ('Sum' ('Prod' ('BVar' s rs)) '[ '[Int, Bool], '[String, Double] ])------ stuff :: 'BP' s '[Baz] a--- stuff = 'withInps' $ \\(baz :< Ø) -\> do---     c <- 'sopVar' barIso baz---     case c of---       'InL' (i :< b :< Ø) -> do---          -- in this branch, baz was made with the A constructor---          -- i and b are the Int and Bool inside it---       'InR' ('InL' (s :< d :< Ø)) -> do---          -- in this branch, baz was made with the B constructor---          -- s and d are the String and Double inside it--- @------ Essentially exists to implement "pattern matching" on multiple--- constructors and fields for the value inside a 'BVar'.------ Note that for a type like @Baz@, @bazIso@ can be generated automatically--- with 'GHC.Generics.Generic' from "GHC.Generics" and--- 'Generics.SOP.Generic' from "Generics.SOP" and /generics-sop/, with--- 'gSOP'.  See 'gSplits' for more information.------ See "Numeric.Backprop#sum" for a mini-tutorial on 'Sum'.-sopVar-    :: forall s rs bss b. Every (Every Num) bss-    => Iso' b (Sum Tuple bss)-    -> BVar s rs b-    -> BP s rs (Sum (Prod (BVar s rs)) bss)-sopVar i r = do-    x <- BP $ resolveVar r-    let xs :: Sum Tuple bss-        xs = view i x-    ifor1 xs $ \ix (ys :: Tuple bs) -> every @_ @(Every Num) ix // do-      let bp :: BPNode s rs '[b] bs-          bp = BPN { _bpnOut       = map1 (const (FRInternal [])) ys-                   , _bpnRes       = ys-                   , _bpnGradFunc  = return . only_-                                   . review i . injectSum ix-                                   . imap1 (\ix' -> every @_ @Num ix' //-                                                      maybe (I 1) I-                                           )-                   , _bpnGradCache = Nothing-                   }-      r' <- BP . liftBase $ newSTRef bp-      registerVar (IRNode IZ r') r-      return $ imap1 (\ix' _ -> BVNode ix' r') ys---- | Using 'GHC.Generics.Generic' from "GHC.Generics" and--- 'Generics.SOP.Generic' from "Generics.SOP", /split/ a 'BVar' containing--- a sum of products (any simple ADT, essentialy) into a 'Sum' of each--- constructor, each containing a tuple ('Prod') of 'BVar's pointing to--- each value inside.------ Building on the example from 'sopVar':------ @--- import qualified Generics.SOP as SOP------ data Baz = A Int    Bool---          | B String Double---   deriving Generic------ instance SOP.Generic Baz------ 'gSplits' :: 'BVar' rs Baz -> 'BP' s rs ('Sum' ('Prod' ('BVar' s rs)) '[ '[Int, Bool], '[String, Double] ])------ stuff :: 'BP' s '[Baz] a--- stuff = 'withInps' $ \\(baz :< Ø) -\> do---     c <- gSplits baz---     case c of---       'InL' (i :< b :< Ø) -> do---          -- in this branch, baz was made with the A constructor---          -- i and b are the Int and Bool inside it---       'InR' ('InL' (s :< d :< Ø)) -> do---          -- in this branch, baz was made with the B constructor---          -- s and d are the String and Double inside it--- @------ Because @Foo@ is a straight up sum-of-products type, 'gSplits' can use--- "GHC.Generics" and take out the items inside.------ Note:------ @--- 'gSplit' = 'splitVars' 'gSOP'--- @------ See "Numeric.Backprop#sum" for a mini-tutorial on 'Sum'.-gSplits-    :: forall s rs b. (SOP.Generic b, Every (Every Num) (SOP.Code b))-    => BVar s rs b-    -> BP s rs (Sum (Prod (BVar s rs)) (SOP.Code b))-gSplits = sopVar gSOP--resolveVar-    :: (MonadReader (Tuple rs) m, MonadBase (ST s) m)-    => BVar s rs a-    -> m a-resolveVar = \case-    BVNode  ix r -> getI . index ix . _bpnRes <$> liftBase (readSTRef r)-    BVInp   ix   -> getI . index ix <$> ask-    BVConst    x -> return x-    BVOp    rs o -> do-      xs <- traverse1 (fmap I . resolveVar) rs-      liftBase $ runOpM o xs--registerVar-    :: forall s rs a. ()-    => BPInpRef s rs a-    -> BVar s rs a-    -> BP s rs ()-registerVar bpir = \case-    BVNode  ix' r' -> BP . liftBase . modifySTRef r' $-                        over (bpnOut . indexP ix' . _FRInternal) (bpir :)-    BVInp   ix'    -> BP $ modifying (bpsSources . indexP ix' . _FRInternal) (bpir :)-    BVConst _      -> return ()-    -- This independently makes a new BPPipe for every usage site of the-    -- BVOp, so it's a bit inefficient.-    BVOp    (rs :: Prod (BVar s rs) ds) (o :: OpM (ST s) ds a) -> do-      xs :: Tuple ds <- traverse1 (fmap I . BP . resolveVar) rs-      (res, gF) <- BP . liftBase $ runOpM' o xs-      let bpp :: BPPipe s rs ds '[a]-          bpp = BPP { _bppOut       = only bpir-                    , _bppRes       = only_ res-                    , _bppGradFunc  = gF . Just . getI . head'-                    , _bppGradCache = Nothing-                    }-      r' <- BP . liftBase $ newSTRef bpp-      ifor1_ rs $ \ix' (bpr :: BVar s rs d) ->-        registerVar (IRPipe ix' r') bpr---- | Infix synonym for 'opVar', which lets you pretend that you're applying--- 'OpB's as if they were functions:------ @--- myOp :: 'OpB' s '[a, b, c] d--- x    :: 'BVar' s rs a--- y    :: 'BVar' s rs b--- z    :: 'BVar' s rs c------ x :< y :< z :< Ø           :: 'Prod' ('BVar' s rs) '[a, b, c]--- myOp '~$' (x :< y :< z :< Ø) :: 'BP' s rs ('BVar' s rs d)--- @------ Note that 'OpB' is a superclass of 'Op', so you can pass in any 'Op'--- here, as well (like those created by 'op1', 'op2', 'constOp', 'op0'--- etc.)------ '~$' can also be thought of as a "binding" version of '.$':------ @--- o '~$' xs = 'bindVar' (o '.$' xs)--- @----infixr 5 ~$-(~$)-    :: Num a-    => OpB s as a-    -> Prod (BVar s rs) as-    -> BP s rs (BVar s rs a)-(~$) = opVar---- | Lets you treat a @'BPOp' s as b@ as an @'Op' as b@, and "apply"--- arguments to it just like you would with an 'Op' and '~$' / 'opVar'.------ Basically a convenient wrapper over 'bpOp' and '~$':------ @--- o '-$' xs = bpOp o '~$' xs--- @------ So for a @'BPOp' s as b@, you can "plug in" 'BVar's to @as@, and get--- a @b@ as a result.------ Useful for running a @'BPOp' s as b@ that you got from a different function, and--- "plugging in" its @as@ inputs with 'BVar's from your current--- environment.-infixr 5 -$-(-$)-    :: (Every Num as, Known Length as, Num a)-    => BPOp s as a-    -> Prod (BVar s rs) as-    -> BPOp s rs a-o -$ xs = bpOp o ~$ xs---- | Create a 'BVar' that represents just a specific value, that doesn't--- depend on any other 'BVar's.-constVar :: a -> BVar s rs a-constVar = BVConst---- | Convenient wrapper over 'opVar' that takes an 'OpB' with one argument--- and a single 'BVar' argument.  Lets you not have to type out the entire--- 'Prod'.------ @--- 'opVar1' o x = 'opVar' o (x ':<' 'Ø')------ myOp :: 'Op' '[a] b--- x    :: 'BVar' s rs a------ 'opVar1' myOp x :: 'BP' s rs ('BVar' s rs b)--- @------ Note that 'OpB' is a superclass of 'Op', so you can pass in an 'Op' here--- (like one made with 'op1') as well.-opVar1-    :: Num b-    => OpB s '[a] b-    -> BVar s rs a-    -> BP s rs (BVar s rs b)-opVar1 o = opVar o . only---- | Convenient wrapper over 'opVar' that takes an 'OpB' with two arguments--- and two 'BVar' arguments.  Lets you not have to type out the entire--- 'Prod'.------ @--- 'opVar2' o x y = 'opVar' o (x ':<' y ':<' 'Ø')------ myOp :: 'Op' '[a, b] c--- x    :: 'BVar' s rs a--- y    :: 'BVar' s rs b------ 'opVar2' myOp x y :: 'BP' s rs ('BVar' s rs c)--- @------ Note that 'OpB' is a superclass of 'Op', so you can pass in an 'Op' here--- (like one made with 'op2') as well.-opVar2-    :: Num c-    => OpB s '[a,b] c-    -> BVar s rs a-    -> BVar s rs b-    -> BP s rs (BVar s rs c)-opVar2 o rx ry = opVar o (rx :< ry :< Ø)---- | Convenient wrapper over 'opVar' that takes an 'OpB' with three arguments--- and three 'BVar' arguments.  Lets you not have to type out the entire--- 'Prod'.------ @--- 'opVar3' o x y z = 'opVar' o (x ':<' y ':<' z ':<' 'Ø')------ myOp :: 'Op' '[a, b, c] d--- x    :: 'BVar' s rs a--- y    :: 'BVar' s rs b--- z    :: 'BVar' s rs c------ 'opVar3' myOp x y z :: 'BP' s rs ('BVar' s rs d)--- @------ Note that 'OpB' is a superclass of 'Op', so you can pass in an 'Op' here--- (like one made with 'op3') as well.-opVar3-    :: Num d-    => OpB s '[a,b,c] d-    -> BVar s rs a-    -> BVar s rs b-    -> BVar s rs c-    -> BP s rs (BVar s rs d)-opVar3 o rx ry rz = opVar o (rx :< ry :< rz :< Ø)---- | Concretizes a delayed 'BVar'.  If you build up a 'BVar' using numeric--- functions like '+' or '*' or using 'liftB', it'll defer the evaluation,--- and all of its usage sites will create a separate graph node.------ Use 'bindVar' if you ever intend to use a 'BVar' in more than one--- location.------ @--- -- bad--- errSquared :: Num a => 'BP' s '[a, a] a--- errSquared = 'withInp' $ \\(r :< t :< Ø) -\> do---     let err = r - t---     'return' (err * err)   -- err is used twice!------ -- good--- errSquared :: Num a => 'BP' s '[a, a] a--- errSquared = 'withInps' $ \\(r :< t :< Ø) -\> do---     let err = r - t---     e <- 'bindVar' err     -- force e, so that it's safe to use twice!---     'return' (e * e)------ -- better--- errSquared :: Num a => 'BP' s '[a, a] a--- errSquared = 'withInps' $ \\(r :< t :< Ø) -\> do---     let err = r - t---     e <- 'bindVar' err---     'bindVar' (e * e)      -- result is forced so user doesn't have to worry--- @------ Note the relation to 'opVar' / '~$' / 'liftB' / '.$':------ @--- 'opVar' o xs    = 'bindVar' ('liftB' o xs)--- o '~$' xs       = 'bindVar' (o '.$' xs)--- 'op2' (*) '~$' (x :< y :< Ø) = 'bindVar' (x * y)--- @------ So you can avoid 'bindVar' altogether if you use the explicitly binding--- '~$' and 'opVar' etc.------ Note that 'bindVar' on 'BVar's that are already forced is a no-op.-bindVar-    :: Num a-    => BVar s rs a-    -> BP s rs (BVar s rs a)-bindVar r = case r of-    BVNode  _  _ -> return r-    BVInp   _    -> return r-    BVConst _    -> return r-    BVOp    rs o -> opVar o rs----backwardPass-    :: forall s rs a. ()-    => BPInpRef s rs a-    -> ST s a-backwardPass = \case-    IRNode  ix r' -> getI . index ix <$> pullNode r'-    IRPipe  ix r' -> getI . index ix <$> pullPipe r'-    IRConst g     -> return g-  where-    pullNode-        :: forall as bs. Every Num bs-        => STRef s (BPNode s rs as bs)-        -> ST s (Tuple as)-    pullNode r = caching bpnGradCache r $ \BPN{..} -> do-        totdervs <- ifor1 _bpnOut $ \ix -> every @_ @Num ix // \case-          FRInternal rs -> Just . sum-            <$> traverse backwardPass rs-          FRTerminal g   -> return g-        g <- _bpnGradFunc totdervs-        return g-    pullPipe-        :: forall as bs. ()-        => STRef s (BPPipe s rs as bs)-        -> ST s (Tuple as)-    pullPipe r = caching bppGradCache r $ \BPP{..} ->-        _bppGradFunc =<< traverse1 (fmap I . backwardPass) _bppOut---- | Perform back-propagation on the given 'BPOp'.  Returns the result of--- the operation it represents, as well as the gradient of the result with--- respect to its inputs.  See module header for "Numeric.Backprop" and--- package documentation for examples and usages.-backprop-    :: Every Num rs-    => (forall s. BPOp s rs a)-    -> Tuple rs-    -> (a, Tuple rs)-backprop bp env = runST $ do-    (res, gFunc) <- backpropWith bp env-    grad <- gFunc Nothing-    return (res, grad)---- | Turn a 'BPOp' into an 'OpB'.  Basically converts a 'BP' taking an @rs@--- and producing an @a@ into an 'Op' taking an @rs@ and returning an @a@,--- with all of the powers and utility of an 'Op', including all of its--- gradient-finding glory.------ Really just reveals the fact that any @'BPOp' s rs a@ is itself an 'Op',--- an @'OpB' s rs a@, which makes it a differentiable function.------ Handy because an 'OpB' can be used with almost all of--- the 'Op'-related functions in this moduel, including 'opVar', '~$', etc.-bpOp-    :: Every Num rs-    => BPOp s rs a-    -> OpB s rs a-bpOp bp = OpM $ backpropWith bp---- | Simply run the 'BPOp' on an input tuple, getting the result without--- bothering with the gradient or with back-propagation.-evalBPOp-    :: (forall s. BPOp s rs a)  -- ^ 'BPOp' to run-    -> Tuple rs                 -- ^ input-    -> a                        -- ^ output-evalBPOp bp env = runST $ do-    r <- evalStateT (runReaderT (bpST bp) env)-                    (BPS (map1 (\_ -> FRInternal []) env))-    runReaderT (resolveVar r) env---- | Run the 'BPOp' on an input tuple and return the gradient of the result--- with respect to the input tuple.-gradBPOp-    :: Every Num rs-    => (forall s. BPOp s rs a)  -- ^ 'BPOp' to differentiate'-    -> Tuple rs                 -- ^ input-    -> Tuple rs                 -- ^ gradient-gradBPOp bp = snd . backprop bp--closeOff-    :: (MonadReader (Tuple rs) m, MonadState (BPState s rs) m, MonadBase (ST s) m)-    => Bool-    -> Maybe a-    -> BVar s rs a-    -> m ()-closeOff isTerminal gOut = \case-    BVNode  ix sr -> liftBase $ modifySTRef sr (over (bpnOut . indexP ix) (<> fr))-    BVInp   ix'   -> modifying (bpsSources . indexP ix') (<> fr)-    BVConst _     -> return ()-    BVOp    rs o  -> do-      xs <- traverse1 (fmap I . resolveVar) rs-      gs <- liftBase $ gradOpWithM' o xs gOut-      for1_ (gs `zipP` rs) $ \(I g :&: r) ->-        closeOff False (Just g) r-  where-    fr | isTerminal = FRTerminal gOut-       | otherwise  = FRInternal (IRConst <$> maybeToList gOut)---- | WARNING: the gradient continuation must only be run ONCE!-backpropWith-    :: Every Num rs-    => BPOp s rs a-    -> Tuple rs-    -> ST s (a, Maybe a -> ST s (Tuple rs))-backpropWith bp env = do-    (r, bps0) <- runStateT (runReaderT (bpST bp) env)-                           (BPS (map1 (\_ -> FRInternal []) env))-    res <- runReaderT (resolveVar r) env-    let gradFunc gradOut = do-          BPS{..} <- execStateT (runReaderT (closeOff True gradOut r) env) bps0-          ifor1 _bpsSources $ \ix rs -> every @_ @Num ix // do-            I <$> case rs of-              FRInternal rs' -> sum <$> traverse backwardPass rs'-              FRTerminal g   -> return $ fromMaybe 1 g-    return (res, gradFunc)---- | A version of 'implicitly' taking explicit 'Length', indicating the--- number of inputs required and their types.------ Requiring an explicit 'Length' is mostly useful for rare "extremely--- polymorphic" situations, where GHC can't infer the type and length of--- the list of inputs.  If you ever actually explicitly write down @rs@ as--- a list of types, you should be able to just use 'implicitly'.-implicitly'-    :: Num a-    => Length rs-    -> BPOpI s rs a-    -> BPOp s rs a-implicitly' l f = withInps' l (bindVar . f)---- | Convert a 'BPOpI' into a 'BPOp'.  That is, convert a function on--- a bundle of 'BVar's (generating an implicit graph) into a fully fledged--- 'BPOp' that you can run 'backprop' on.  See 'BPOpI' for more--- information.------ If you are going to write exclusively using implicit 'BVar' operations,--- it might be more convenient to use "Numeric.Backprop.Implicit" instead,--- which is geared around that use case.-implicitly-    :: (Known Length rs, Num a)-    => BPOpI s rs a-    -> BPOp s rs a-implicitly = implicitly' known---- | Create a 'BVar' given an index into the input environment.  For an--- example,------ @--- 'inpVar' 'IZ'--- @------ would refer to the /first/ input variable (the 'Int' in a--- @'BP' s '[Int, Bool]@), and------ @--- 'inpVar' ('IS' 'IZ')--- @------ Would refer to the /second/ input variable (the 'Bool' in a--- @'BP' s '[Int, Bool]@)------ Typically, there shouldn't be any reason to use 'inpVar' directly.  It's--- cleaner to get all of your input 'BVar's together using 'withInps' or--- 'inpVars'.-inpVar-    :: Index rs a-    -> BVar s rs a-inpVar = BVInp---- | Get a 'Prod' (tupling) of 'BVar's for all of the input environment--- (@rs@) of the @'BP' s rs@------ For example, if your 'BP' has an 'Int' and 'Double' in its input--- environment (a @'BP' s '[Int, Double]@), this would return a 'BVar'--- pointing to the 'Int' and a 'BVar' pointing to the 'Double'.------ @--- case ('inpVars' :: 'Prod' ('BVar' s '[Int, Double]) '[Int, Double]) of---   x :\< y :\< Ø -\> do---     -- the first item, x, is a var to the input 'Int'---     -- x :: 'BVar' s '[Int, Double] Int---     -- the second item, y, is a var to the input 'Double'---     -- y :: 'BVar' s '[Int, Double] Double--- @-inpVars-    :: Known Length rs-    => Prod (BVar s rs) rs-inpVars = inpVars' known---- | A version of 'inpVars' taking explicit 'Length', indicating the--- number of inputs required and their types.------ Mostly useful for rare "extremely polymorphic" situations, where GHC--- can't infer the type and length of the list of inputs.  If you ever--- actually explicitly write down @rs@ as a list of types, you should be--- able to just use 'inpVars'.-inpVars'-    :: Length rs-    -> Prod (BVar s rs) rs-inpVars' = map1 inpVar . indices'---- | A version of 'withInps' taking explicit 'Length', indicating the--- number of inputs required and their types.------ Mostly useful for rare "extremely polymorphic" situations, where GHC--- can't infer the type and length of the list of inputs.  If you ever--- actually explicitly write down @rs@ as a list of types, you should be--- able to just use 'withInps'.-withInps'-    :: Length rs-    -> (Prod (BVar s rs) rs -> BP s rs a)-    -> BP s rs a-withInps' l f = f (inpVars' l)---- | Runs a continuation on a 'Prod' of all of the input 'BVar's.------ Handy for bringing the environment into scope and doing stuff with it:------ @--- foo :: 'BPOp' '[Double, Int] a--- foo = 'withInps' $ \\(x :< y :< Ø) -\> do---     -- do stuff with inputs--- @------ Looks kinda like @foo (x :< y :< Ø) = -- ...@, don't it?------ Note that the above is the same as------ @--- foo :: 'BPOp' '[Double, Int] a--- foo = do---     case 'inpVars' of---       x :< y :< Ø -> do---         -- do stuff with inputs--- @------ But just a little nicer!-withInps-    :: Known Length rs-    => (Prod (BVar s rs) rs -> BP s rs a)-    -> BP s rs a-withInps = withInps' known---- | Apply 'OpB' over a 'Prod' of 'BVar's, as inputs. Provides--- "implicit-graph" back-propagation, with deferred evaluation.------ If you had an @'OpB' s '[a, b, c] d@, this function will expect a 3-Prod--- of a @'BVar' s rs a@, a @'BVar' s rs b@, and a @'BVar' s rs c@, and the--- result will be a @'BVar' s rs d@:------ @--- myOp :: 'OpB' s '[a, b, c] d--- x    :: 'BVar' s rs a--- y    :: 'BVar' s rs b--- z    :: 'BVar' s rs c------ x :< y :< z :< Ø              :: 'Prod' ('BVar' s rs) '[a, b, c]--- 'liftB' myOp (x :< y :< z :< Ø) :: 'BVar' s rs d--- @------ Note that 'OpB' is a superclass of 'Op', so you can provide any 'Op'--- here, as well (like those created by 'op1', 'op2', 'constOp', 'op0'--- etc.)------ 'liftB' has an infix alias, '.$', so the above example can also be--- written as:------ @--- myOp '.$' (x :< y :< z :< Ø) :: 'BVar' s rs d--- @------ to let you pretend that you're applying the 'myOp' function to three--- inputs.------ The result is a new /deferred/ 'BVar'.  This should be fine in most--- cases, unless you use the result in more than one location.  This will--- cause evaluation to be duplicated and multiple redundant graph nodes to--- be created.  If you need to use it in two locations, you should use--- 'opVar' instead of 'liftB', or use 'bindVar':------ @--- 'opVar' o xs = 'bindVar' ('liftB' o xs)--- @------ 'liftB' can be thought of as a "deferred evaluation" version of 'opVar'.-liftB-    :: OpB s as a-    -> Prod (BVar s rs) as-    -> BVar s rs a-liftB = flip BVOp----- | Infix synonym for 'liftB', which lets you pretend that you're applying--- 'OpB's as if they were functions:------ @--- myOp :: 'OpB' s '[a, b, c] d--- x    :: 'BVar' s rs a--- y    :: 'BVar' s rs b--- z    :: 'BVar' s rs c------ x :< y :< z :< Ø           :: 'Prod' ('BVar' s rs) '[a, b, c]--- myOp '.$' (x :< y :< z :< Ø) :: 'BVar' s rs d--- @------ Note that 'OpB' is a superclass of 'Op', so you can pass in any 'Op'--- here, as well (like those created by 'op1', 'op2', 'constOp', 'op0'--- etc.)------ See the documentation for 'liftB' for all the caveats of this usage.------ '.$' can also be thought of as a "deferred evaluation" version of '~$':------ @--- o '~$' xs = 'bindVar' (o '.$' xs)--- @----infixr 5 .$-(.$)-    :: OpB s as a-    -> Prod (BVar s rs) as-    -> BVar s rs a-(.$) = liftB----- | Convenient wrapper over 'liftB' that takes an 'OpB' with one argument--- and a single 'BVar' argument.  Lets you not have to type out the entire--- 'Prod'.------ @--- 'liftB1' o x = 'liftB' o (x ':<' 'Ø')------ myOp :: 'Op' '[a] b--- x    :: 'BVar' s rs a------ 'liftB1' myOp x :: 'BVar' s rs b--- @------ Note that 'OpB' is a superclass of 'Op', so you can pass in an 'Op' here--- (like one made with 'op1') as well.------ See the documentation for 'liftB' for caveats and potential problematic--- situations with this.-liftB1-    :: OpB s '[a] b-    -> BVar s rs a-    -> BVar s rs b-liftB1 o = liftB o . only---- | Convenient wrapper over 'liftB' that takes an 'OpB' with two arguments--- and two 'BVar' arguments.  Lets you not have to type out the entire--- 'Prod'.------ @--- 'liftB2' o x y = 'liftB' o (x ':<' y ':<' 'Ø')------ myOp :: 'Op' '[a, b] c--- x    :: 'BVar' s rs a--- y    :: 'BVar' s rs b------ 'liftB2' myOp x y :: 'BVar' s rs c--- @------ Note that 'OpB' is a superclass of 'Op', so you can pass in an 'Op' here--- (like one made with 'op2') as well.------ See the documentation for 'liftB' for caveats and potential problematic--- situations with this.-liftB2-    :: OpB s '[a,b] c-    -> BVar s rs a-    -> BVar s rs b-    -> BVar s rs c-liftB2 o x y = liftB o (x :< y :< Ø)---- | Convenient wrapper over 'liftB' that takes an 'OpB' with three arguments--- and three 'BVar' arguments.  Lets you not have to type out the entire--- 'Prod'.------ @--- 'liftB3' o x y z = 'liftB' o (x ':<' y ':<' z ':<' 'Ø')------ myOp :: 'Op' '[a, b, c] d--- x    :: 'BVar' s rs a--- y    :: 'BVar' s rs b--- z    :: 'BVar' s rs c------ 'liftB3' myOp x y z :: 'BVar' s rs d--- @------ Note that 'OpB' is a superclass of 'Op', so you can pass in an 'Op' here--- (like one made with 'op3') as well.------ See the documentation for 'liftB' for caveats and potential problematic--- situations with this.-liftB3-    :: OpB s '[a,b,c] d-    -> BVar s rs a-    -> BVar s rs b-    -> BVar s rs c-    -> BVar s rs d-liftB3 o x y z = liftB o (x :< y :< z :< Ø)---- | For usage with 'withGADT', to handle constructors of a GADT.  See--- documentation for 'withGADT' for more information.-data BPCont :: Type -> [Type] -> Type -> Type -> Type where-    BPC :: Every Num as-        => Tuple as-        -> (Tuple as -> a)-        -> (Prod (BVar s rs) as -> BP s rs b)-        -> BPCont s rs a b---- | Special __unsafe__ combinator that lets you pattern match and work on--- GADTs.------ @--- data MyGADT :: Bool -> Type where---     A :: String -> Int    -> MyGADT 'True---     B :: Bool   -> Double -> MyGADT 'False--------- foo :: BP s '[ MyGADT b ] a--- foo = 'withInps' $ \\( gVar :< Ø ) -\>---     withGADT gVar $ \\case---       A s i -\> BPC (s ::< i ::< Ø) (\\(s' ::< i' ::< Ø) -\> A s i) $---         \\(sVar :< iVar) -> do---           -- .. in this 'BP' action, sVar and iVar are 'BPVar's that---           -- refer to the String and Int inside the A constructor in---           -- gVar---       B b d -\> BPC (b ::< d ::< Ø) (\\(b' ::< d' ::< Ø) -\> B b d) $---         \\(bVar :< dVar) -> do---           -- .. in this 'BP' action, bVar and dVar are 'BPVar's that---           -- refer to the Bool and DOuble inside the B constructor in---           -- gVar--- @------ 'withGADT' lets to directly pattern match on the GADT, but as soon as--- you pattern match, you must handle the results with a 'BPCont'--- containing:------ 1.   /All/ of the items inside the GADT constructor, in a 'Tuple'--- 2.   A function from a 'Tuple' of items inside the GADT constructor that---      assembles them back into the original /same/ constructor.--- 3.   A function from a 'Prod' of 'BVar's (that contain the items inside---      the constructor) and doing whatever you wanted to do with it,---      inside 'BP'.------ If you don't provide all of the items inside the GADT into the 'BPC', or--- if your "re-assembling" function doesn't properly reassemble things--- correctly or changes some of the values, this will not work.----withGADT-    :: forall s rs a b. ()-    => BVar s rs a-    -> (a -> BPCont s rs a b)-    -> BP s rs b-withGADT v f = do-    x <- BP (resolveVar v)-    case f x of-      BPC (xs :: Tuple as) g h -> do-        let bp :: BPNode s rs '[a] as-            bp = BPN { _bpnOut       = map1 (const (FRInternal [])) xs-                     , _bpnRes       = xs-                     , _bpnGradFunc  = return . only_ . g-                                     . imap1 (\ix -> every @_ @Num ix // maybe (I 1) I)-                     , _bpnGradCache = Nothing-                     }-        r <- BP . liftBase $ newSTRef bp-        registerVar (IRNode IZ r) v-        h $ imap1 (\ix _ -> BVNode ix r) xs---- | Apply a function to the contents of an STRef, and cache the results--- using the given lens.  If already calculated, simply returned the cached--- result.-caching-    :: Lens' a (Maybe b)-    -> STRef s a-    -> (a -> ST s b)-    -> ST s b-caching l r f = do-    x <- readSTRef r-    let y = view l x-    case y of-      Just z ->-        return z-      Nothing -> do-        z <- f x-        modifySTRef r (set l (Just z))-        return z-{-# INLINE caching #-}-+{-# LANGUAGE FlexibleContexts #-}+{-# LANGUAGE GADTs            #-}+{-# LANGUAGE PatternSynonyms  #-}+{-# LANGUAGE RankNTypes       #-}++-- |+-- Module      : Numeric.Backprop+-- Copyright   : (c) Justin Le 2018+-- License     : BSD3+--+-- Maintainer  : justin@jle.im+-- Stability   : experimental+-- Portability : non-portable+--+-- Automatic differentation and backpropagation.+--+-- Main idea: Write a function computing what you want, and the library+-- automatically provies the gradient of that function as well, for usage+-- with gradient descent and other training methods.+--+-- In more detail: instead of working directly with values to produce your+-- result, you work with 'BVar's containing those values.  Working with+-- these 'BVar's is made smooth with the usage of lenses and other+-- combinators, and libraries can offer operatons on 'BVar's instead of+-- those on normal types directly.+--+-- Then, you can use:+--+-- @+-- 'evalBP' :: (forall s. 'Reifies' s 'W'. 'BVar' s a -> BVar s b) -> (a -> b)+-- @+--+-- to turn a 'BVar' function into the function on actual values @a -> b@.+-- This has virtually zero overhead over writing the actual function+-- directly.+--+-- Then, there's:+--+-- @+-- 'gradBP' :: (forall s. 'Reifies' s 'W'. 'BVar' s a -> BVar s b) -> (a -> a)+-- @+--+-- to automatically get the /gradient/, as well, for a given input.+--+-- See the <https://github.com/mstksg/backprop README> for more information+-- and links to demonstrations and tutorials, or dive striaght in by+-- reading the docs for 'BVar'.+--++module Numeric.Backprop (+    -- * Types+    BVar, W+    -- * Running+  , backprop, evalBP, gradBP+    -- ** Multiple inputs+  , backprop2, evalBP2, gradBP2+  , backpropN, evalBPN, gradBPN, Every+    -- * Manipulating 'BVar'+  , constVar+  , (^^.), (.~~), (^^?), (^^..)+  , viewVar, setVar+  , sequenceVar, collectVar+  , previewVar, toListOfVar+    -- ** With 'Op's#liftops#+    -- $liftops+  , liftOp+  , liftOp1, liftOp2, liftOp3+    -- * 'Op'+  , Op(..)+    -- ** Creation+  , op0, opConst, idOp+  , opConst'+    -- *** Giving gradients directly+  , op1, op2, op3+    -- *** From Isomorphisms+  , opCoerce, opTup, opIso, opLens+    -- * Utility+    -- ** Inductive tuples/heterogeneous lists+  , Prod(..), pattern (:>), only, head'+  , Tuple, pattern (::<), only_+  , I(..)+    -- ** Misc+  , Reifies+  ) where++import           Data.Bifunctor+import           Data.Reflection+import           Data.Type.Index+import           Lens.Micro+import           Numeric.Backprop.Internal+import           Numeric.Backprop.Op++-- $liftops+--+-- This library provides a few primitive actions for manipulating 'BVar's+-- and the values inside them, including its 'Num', 'Fractional', and+-- 'Floating' instances, and lens-based operations like '^^.', '.~~' '^^?',+-- and '^^..'.+--+-- However, the power of this library comes from manipulating many+-- different types from libraries, like matrices and vectors.  Libraries+-- can provide their own @'BVar' s a -> 'BVar' s b@ functions, alongside+-- (or in lieu of) @a -> b@ functions for their types.+--+-- The easiest way to create a 'BVar' function is to use 'liftOp' with an+-- 'Op' constructor.  For example, imagine a vector library providing a dot+-- product function.  We can write this using 'liftOp2' and 'op2':+--+-- @+-- dot :: 'BVar' s Vec -> BVar s Vec -> BVar s Double+-- dot = 'liftOp2' . op2 $ \\xs ys ->+--         ( sum (zipWith (*) xs ys)+--         , \\g -> (map (*g) ys, map (*g) xs)+--         )+-- @+--+-- We provide a function that, given the two inputs, returns:+--+--     (1) The result of the function on those two inputs+--     (2) A function taking the "total derivative", and returning the+--     gradient with respect to each of the inputs.+--+-- See documentation in "Numeric.Backprop.Op" for more information on the+-- second part (the gradient).+--+-- Nice 'Op's are how /backprop/ links together 'BVar's and tracks them to+-- determine their gradient.  Ideally, users would never have to deal with+-- these when backpropagating their own functions, and library authors+-- providing their matrix and vector operations, etc. would provide 'BVar'+-- variants of their normal operations.+--+-- In fact, 'BVar' operations could even be defined /instead/ of normal+-- operations, since it is easy to go from @'BVar' s a -> 'BVar' s b@ to @a+-- -> b@, using 'evalBP', and this carries virtually zero overhead, so some+-- libraries might even provide 'BVar' versions by default.++-- | Turn a function @'BVar' s a -> 'BVar' s b@ into the function @a -> b@+-- that it represents, also computing its gradient @a@ as well.+--+-- The Rank-N type @forall s. 'Reifies' s 'W' => ...@ is used to ensure+-- that 'BVar's do not leak out of the context (similar to how it is used+-- in "Control.Monad.ST"), and also as a reference to an ephemeral Wengert+-- tape used to track the graph of references.+--+-- Note that every type involved has to be an instance of 'Num'.  This is+-- because gradients all need to be "summable" (which is implemented using+-- 'sum' and '+'), and we also need to able to generate gradients of 1+-- and 0.  Really, only '+' and 'fromInteger' methods are used from the+-- 'Num' typeclass.+--+-- This might change in the future, to allow easier integration with tuples+-- (which typically do not have a 'Num' instance), and potentially make+-- types easier to use (by only requiring '+', 0, and 1, and not the rest+-- of the 'Num' class).+--+-- See the <https://github.com/mstksg/backprop README> for a more detailed+-- discussion on this issue.+--+-- If you need a 'Num' instance for tuples, consider the+-- <https://hackage.haskell.org/package/NumInstances NumInstances> package+-- (in particular, "Data.NumInstances.Tuple"), or else using a named+-- product type instead.+backprop+    :: forall a b. (Num a, Num b)+    => (forall s. Reifies s W => BVar s a -> BVar s b)+    -> a+    -> (b, a)+backprop f = second (getI . head')+           . backpropN (f . head')+           . only_+{-# INLINE backprop #-}++-- | Turn a function @'BVar' s a -> 'BVar' s b@ into the function @a -> b@+-- that it represents.+--+-- Benchmarks show that this should have virtually no overhead over+-- directly writing a @a -> b@. 'BVar' is, in this situation, a zero-cost+-- abstraction, performance-wise.+--+-- Has a nice advantage over using 'backprop' in that it doesn't require+-- 'Num' constraints on the input and output.+--+-- See documentation of 'backprop' for more information.+--+evalBP :: (forall s. Reifies s W => BVar s a -> BVar s b) -> a -> b+evalBP f = evalBPN (f . head') . only_+{-# INLINE evalBP #-}++-- | Take a function @'BVar' s a -> 'BVar' s b@, interpreted as a function+-- @a -> b@, and compute its gradient with respect to its input.+--+-- The resulting @a -> a@ tells how the input (and its components) affects+-- the output.  Positive numbers indicate that the result will vary in the+-- same direction as any adjustment in the input.  Negative numbers+-- indicate that the result will vary in the opposite direction as any+-- adjustment in the input.  Larger numbers indicate a greater sensitivity+-- of change, and small numbers indicate lower sensitivity.+--+-- See documentation of 'backprop' for more information.+--+gradBP+    :: forall a b. (Num a, Num b)+    => (forall s. Reifies s W => BVar s a -> BVar s b)+    -> a+    -> a+gradBP f = snd . backprop f+{-# INLINE gradBP #-}++-- | 'gradBP' generalized to multiple inputs of different types.  See+-- documentation for 'backpropN' for more details.+gradBPN+    :: forall as b. (Every Num as, Num b)+    => (forall s. Reifies s W => Prod (BVar s) as -> BVar s b)+    -> Tuple as+    -> Tuple as+gradBPN f = snd . backpropN f+{-# INLINE gradBPN #-}++-- | 'backprop' for a two-argument function.+--+-- Not strictly necessary, because you can always uncurry a function by+-- putting the arguments inside a data type, or using a tuple with+-- <https://hackage.haskell.org/package/NumInstances NumInstances>.+-- However, this can be convenient if you don't want to make a custom tuple+-- type or pull in orphan instances.  This could potentially also be more+-- performant.+--+-- For 3 and more arguments, consider using 'backpropN'.+backprop2+    :: forall a b c. (Num a, Num b, Num c)+    => (forall s. Reifies s W => BVar s a -> BVar s b -> BVar s c)+    -> a+    -> b+    -> (c, (a, b))+backprop2 f x y = second (\(dx ::< dy ::< Ø) -> (dx, dy))+                $ backpropN (\(x' :< y' :< Ø) -> f x' y') (x ::< y ::< Ø)+{-# INLINE backprop2 #-}++-- | 'evalBP' for a two-argument function.  See 'backprop2' for notes.+evalBP2+    :: (forall s. Reifies s W => BVar s a -> BVar s b -> BVar s c)+    -> a+    -> b+    -> c+evalBP2 f x y = evalBPN (\(x' :< y' :< Ø) -> f x' y') (x ::< y ::< Ø)+{-# INLINE evalBP2 #-}++-- | 'gradBP' for a two-argument function.  See 'backprop2' for notes.+gradBP2+    :: (Num a, Num b, Num c)+    => (forall s. Reifies s W => BVar s a -> BVar s b -> BVar s c)+    -> a+    -> b+    -> (a, b)+gradBP2 f x = snd . backprop2 f x+{-# INLINE gradBP2 #-}++-- | An infix version of 'viewVar', meant to evoke parallels to '^.' from+-- lens.+--+-- With normal values, you can extract something from that value with+-- a lens:+--+-- @+-- x '^.' myLens+-- @+--+-- would extract a piece of @x :: b@, specified by @myLens :: 'Lens'' b a@.+-- The result has type @a@.+--+-- @+-- xVar '^^.' myLens+-- @+--+-- would extract a piece out of @xVar :: 'BVar' s b@ (a 'BVar' holding a+-- @b@), specified by @myLens :: Lens' b a@.   The result has type @'BVar'+-- s a@ (a 'BVar' holding a @a@)+--+-- This is the main way to pull out values from 'BVar' of container types.+--+(^^.)+    :: forall a b s. (Reifies s W, Num a)+    => BVar s b+    -> Lens' b a+    -> BVar s a+x ^^. l = viewVar l x+infixl 8 ^^.+{-# INLINE (^^.) #-}++-- | An infix version of 'setVar', meant to evoke parallels to '.~' from+-- lens.+--+-- With normal values, you can set something in a value with a lens:+-- a lens:+--+-- @+-- x '&' myLens '.~' 'y'+-- @+--+-- would "set" a part of @x :: b@, specified by @myLens :: 'Lens'' a b@, to+-- a new value @y :: a@.+--+-- @+-- xVar '&' myLens '.~~' yVar+-- @+--+-- would "set" a part of @xVar :: 'BVar' s b@ (a 'BVar' holding a @b@),+-- specified by @myLens :: 'Lens'' a b@, to a new value given by @yVar ::+-- 'BVar' s a@.  The result is a new (updated) value of type @'BVar' s b@.+--+-- This is the main way to set values inside 'BVar's of container types.+--+(.~~)+    :: forall a b s. (Reifies s W, Num a, Num b)+    => Lens' b a+    -> BVar s a+    -> BVar s b+    -> BVar s b+l .~~ x = setVar l x+infixl 8 .~~+{-# INLINE (.~~) #-}++-- | An infix version of 'previewVar', meant to evoke parallels to '^?'+-- from lens.+--+-- With normal values, you can (potentially) extract something from that+-- value with a lens:+--+-- @+-- x '^?' myPrism+-- @+--+-- would (potentially) extract a piece of @x :: b@, specified by+-- @myPrism :: 'Traversal'' b a@. The result has type @'Maybe' a@.+--+-- @+-- xVar '^^?' myPrism+-- @+--+-- would (potentially) extract a piece out of @xVar :: 'BVar' s b@ (a+-- 'BVar' holding a @b@), specified by @myPrism :: Prism' b a@.+-- The result has type @'Maybe' ('BVar' s a)@ ('Maybe' a 'BVar' holding+-- a @a@).+--+-- This is intended to be used with 'Prism''s (which hits at most one+-- target), but will actually work with /any/ 'Traversal''.  If the+-- traversal hits more than one target, the first one found will be+-- extracted.+--+-- This can be used to "pattern match" on 'BVar's, by using prisms on+-- constructors.+--+-- Note that many automatically-generated prisms by the /lens/ package use+-- tuples, which cannot normally be backpropagated (because they do not+-- have a 'Num' instance).  However, you can pull in orphan instances from+-- <https://hackage.haskell.org/package/NumInstances NumInstances>, or also+-- chain those prisms with functions to convert tuples to your own custom+-- product types (or tuple types with 'Num' instances).+(^^?)+    :: forall b a s. (Num a, Reifies s W)+    => BVar s b+    -> Traversal' b a+    -> Maybe (BVar s a)+v ^^? t = previewVar t v+{-# INLINE (^^?) #-}++-- | An infix version of 'toListOfVar', meant to evoke parallels to '^..'+-- from lens.+--+-- With normal values, you can extract all targets of a 'Traversal' from+-- that value with a:+--+-- @+-- x '^..' myTraversal+-- @+--+-- would extract all targets inside of @x :: b@, specified by @myTraversal+-- :: 'Traversal'' b a@. The result has type @[a]@.+--+-- @+-- xVar '^^..' myTraversal+-- @+--+-- would extract all targets inside of @xVar :: 'BVar' s b@ (a 'BVar'+-- holding a @b@), specified by @myTraversal :: Traversal' b a@.   The result+-- has type @['BVar' s a]@ (A list of 'BVar's holding @a@s).+--+(^^..)+    :: forall b a s. (Num a, Reifies s W)+    => BVar s b+    -> Traversal' b a+    -> [BVar s a]+v ^^.. t = toListOfVar t v+{-# INLINE (^^..) #-}

@@ -1,393 +0,0 @@-{-# LANGUAGE DataKinds           #-}-{-# LANGUAGE FlexibleContexts    #-}-{-# LANGUAGE PatternSynonyms     #-}-{-# LANGUAGE RankNTypes          #-}-{-# LANGUAGE ScopedTypeVariables #-}-{-# LANGUAGE TypeApplications    #-}-{-# LANGUAGE TypeFamilies        #-}-{-# LANGUAGE TypeOperators       #-}---- |--- Module      : Numeric.Backprop.Implicit--- Copyright   : (c) Justin Le 2017--- License     : BSD3------ Maintainer  : justin@jle.im--- Stability   : experimental--- Portability : non-portable------ Offers full functionality for implicit-graph back-propagation.  The--- intended usage is to write a 'BPOp', which is a normal Haskell--- function from 'BVar's to a result 'BVar'. These 'BVar's can be--- manipulated using their 'Num' \/ 'Fractional' \/ 'Floating' instances.------ The library can then perform back-propagation on the function (using--- 'backprop' or 'grad') by using an implicitly built graph.------ This should actually be powerful enough for most use cases, but falls--- short for a couple of situations:------ 1. If the result of a function on 'BVar's is used twice--- (like @z@ in @let z = x * y in z + z@), this will allocate a new--- redundant graph node for every usage site of @z@.  You can explicitly--- /force/ @z@, but only using an explicit graph description using--- "Numeric.Backprop".------ 2. This can't handle sum types, like "Numeric.Backprop" can.  You can--- never pattern match on the constructors of a value inside a 'BVar'.  I'm--- not sure if this is a fundamental limitation (I suspect it might be) or--- if I just can't figure out how to implement it.  Suggestions welcome!------ As a comparison, this module offers functionality and an API very--- similar to "Numeric.AD.Mode.Reverse" from the /ad/ library, except for--- the fact that it can handle /heterogeneous/ values.------ Note that every type involved has to be an instance of 'Num'.  This is--- because gradients all need to be "summable" (which is implemented using--- 'sum' and '+'), and we also need to able to generate gradients of '1'--- and '0'.---module Numeric.Backprop.Implicit (-  -- * Types-  -- ** Backprop types-    BPOp, BVar, Op, OpB-  -- ** Tuple types-  -- | See "Numeric.Backprop#prod" for a mini-tutorial on 'Prod' and-  -- 'Tuple'-  , Prod(..), Tuple, I(..)-  -- * back-propagation-  , backprop, grad, eval-  -- * Var manipulation-  , BP.constVar, BP.liftB, (BP..$), BP.liftB1, BP.liftB2, BP.liftB3-  -- ** As Parts-  , partsVar, withParts-  , splitVars, gSplit, gTuple-  , partsVar', withParts'-  , splitVars', gSplit'-  -- * Op-  , BP.op1, BP.op2, BP.op3, BP.opN-  , BP.op1', BP.op2', BP.op3'-  -- * Utility-  , pattern (:>), only, head'-  , pattern (::<), only_-  -- ** Numeric Ops-  -- | Optimized ops for numeric functions.  See-  -- "Numeric.Backprop.Op#numops" for more information.-  , (+.), (-.), (*.), negateOp, absOp, signumOp-  , (/.), recipOp-  , expOp, logOp, sqrtOp, (**.), logBaseOp-  , sinOp, cosOp, tanOp, asinOp, acosOp, atanOp-  , sinhOp, coshOp, tanhOp, asinhOp, acoshOp, atanhOp-  ) where--import           Data.Type.Combinator-import           Data.Type.Index-import           Data.Type.Length-import           Data.Type.Product-import           Data.Type.Util-import           Lens.Micro hiding         (ix)-import           Lens.Micro.Extras-import           Numeric.Backprop.Internal-import           Numeric.Backprop.Iso-import           Numeric.Backprop.Op-import           Type.Class.Higher-import           Type.Class.Known-import           Type.Class.Witness-import qualified Generics.SOP              as SOP-import qualified Numeric.Backprop          as BP---- | An operation on 'BVar's that can be backpropagated. A value of type:------ @--- 'BPOp' rs a--- @------ takes a bunch of 'BVar's containg @rs@ and uses them to (purely) produce--- a 'BVar' containing an @a@.------ @--- foo :: 'BPOp' '[ Double, Double ] Double--- foo (x ':<' y ':<' 'Ø') = x + sqrt y--- @------ 'BPOp' here is related to 'Numeric.Backprop.BPOpI' from the normal--- explicit-graph backprop module "Numeric.Backprop".-type BPOp rs a = forall s. Prod (BVar s rs) rs -> BVar s rs a---- | Run back-propagation on a 'BPOp' function, getting both the result and--- the gradient of the result with respect to the inputs.------ @--- foo :: 'BPOp' '[Double, Double] Double--- foo (x :< y :< Ø) =---   let z = x * sqrt y---   in  z + x ** y--- @------ >>> 'backprop' foo (2 ::< 3 ::< Ø)--- (11.46, 13.73 ::< 6.12 ::< Ø)-backprop-    :: Every Num rs-    => BPOp rs a-    -> Tuple rs-    -> (a, Tuple rs)-backprop f xs = BP.backprop (BP.withInps' (prodLength xs) (return . f)) xs---- | Run the 'BPOp' on an input tuple and return the gradient of the result--- with respect to the input tuple.------ @--- foo :: 'BPOp' '[Double, Double] Double--- foo (x :< y :< Ø) =---   let z = x * sqrt y---   in  z + x ** y--- @------ >>> grad foo (2 ::< 3 ::< Ø)--- 13.73 ::< 6.12 ::< Ø-grad-    :: Every Num rs-    => BPOp rs a-    -> Tuple rs-    -> Tuple rs-grad f = snd . backprop f---- | Simply run the 'BPOp' on an input tuple, getting the result without--- bothering with the gradient or with back-propagation.------ @--- foo :: 'BPOp' '[Double, Double] Double--- foo (x :< y :< Ø) =---   let z = x * sqrt y---   in  z + x ** y--- @------ >>> eval foo (2 ::< 3 ::< Ø)--- 11.46-eval-    :: (Known Length rs, Num a)-    => BPOp rs a-    -> Tuple rs-    -> a-eval f = BP.evalBPOp $ BP.implicitly f---- | A version of 'partsVar' taking explicit 'Length', indicating the--- number of items in the input tuple and their types.------ Requiring an explicit 'Length' is mostly useful for rare "extremely--- polymorphic" situations, where GHC can't infer the type and length of--- the internal tuple.  If you ever actually explicitly write down @bs@ as--- a list of types, you should be able to just use 'partsVar'.-partsVar'-    :: forall s rs bs a. Every Num bs-    => Length bs-    -> Iso' a (Tuple bs)-    -> BVar s rs a-    -> Prod (BVar s rs) bs-partsVar' l i r = map1 (\ix -> every @_ @Num ix //-                                 BP.liftB1 (BP.op1' (f ix)) r-                       ) ixes-  where-    f :: Num b-      => Index bs b-      -> a-      -> (b, Maybe b -> a)-    f ix x = ( getI . index ix . view i $ x-             , review i-             . flip (set (indexP ix)) zeroes-             . maybe (I 1) I-             )-    zeroes :: Tuple bs-    zeroes = map1 (\ix -> I 0 \\ every @_ @Num ix) ixes-    ixes :: Prod (Index bs) bs-    ixes = indices' l---- | Use an 'Iso' (or compatible 'Control.Lens.Iso.Iso' from the lens--- library) to "pull out" the parts of a data type and work with each part--- as a 'BVar'.------ If there is an isomorphism between a @b@ and a @'Tuple' as@ (that is, if--- an @a@ is just a container for a bunch of @as@), then it lets you break--- out the @as@ inside and work with those.------ @--- data Foo = F Int Bool------ fooIso :: 'Iso'' Foo (Tuple '[Int, Bool])--- fooIso = 'iso' (\\(F i b)         -\> i ::\< b ::\< Ø)---              (\\(i ::\< b ::\< Ø) -\> F i b        )------ 'partsVar' fooIso :: 'BVar' rs Foo -> 'Prod' ('BVar' s rs) '[Int, Bool]------ stuff :: 'BPOp' s '[Foo] a--- stuff (foo :< Ø) =---     case 'partsVar' fooIso foo of---       i :< b :< Ø ->---         -- now, i is a 'BVar' pointing to the 'Int' inside foo---         -- and b is a 'BVar' pointing to the 'Bool' inside foo---         -- you can do stuff with the i and b here--- @------ You can use this to pass in product types as the environment to a 'BP',--- and then break out the type into its constituent products.------ Note that for a type like @Foo@, @fooIso@ can be generated automatically--- with 'GHC.Generics.Generic' from "GHC.Generics" and--- 'Generics.SOP.Generic' from "Generics.SOP" and /generics-sop/, using the--- 'gTuple' iso.  See 'gSplit' for more information.------ Also, if you are literally passing a tuple (like--- @'BP' s '[Tuple '[Int, Bool]@) then you can give in the identity--- isomorphism ('id') or use 'splitVars'.------ At the moment, this implicit 'partsVar' is less efficient than the--- explicit 'Numeric.Backprop.partsVar', but this might change in the--- future.-partsVar-    :: forall s rs bs a. (Every Num bs, Known Length bs)-    => Iso' a (Tuple bs)-    -> BVar s rs a-    -> Prod (BVar s rs) bs-partsVar = partsVar' known---- | A version of 'withParts' taking explicit 'Length', indicating the--- number of internal items and their types.------ Requiring an explicit 'Length' is mostly useful for rare "extremely--- polymorphic" situations, where GHC can't infer the type and length of--- the internal tuple.  If you ever actually explicitly write down @bs@ as--- a list of types, you should be able to just use 'withParts'.-withParts'-    :: forall s rs bs a r. Every Num bs-    => Length bs-    -> Iso' a (Tuple bs)-    -> BVar s rs a-    -> (Prod (BVar s rs) bs -> r)-    -> r-withParts' l i r f = f (partsVar' l i r)---- | A continuation-based version of 'partsVar'.  Instead of binding the--- parts and using it in the rest of the block, provide a continuation to--- handle do stuff with the parts inside.------ Building on the example from 'partsVar':------ @--- data Foo = F Int Bool------ fooIso :: 'Iso'' Foo (Tuple '[Int, Bool])--- fooIso = 'iso' (\\(F i b)         -\> i ::\< b ::\< Ø)---              (\\(i ::\< b ::\< Ø) -\> F i b        )------ stuff :: 'BPOp' s '[Foo] a--- stuff (foo :< Ø) = 'withParts' fooIso foo $ \\case---     i :\< b :< Ø -\>---       -- now, i is a 'BVar' pointing to the 'Int' inside foo---       -- and b is a 'BVar' pointing to the 'Bool' inside foo---       -- you can do stuff with the i and b here--- @------ Mostly just a stylistic alternative to 'partsVar'.-withParts-    :: forall s rs bs a r. (Every Num bs, Known Length bs)-    => Iso' a (Tuple bs)-    -> BVar s rs a-    -> (Prod (BVar s rs) bs -> r)-    -> r-withParts = withParts' known---- | A version of 'splitVars' taking explicit 'Length', indicating the--- number of internal items and their types.------ Requiring an explicit 'Length' is mostly useful for rare "extremely--- polymorphic" situations, where GHC can't infer the type and length of--- the internal tuple.  If you ever actually explicitly write down @as@ as--- a list of types, you should be able to just use 'splitVars'.-splitVars'-    :: forall s rs as. Every Num as-    => Length as-    -> BVar s rs (Tuple as)-    -> Prod (BVar s rs) as-splitVars' l = partsVar' l id---- | Split out a 'BVar' of a tuple into a tuple ('Prod') of 'BVar's.------ @--- -- the environment is a single Int-Bool tuple, tup--- stuff :: 'BPOp' s '[ Tuple '[Int, Bool] ] a--- stuff (tup :< Ø) =---   case 'splitVar' tup of---     i :< b :< Ø <- 'splitVars' tup---     -- now, i is a 'BVar' pointing to the 'Int' inside tup---     -- and b is a 'BVar' pointing to the 'Bool' inside tup---     -- you can do stuff with the i and b here--- @------ Note that------ @--- 'splitVars' = 'partsVar' 'id'--- @-splitVars-    :: forall s rs as. (Every Num as, Known Length as)-    => BVar s rs (Tuple as)-    -> Prod (BVar s rs) as-splitVars = splitVars' known---- | A version of 'gSplit' taking explicit 'Length', indicating the--- number of internal items and their types.------ Requiring an explicit 'Length' is mostly useful for rare "extremely--- polymorphic" situations, where GHC can't infer the type and length of--- the internal tuple.  If you ever actually explicitly write down @as@ as--- a list of types, you should be able to just use 'gSplit'.-gSplit'-    :: forall s rs as a. (SOP.Generic a, SOP.Code a ~ '[as], Every Num as)-    => Length as-    -> BVar s rs a-    -> Prod (BVar s rs) as-gSplit' l = partsVar' l gTuple---- | Using 'GHC.Generics.Generic' from "GHC.Generics" and--- 'Generics.SOP.Generic' from "Generics.SOP", /split/ a 'BVar' containing--- a product type into a tuple ('Prod') of 'BVar's pointing to each value--- inside.------ Building on the example from 'partsVar':------ @--- import qualified Generics.SOP as SOP------ data Foo = F Int Bool---   deriving Generic------ instance SOP.Generic Foo------ 'gSplit' :: 'BVar' rs Foo -> 'Prod' ('BVar' s rs) '[Int, Bool]------ stuff :: 'BPOp' s '[Foo] a--- stuff (foo :< Ø) =---     case 'gSplit' foo of---       i :< b :< Ø ->---         -- now, i is a 'BVar' pointing to the 'Int' inside foo---         -- and b is a 'BVar' pointing to the 'Bool' inside foo---         -- you can do stuff with the i and b here--- @------ Because @Foo@ is a straight up product type, 'gSplit' can use--- "GHC.Generics" and take out the items inside.------ Note that------ @--- 'gSplit' = 'splitVars' 'gTuple'--- @-gSplit-    :: forall s rs as a. (SOP.Generic a, SOP.Code a ~ '[as], Every Num as, Known Length as)-    => BVar s rs a-    -> Prod (BVar s rs) as-gSplit = gSplit' known---- TODO: figure out how to split sums--- TODO: refactor these out to not need Known Length

@@ -1,23 +1,19 @@-{-# LANGUAGE DeriveFunctor              #-}-{-# LANGUAGE FlexibleContexts           #-}-{-# LANGUAGE FlexibleInstances          #-}-{-# LANGUAGE GADTs                      #-}-{-# LANGUAGE GeneralizedNewtypeDeriving #-}-{-# LANGUAGE LambdaCase                 #-}-{-# LANGUAGE MultiParamTypeClasses      #-}-{-# LANGUAGE PolyKinds                  #-}-{-# LANGUAGE RankNTypes                 #-}-{-# LANGUAGE ScopedTypeVariables        #-}-{-# LANGUAGE TemplateHaskell            #-}-{-# LANGUAGE TypeApplications           #-}-{-# LANGUAGE TypeFamilies               #-}-{-# LANGUAGE TypeInType                 #-}-{-# LANGUAGE TypeOperators              #-}-{-# LANGUAGE UndecidableInstances       #-}+{-# LANGUAGE BangPatterns        #-}+{-# LANGUAGE DeriveGeneric       #-}+{-# LANGUAGE FlexibleContexts    #-}+{-# LANGUAGE GADTs               #-}+{-# LANGUAGE LambdaCase          #-}+{-# LANGUAGE RankNTypes          #-}+{-# LANGUAGE RecordWildCards     #-}+{-# LANGUAGE ScopedTypeVariables #-}+{-# LANGUAGE TupleSections       #-}+{-# LANGUAGE TypeApplications    #-}+{-# LANGUAGE TypeInType          #-}+{-# LANGUAGE ViewPatterns        #-}  -- | -- Module      : Numeric.Backprop.Internal--- Copyright   : (c) Justin Le 2017+-- Copyright   : (c) Justin Le 2018 -- License     : BSD3 -- -- Maintainer  : justin@jle.im@@ -27,291 +23,659 @@ -- Provides the types and instances used for the graph -- building/back-propagation for the library. -module Numeric.Backprop.Internal-  ( OpB-  , BPState(..), bpsSources-  , BP(..)-  , BPInpRef(..)-  , BPNode(..), bpnOut, bpnRes, bpnGradFunc, bpnGradCache-  , BPPipe(..), bppOut, bppRes, bppGradFunc, bppGradCache-  , BVar(..)-  , ForwardRefs(..), _FRInternal+module Numeric.Backprop.Internal (+    BVar+  , W+  , backpropN, evalBPN+  , constVar+  , liftOp, liftOp1, liftOp2, liftOp3+  , viewVar, setVar, sequenceVar, collectVar, previewVar, toListOfVar+  -- * Debug+  , debugSTN+  , debugIR   ) where -import           Control.Monad.Reader+import           Control.DeepSeq+import           Control.Exception+import           Control.Monad+import           Control.Monad.Primitive import           Control.Monad.ST-import           Control.Monad.State+import           Control.Monad.Trans.State+import           Data.Bifunctor+import           Data.Foldable+import           Data.IORef import           Data.Kind-import           Data.STRef+import           Data.Maybe+import           Data.Monoid+import           Data.Proxy+import           Data.Reflection import           Data.Type.Index-import           Data.Type.Product-import           Lens.Micro hiding    (ix)-import           Lens.Micro.TH+import           Data.Type.Product hiding   (toList)+import           Data.Type.Util+import           Data.Type.Vector hiding    (itraverse, head')+import           GHC.Generics+import           Lens.Micro import           Numeric.Backprop.Op+import           System.IO.Unsafe+import           Type.Class.Higher+import           Type.Class.Witness+import           Unsafe.Coerce+import qualified Data.Vector                as V+import qualified Data.Vector.Mutable        as MV --- | A subclass of 'OpM' (and superclass of 'Op'), representing 'Op's that--- the /backprop/ library uses to perform backpropation.+-- | A @'BVar' s a@ is a value of type @a@ that can be "backpropagated". ----- An+-- Functions referring to 'BVar's are tracked by the library and can be+-- automatically differentiated to get their gradients and results. --+-- For simple numeric values, you can use its 'Num', 'Fractional', and+-- 'Floating' instances to manipulate them as if they were the numbers they+-- represent.+--+-- If @a@ contains items, the items can be accessed and extracted using+-- lenses. A @'Lens'' b a@ can be used to access an @a@ inside a @b@, using+-- '^^.' ('viewVar'):+-- -- @--- 'OpB' s rs a+-- ('^.')  ::        a -> 'Lens'' a b ->        b+-- ('^^.') :: 'BVar' s a -> 'Lens'' a b -> 'BVar' s b -- @ ----- represents a differentiable function that takes a tuple of @rs@ and--- produces an a @a@, which can be run on @'BVar' s@s and also inside @'BP'--- s@s.  For example, an @'OpB' s '[ Int, Double ] Bool@ takes an 'Int' and--- a 'Double' and produces a 'Bool', and does it in a differentiable way.+-- There is also '^^?' ('previewVar'), to use a 'Prism'' or 'Traversal'' to+-- extract a target that may or may not be present (which can implement+-- pattern matching), '^^..' ('toListOfVar') to use a 'Traversal'' to+-- extract /all/ targets inside a 'BVar', and '.~~' ('setVar') to set and+-- update values inside a 'BVar'. ----- 'OpB' is a /superset/ of 'Op', so, if you see any function--- that expects an 'OpB' (like 'Numeric.Backprop.opVar'' and--- 'Numeric.Backprop.~$', for example), you can give them an 'Op', as well.+-- For more complex operations, libraries can provide functions on 'BVar's+-- using 'liftOp' and related functions.  This is how you can create+-- primitive functions that users can use to manipulate your library's+-- values. ----- You can think of 'OpB' as a superclass/parent class of 'Op' in this--- sense, and of 'Op' as a subclass of 'OpB'.-type OpB s as a = OpM (ST s) as a+-- For example, the /hmatrix/ library has a matrix-vector multiplication+-- function, @#> :: L m n -> R n -> L m@.+--+-- A library could instead provide a function @#> :: 'BVar' (L m n) -> BVar+-- (R n) -> BVar (R m)@, which the user can then use to manipulate their+-- 'BVar's of @L m n@s and @R n@s, etc.+--+-- See "Numeric.Backprop#liftops" and documentation for 'liftOp' for more+-- information.+--+data BVar s a = BV { _bvRef :: !(BRef s)+                   , _bvVal :: !a+                   } --- | Reference to /usage sites/ for a given entity, used to get partial or--- total derivatives.-data ForwardRefs s rs a-    -- | A list of 'BPInpRef's pointing to places that use the entity, to-    -- provide partial derivatives.-    = FRInternal ![BPInpRef s rs a]-    -- | The entity is the terminal result of a BP, so its total derivative-    -- is fixed.-    | FRTerminal !(Maybe a)+data BRef (s :: Type) = BRInp !Int+                      | BRIx !Int+                      | BRC+  deriving (Generic, Show) --- | Combines two 'FRInternal' lists.  If either input is an 'FRTerminal',--- then throws away the other result and keeps the new terminal forced--- total derivative.  (Biases to the left)-instance Monoid (ForwardRefs s rs a) where-    mempty  = FRInternal []-    mappend = \case-        FRInternal rs -> \case-          FRInternal rs'   -> FRInternal (rs ++ rs')-          t@(FRTerminal _) -> t-        FRTerminal _  -> id+instance NFData (BRef s) --- | The "state" of a 'BP' action, which keeps track of what nodes, if any,--- refer to any of the inputs.-data BPState :: Type -> [Type] -> Type where-    BPS :: { _bpsSources :: !(Prod (ForwardRefs s rs) rs)-           }-        -> BPState s rs+-- | This will force the value inside, as well.+instance NFData a => NFData (BVar s a) where+    rnf (BV r v) = force r `seq` force v `seq` () --- | A Monad allowing you to explicitly build hetereogeneous data--- dependency graphs and that the library can perform back-propagation on.+-- | Project out a constant value if the 'BVar' refers to one.+bvConst :: BVar s a -> Maybe a+bvConst (BV BRC !x) = Just x+bvConst _           = Nothing+{-# INLINE bvConst #-}++forceBVar :: BVar s a -> ()+forceBVar (BV !r !_) = force r `seq` ()+{-# INLINE forceBVar #-}++data InpRef :: Type -> Type where+    IR :: Num a+       => { _irIx  :: !(BVar s b)+          , _irUpd :: !(Lens' b a)+          }+       -> InpRef a++forceInpRef :: InpRef a -> ()+forceInpRef (IR !v !_) = forceBVar v `seq` ()+{-# INLINE forceInpRef #-}++-- | Debugging string for an 'InpRef'.+debugIR :: InpRef a -> String+debugIR IR{..} = show (_bvRef _irIx)++data TapeNode :: Type -> Type where+    TN :: { _tnInputs :: !(Prod InpRef as)+          , _tnGrad   :: !(a -> Tuple as)+          }+       -> TapeNode a++forceTapeNode :: TapeNode a -> ()+forceTapeNode (TN !inps !_) = foldMap1 forceInpRef inps `seq` ()+{-# INLINE forceTapeNode #-}++data SomeTapeNode :: Type where+    STN :: forall a. Num a+        => !(TapeNode a)+        -> SomeTapeNode++forceSomeTapeNode :: SomeTapeNode -> ()+forceSomeTapeNode (STN !tn) = forceTapeNode tn `seq` ()+{-# INLINE forceSomeTapeNode #-}++-- | Debugging string for a 'SomeTapeMode'.+debugSTN :: SomeTapeNode -> String+debugSTN (STN TN{..}) = show . foldMap1 ((:[]) . debugIR) $ _tnInputs++-- | An ephemeral Wengert Tape in the environment.  Used internally to+-- track of the computational graph of variables. ----- A @'BP' s rs a@ is a 'BP' action that uses an environment of @rs@--- returning a @a@.  When "run", it will compute a gradient that is a tuple--- of @rs@.  (The phantom parameter @s@ is used to ensure that any 'BVar's--- aren't leaked out of the monad)+-- For the end user, one can just imagine @'Reifies' s 'W'@ as a required+-- constraint on @s@ that allows backpropagation to work.+newtype W = W { wRef :: IORef (Int, [SomeTapeNode]) }++initWengert :: IO W+initWengert = W <$> newIORef (0,[])+{-# INLINE initWengert #-}++insertNode+    :: Num a+    => TapeNode a+    -> a+    -> W+    -> IO (BVar s a)+insertNode !tn !x !w = fmap ((`BV` x) . BRIx) . atomicModifyIORef' (wRef w) $ \(!(!n,!t)) ->+    let n' = n + 1+        t' = STN tn:t+    in  forceTapeNode tn `seq` n' `seq` t' `seq` ((n', t'), n)+{-# INLINE insertNode #-}++-- | Lift a value into a 'BVar' representing a constant value. ----- Note that you can only "run" a @'BP' s rs@ that produces a 'BVar' ----- that is, things of the form+-- This value will not be considered an input, and its gradients will not+-- be backpropagated.+constVar :: a -> BVar s a+constVar = BV BRC+{-# INLINE constVar #-}++liftOp_+    :: forall s as b. (Reifies s W, Num b, Every Num as)+    => Op as b+    -> Prod (BVar s) as+    -> IO (BVar s b)+liftOp_ o !vs = case traverse1 (fmap I . bvConst) vs of+                   Just xs -> return $ constVar (evalOp o xs)+                   Nothing -> insertNode tn y (reflect (Proxy @s))+  where+    (y,g) = runOpWith o (map1 (I . _bvVal) vs)+    tn = TN { _tnInputs = imap1 go vs+            , _tnGrad   = g+            }+    go :: forall a. Index as a -> BVar s a -> InpRef a+    go i !v = forceBVar v `seq` (IR v id \\ every @_ @Num i)+{-# INLINE liftOp_ #-}++-- | Lift an 'Op' with an arbitrary number of inputs to a function on the+-- appropriate number of 'BVar's. ----- @--- 'BP' s rs ('BVar' s rs a)--- @+-- Should preferably be used only by libraries to provide primitive 'BVar'+-- functions for their types for users. ----- The above is a 'BP' action that returns a 'BVar' containing an @a@.--- When this is run, it'll produce a result of type @a@ and a gradient of--- that is a tuple of @rs@.  (This form has a type synonym,--- 'Numeric.Backprop.BPOp', for convenience)+-- See "Numeric.Backprop#liftops" and documentation for 'liftOp' for more+-- information, and "Numeric.Backprop.Op#prod" for a mini-tutorial on using+-- 'Prod' and 'Tuple'.+liftOp+    :: forall s as b. (Reifies s W, Num b, Every Num as)+    => Op as b+    -> Prod (BVar s) as+    -> BVar s b+liftOp o !vs = unsafePerformIO $ liftOp_ o vs+{-# INLINE liftOp #-}++liftOp1_+    :: forall s a b. (Reifies s W, Num a, Num b)+    => Op '[a] b+    -> BVar s a+    -> IO (BVar s b)+liftOp1_ o (bvConst->Just x) = return . constVar . evalOp o $ (x ::< Ø)+liftOp1_ o !v = forceBVar v `seq` insertNode tn y (reflect (Proxy @s))+  where+    (y,g) = runOpWith o (_bvVal v ::< Ø)+    tn = TN { _tnInputs = IR v id :< Ø+            , _tnGrad   = g+            }+{-# INLINE liftOp1_ #-}++-- | Lift an 'Op' with a single input to be a function on a single 'BVar'. ----- For example, a @'BP' s '[ Int, Double, Double ]@ is a monad that--- represents a computation with an 'Int', 'Double', and 'Double' as--- inputs.   And, if you ran a+-- Should preferably be used only by libraries to provide primitive 'BVar'+-- functions for their types for users. ----- @--- 'BP' s '[ Int, Double, Double ] ('BVar' s '[ Int, Double, Double ] Double)--- @+-- See "Numeric.Backprop#liftops" and documentation for 'liftOp' for more+-- information.+liftOp1+    :: forall s a b. (Reifies s W, Num a, Num b)+    => Op '[a] b+    -> BVar s a+    -> BVar s b+liftOp1 o !v = unsafePerformIO $ liftOp1_ o v+{-# INLINE liftOp1 #-}++liftOp2_+    :: forall s a b c. (Reifies s W, Num a, Num b, Num c)+    => Op '[a,b] c+    -> BVar s a+    -> BVar s b+    -> IO (BVar s c)+liftOp2_ o (bvConst->Just x) (bvConst->Just y) = return . constVar . evalOp o $ x ::< y ::< Ø+liftOp2_ o !v !u = forceBVar v+             `seq` forceBVar u+             `seq` insertNode tn y (reflect (Proxy @s))+  where+    (y,g) = runOpWith o (_bvVal v ::< _bvVal u ::< Ø)+    tn = TN { _tnInputs = IR v id :< IR u id :< Ø+            , _tnGrad   = g+            }+{-# INLINE liftOp2_ #-}++-- | Lift an 'Op' with two inputs to be a function on a two 'BVar's. ----- Or, using the 'BPOp' type synonym:+-- Should preferably be used only by libraries to provide primitive 'BVar'+-- functions for their types for users. ----- @--- 'Numeric.Backprop.BPOp' s '[ Int, Double, Double ] Double--- @+-- See "Numeric.Backprop#liftops" and documentation for 'liftOp' for more+-- information.+liftOp2+    :: forall s a b c. (Reifies s W, Num a, Num b, Num c)+    => Op '[a,b] c+    -> BVar s a+    -> BVar s b+    -> BVar s c+liftOp2 o !v !u = unsafePerformIO $ liftOp2_ o v u+{-# INLINE liftOp2 #-}++liftOp3_+    :: forall s a b c d. (Reifies s W, Num a, Num b, Num c, Num d)+    => Op '[a,b,c] d+    -> BVar s a+    -> BVar s b+    -> BVar s c+    -> IO (BVar s d)+liftOp3_ o (bvConst->Just x) (bvConst->Just y) (bvConst->Just z)+    = return . constVar . evalOp o $ x ::< y ::< z ::< Ø+liftOp3_ o !v !u !w = forceBVar v+                `seq` forceBVar u+                `seq` forceBVar w+                `seq` insertNode tn y (reflect (Proxy @s))+  where+    (y, g) = runOpWith o (_bvVal v ::< _bvVal u ::< _bvVal w ::< Ø)+    tn = TN { _tnInputs = IR v id :< IR u id :< IR w id :< Ø+            , _tnGrad   = g+            }+{-# INLINE liftOp3_ #-}++-- | Lift an 'Op' with three inputs to be a function on a three 'BVar's. ----- with 'Numeric.Backprop.backprop' or 'Numeric.Backprop.gradBPOp', it'll--- return a gradient on the inputs ('Int', 'Double', and 'Double') and--- produce a value of type 'Double'.+-- Should preferably be used only by libraries to provide primitive 'BVar'+-- functions for their types for users. ----- Now, one powerful thing about this type is that a 'BP' is itself an--- 'Op' (or more precisely, an 'Numeric.Backprop.OpB', which is a subtype of--- 'OpM').  So, once you create your fancy 'BP' computation, you can--- transform it into an 'OpM' using 'Numeric.Backprop.bpOp'.-newtype BP s rs a-    = BP { bpST :: ReaderT (Tuple rs) (StateT (BPState s rs) (ST s)) a }-    deriving (Functor, Applicative, Monad)+-- See "Numeric.Backprop#liftops" and documentation for 'liftOp' for more+-- information.+liftOp3+    :: forall s a b c d. (Reifies s W, Num a, Num b, Num c, Num d)+    => Op '[a,b,c] d+    -> BVar s a+    -> BVar s b+    -> BVar s c+    -> BVar s d+liftOp3 o !v !u !w = unsafePerformIO $ liftOp3_ o v u w+{-# INLINE liftOp3 #-} --- | The basic unit of manipulation inside 'BP' (or inside an--- implicit-graph backprop function).  Instead of directly working with--- values, you work with 'BVar's contating those values.  When you work--- with a 'BVar', the /backprop/ library can keep track of what values--- refer to which other values, and so can perform back-propagation to--- compute gradients.+viewVar_+    :: forall a b s. (Reifies s W, Num a)+    => Lens' b a+    -> BVar s b+    -> IO (BVar s a)+viewVar_ l !v = forceBVar v `seq` insertNode tn y (reflect (Proxy @s))+  where+    y = _bvVal v ^. l+    tn = TN { _tnInputs = IR v l :< Ø+            , _tnGrad   = only_+            }+{-# INLINE viewVar_ #-}++-- | Using a 'Lens'', extract a value /inside/ a 'BVar'.  Meant to evoke+-- parallels to 'view' from lens. ----- A @'BVar' s rs a@ refers to a value of type @a@, with an environment--- of values of the types @rs@.  The phantom parameter @s@ is used to--- ensure that stray 'BVar's don't leak outside of the backprop process.+-- See documentation for '^^.' for more information.+viewVar+    :: forall a b s. (Reifies s W, Num a)+    => Lens' b a+    -> BVar s b+    -> BVar s a+viewVar l !v = unsafePerformIO $ viewVar_ l v+{-# INLINE viewVar #-}++setVar_+    :: forall a b s. (Reifies s W, Num a, Num b)+    => Lens' b a+    -> BVar s a+    -> BVar s b+    -> IO (BVar s b)+setVar_ l !w !v = forceBVar v+            `seq` forceBVar w+            `seq` insertNode tn y (reflect (Proxy @s))+  where+    y = _bvVal v & l .~ _bvVal w+    tn = TN { _tnInputs = IR w id :< IR v id :< Ø+            , _tnGrad   = \d -> let (dw,dv) = l (,0) d+                                in  dw ::< dv ::< Ø+            }+{-# INLINE setVar_ #-}++-- | Using a 'Lens'', set a value /inside/ a 'BVar'.  Meant to evoke+-- parallels to "set" from lens. ----- (That is, if you're using implicit backprop, it ensures that you interact--- with 'BVar's in a polymorphic way.  And, if you're using explicit--- backprop, it ensures that a @'BVar' s rs a@ never leaves the @'BP' s rs@--- that it was created in.)+-- See documentation for '.~~' for more information.+setVar+    :: forall a b s. (Reifies s W, Num a, Num b)+    => Lens' b a+    -> BVar s a+    -> BVar s b+    -> BVar s b+setVar l !w !v = unsafePerformIO $ setVar_ l w v+{-# INLINE setVar #-}++-- | Extract all of the 'BVar's out of a 'Traversable' container of+-- 'BVar's.+sequenceVar+    :: forall t a s. (Reifies s W, Traversable t, Num a)+    => BVar s (t a)+    -> t (BVar s a)+sequenceVar !v = unsafePerformIO $ traverseVar' id traverse v+{-# INLINE sequenceVar #-}++collectVar_+    :: forall a t s. (Reifies s W, Foldable t, Functor t, Num (t a), Num a)+    => t (BVar s a)+    -> IO (BVar s (t a))+collectVar_ !vs = withV (toList vs) $ \(vVec :: Vec n (BVar s a)) -> do+    let tn :: TapeNode (t a)+        tn = TN { _tnInputs = vecToProd (vmap ((`IR` id) . getI) vVec)+                , _tnGrad   = maybe (error "distributeVar") vecToProd+                            . listToVec (vecLen vVec)+                            . map I . toList+                }+    traverse_ (evaluate . forceBVar) vs+    insertNode tn (_bvVal <$> vs) (reflect (Proxy @s))+{-# INLINE collectVar_ #-}++-- | Collect all of the 'BVar's in a container into a 'BVar' of that+-- container's contents.+collectVar+    :: forall a t s. (Reifies s W, Foldable t, Functor t, Num (t a), Num a)+    => t (BVar s a)+    -> BVar s (t a)+collectVar !vs = unsafePerformIO $ collectVar_ vs+{-# INLINE collectVar #-}++traverseVar'+    :: forall b a f s. (Num a, Reifies s W, Traversable f)+    => (b -> f a)+    -> Traversal' b a+    -> BVar s b+    -> IO (f (BVar s a))+traverseVar' f t !v = forceBVar v+              `seq` itraverse go (f (_bvVal v))+  where+    go :: Int -> a -> IO (BVar s a)+    go i y = insertNode tn y (reflect (Proxy @s))+      where+        tn = TN { _tnInputs = IR v (ixt t i) :< Ø+                , _tnGrad   = only_+                }+{-# INLINE traverseVar' #-}++-- | Using a 'Traversal'', extract a single value /inside/ a 'BVar', if it+-- exists.  If more than one traversal target exists, returns te first.+-- Meant to evoke parallels to 'preview' from lens.  Really only intended+-- to be used wth 'Prism''s, or up-to-one target traversals. ----- 'BVar's have 'Num', 'Fractional', 'Floating', etc. instances, so they--- can be manipulated using polymorphic functions and numeric functions in--- Haskell.  You can add them, subtract them, etc., in "implicit" backprop--- style.+-- See documentation for '^^?' for more information.+previewVar+    :: forall b a s. (Num a, Reifies s W)+    => Traversal' b a+    -> BVar s b+    -> Maybe (BVar s a)+previewVar t !v = unsafePerformIO $ traverseVar' (listToMaybe . toListOf t) t v+{-# INLINE previewVar #-}++-- | Using a 'Traversal'', extract all targeted values /inside/ a 'BVar'.+-- Meant to evoke parallels to 'toListOf' from lens. ----- (However, note that if you directly manipulate 'BVar's using those--- instances or using 'Numeric.Backprop.liftB', it delays evaluation, so every usage site--- has to re-compute the result/create a new node.  If you want to re-use--- a 'BVar' you created using '+' or '-' or 'Numeric.Backprop.liftB', use--- 'Numeric.Backprop.bindVar' to force it first.  See documentation for--- 'Numeric.Backprop.bindVar' for more details.)-data BVar :: Type -> [Type] -> Type -> Type where-    -- | A BVar referring to a 'BPNode'-    BVNode  :: !(Index bs a)-            -> !(STRef s (BPNode s rs as bs))-            -> BVar s rs a-    -- | A BVar referring to an environment input variable-    BVInp   :: !(Index rs a)-            -> BVar s rs a-    -- | A constant BVar that refers to a specific Haskell value-    BVConst :: !a-            -> BVar s rs a-    -- | A BVar that combines several other BVars using a function (an-    -- 'Op').  Essentially a branch of a tree.-    BVOp    :: !(Prod (BVar s rs) as)-            -> !(OpB s as a)-            -> BVar s rs a+-- See documentation for '^^..' for more information.+toListOfVar+    :: forall b a s. (Num a, Reifies s W)+    => Traversal' b a+    -> BVar s b+    -> [BVar s a]+toListOfVar t !v = unsafePerformIO $ traverseVar' (toListOf t) t v+{-# INLINE toListOfVar #-} --- | Used exclusively by 'ForwardRefs' to specify "where" and "how" to look--- for partial derivatives at usage sites of a given entity.-data BPInpRef :: Type -> [Type] -> Type -> Type where-    -- | The entity is used in a 'BPNode', and as an Nth input-    IRNode  :: Every Num cs-            => !(Index bs a)-            -> !(STRef s (BPNode s rs bs cs))-            -> BPInpRef s rs a-    -- | The entity is used in a 'BPPipe', and as an Nth input-    IRPipe  :: !(Index bs a)-            -> !(STRef s (BPPipe s rs bs cs))-            -> BPInpRef s rs a-    -- | The entity is used somehow in the terminal result of a 'BP', and-    -- so therefore has a fixed partial derivative contribution.-    IRConst :: !a-            -> BPInpRef s rs a+data SomeNum :: Type where+    SN  :: Num a+        => Proxy a+        -> a+        -> SomeNum --- | A (stateful) node in the graph of operations/data dependencies in 'BP'--- that the library uses.  'BVar's can refer to these to get results from--- them, and 'BPInpRef's can refer to these to get partial derivatives from--- them.-data BPNode :: Type -> [Type] -> [Type] -> [Type] -> Type where-    BPN :: { _bpnOut       :: !(Prod (ForwardRefs s rs) bs)-           , _bpnRes       :: !(Tuple bs)-           , _bpnGradFunc  :: !(Prod Maybe bs -> ST s (Tuple as))-           , _bpnGradCache :: !(Maybe (Tuple as))  -- nothing if is the "final output"-           }-        -> BPNode s rs as bs+data Runner s = R { _rDelta  :: MV.MVector s SomeNum+                  , _rInputs :: MV.MVector s SomeNum+                  } --- | Essentially a "single-usage" 'BPNode'.  It's a stateful node, but only--- ever has a single consumer (and so its total derivative comes from--- a single partial derivative).  Used when keeping track of 'BVOp's.-data BPPipe :: Type -> [Type] -> [Type] -> [Type] -> Type where-    BPP :: { _bppOut       :: !(Prod (BPInpRef s rs) bs)-           , _bppRes       :: !(Tuple bs)-           , _bppGradFunc  :: !(Tuple bs -> ST s (Tuple as))-           , _bppGradCache :: !(Maybe (Tuple as))-           }-        -> BPPipe s rs as bs+initRunner+    :: (PrimMonad m, PrimState m ~ s)+    => (Int, [SomeTapeNode])+    -> (Int, [Some (Wit1 Num)])+    -> m (Runner s)+initRunner (n, stns) (nx,xs) = do+    delts <- MV.new n+    for_ (zip [n-1,n-2..] stns) $ \(i, STN (TN{..} :: TapeNode c)) -> do+      MV.write delts i $ SN (Proxy @c) 0+    inps <- MV.new nx+    for_ (zip [0..] xs) $ \(i, Some (Wit1 :: Wit1 Num c)) -> do+      MV.write inps i $ SN (Proxy @c) 0+    return $ R delts inps+{-# INLINE initRunner #-} -makeLenses ''BPState-makeLenses ''BPNode-makeLenses ''BPPipe+gradRunner+    :: forall m b s p. (PrimMonad m, PrimState m ~ s, Num b)+    => p b+    -> Runner s+    -> (Int, [SomeTapeNode])+    -> m ()+gradRunner _ R{..} (n,stns) = do+    when (n > 0) $+      MV.write _rDelta (n - 1) (SN (Proxy @b) 1)+    zipWithM_ go [n-1,n-2..] stns+  where+    go :: Int -> SomeTapeNode -> m ()+    go i (STN TN{..}) = do+      SN _ delt  <- MV.read _rDelta i+      let gs = _tnGrad (unsafeCoerce delt)+      zipWithPM_ propagate _tnInputs gs+    propagate :: forall x. InpRef x -> I x -> m ()+    propagate (IR v ln) (I !d) = case _bvRef v of+      BRInp !i -> flip (MV.modify _rInputs) i $ \case+        SN p !y -> let y' = unsafeCoerce y & ln %~ (+d)+                   in  y' `seq` SN p (unsafeCoerce y')+      BRIx !i -> flip (MV.modify _rDelta) i $ \case+        SN p !y -> let y' = unsafeCoerce y & ln %~ (+d)+                   in  y' `seq` SN p (unsafeCoerce y')+      BRC     -> return ()+{-# INLINE gradRunner #-} --- | Traversal (fake prism) to refer to the list of internal refs if the--- 'ForwardRef' isn't associated with a terminal entity.-_FRInternal-    :: Traversal (ForwardRefs s as a) (ForwardRefs t bs a)-                 [BPInpRef s as a]    [BPInpRef t bs a]-_FRInternal f = \case-    FRInternal xs -> FRInternal <$> f xs-    FRTerminal g  -> pure (FRTerminal g)+-- | 'backprop' generalized to multiple inputs of different types.  See the+-- "Numeric.Backprop.Op#prod" for a mini-tutorial on heterogeneous lists.+--+-- Not strictly necessary, because you can always uncurry a function by+-- passing in all of the inputs in a data type containing all of the+-- arguments.   You could also pass in a giant tuple with+-- <https://hackage.haskell.org/package/NumInstances NumInstances>.+-- However, this can be convenient if you don't want to make a custom tuple+-- type or pull in orphan instances.  This could potentially also be more+-- performant.+--+-- A @'Prod' ('BVar' s) '[Double, Float, Double]@, for instance, is a tuple+-- of @'BVar' s 'Double'@, @'BVar' s 'Float'@, and @'BVar' s 'Double'@, and+-- can be pattern matched on using ':<' (cons) and 'Ø' (nil).+--+-- Tuples can be built and pattern matched on using '::<' (cons) and 'Ø'+-- (nil), as well.+--+-- The @'Every' 'Num' as@ in the constraint says that every value in the+-- type-level list @as@ must have a 'Num' instance.  This means you can+-- use, say, @'[Double, Float, Int]@, but not @'[Double, Bool, String]@.+--+-- If you stick to /concerete/, monomorphic usage of this (with specific+-- types, typed into source code, known at compile-time), then @'Every'+-- 'Num' as@ should be fulfilled automatically.+--+backpropN+    :: forall as b. (Every Num as, Num b)+    => (forall s. Reifies s W => Prod (BVar s) as -> BVar s b)+    -> Tuple as+    -> (b, Tuple as)+backpropN f xs = (y, g)+  where+    !(!tp@(!_,!_),!y) = unsafePerformIO $ fillWengert f xs+    g :: Tuple as+    g = runST $ do+        r <- initRunner tp (getSum `first` ifoldMap1 go xs)+        gradRunner (Proxy @b) r tp+        delts <- toList <$> V.freeze (_rInputs r)+        return . fromMaybe (error "backpropN") $+          fillProd (\_ (SN _ d) -> I (unsafeCoerce d)) xs delts+      where+        go :: forall a. Index as a -> I a -> (Sum Int, [Some (Wit1 Num)])+        go i (I _) = (1, [Some (Wit1 :: Wit1 Num a)]) \\ every @_ @Num i+{-# INLINE backpropN #-} +-- | 'evalBP' generalized to multiple inputs of different types.  See+-- documentation for 'backpropN' for more details.+evalBPN+    :: forall as b. ()+    => (forall s. Reifies s W => Prod (BVar s) as -> BVar s b)+    -> Tuple as+    -> b+evalBPN f = snd . unsafePerformIO . fillWengert f+{-# INLINE evalBPN #-} +fillWengert+    :: forall as b. ()+    => (forall s. Reifies s W => Prod (BVar s) as -> BVar s b)+    -> Tuple as+    -> IO ((Int, [SomeTapeNode]), b)+fillWengert f xs = do+    w <- initWengert+    o <- reify w $ \(Proxy :: Proxy s) -> do+      let oVar = f (inpProd @s)+      evaluate (forceBVar oVar)+      return (_bvVal oVar)+    t <- readIORef (wRef w)+    traverse_ (evaluate . forceSomeTapeNode) (snd t)+    return (t, o)+  where+    inpProd :: forall s. Prod (BVar s) as+    inpProd = evalState (traverse1 (state . go . getI) xs) 0+      where+        go :: a -> Int -> (BVar s a, Int)+        go x i = (BV (BRInp i) x, i + 1)+{-# INLINE fillWengert #-}  --- | Note that if you use the 'Num' instance to create 'BVar's, the--- resulting 'BVar' is deferred/delayed.  At every location you use it, it--- will be recomputed, and a separate graph node will be created.  If you--- are using a 'BVar' you made with the 'Num' instance in multiple--- locations, use 'Numeric.Backprop.bindVar' first to force it and prevent--- recomputation.-instance Num a => Num (BVar s rs a) where-    r1 + r2       = BVOp (r1 :< r2 :< Ø) (+.)+instance (Num a, Reifies s W) => Num (BVar s a) where+    (+)         = liftOp2 (+.)     {-# INLINE (+) #-}-    r1 - r2       = BVOp (r1 :< r2 :< Ø) (-.)+    (-)         = liftOp2 (-.)     {-# INLINE (-) #-}-    r1 * r2       = BVOp (r1 :< r2 :< Ø) (*.)+    (*)         = liftOp2 (*.)     {-# INLINE (*) #-}-    negate r      = BVOp (r  :< Ø)       negateOp+    negate      = liftOp1 negateOp     {-# INLINE negate #-}-    signum r      = BVOp (r  :< Ø)       signumOp+    signum      = liftOp1 signumOp     {-# INLINE signum #-}-    abs r         = BVOp (r  :< Ø)       absOp+    abs         = liftOp1 absOp     {-# INLINE abs #-}-    fromInteger x = BVConst (fromInteger x)+    fromInteger = constVar . fromInteger     {-# INLINE fromInteger #-} --- | See note for 'Num' instance.-instance Fractional a => Fractional (BVar s rs a) where-    r1 / r2        = BVOp (r1 :< r2 :< Ø) (/.)+instance (Fractional a, Reifies s W) => Fractional (BVar s a) where+    (/)          = liftOp2 (/.)     {-# INLINE (/) #-}-    recip r        = BVOp (r  :< Ø)       recipOp+    recip        = liftOp1 recipOp     {-# INLINE recip #-}-    fromRational x = BVConst (fromRational x)+    fromRational = constVar . fromRational     {-# INLINE fromRational #-} --- | See note for 'Num' instance.-instance Floating a => Floating (BVar s rs a) where-    pi            = BVConst pi+instance (Floating a, Reifies s W) => Floating (BVar s a) where+    pi      = constVar pi     {-# INLINE pi #-}-    exp   r       = BVOp (r :< Ø)        expOp+    exp     = liftOp1 expOp     {-# INLINE exp #-}-    log   r       = BVOp (r :< Ø)        logOp+    log     = liftOp1 logOp     {-# INLINE log #-}-    sqrt  r       = BVOp (r :< Ø)        sqrtOp+    sqrt    = liftOp1 sqrtOp     {-# INLINE sqrt #-}-    r1 ** r2      = BVOp (r1 :< r2 :< Ø) (**.)+    (**)    = liftOp2 (**.)     {-# INLINE (**) #-}-    logBase r1 r2 = BVOp (r1 :< r2 :< Ø) logBaseOp+    logBase = liftOp2 logBaseOp     {-# INLINE logBase #-}-    sin   r       = BVOp (r :< Ø)        sinOp+    sin     = liftOp1 sinOp     {-# INLINE sin #-}-    cos   r       = BVOp (r :< Ø)        cosOp+    cos     = liftOp1 cosOp     {-# INLINE cos #-}-    tan   r       = BVOp (r :< Ø)        tanOp+    tan     =  liftOp1 tanOp     {-# INLINE tan  #-}-    asin  r       = BVOp (r :< Ø)        asinOp+    asin    = liftOp1 asinOp     {-# INLINE asin #-}-    acos  r       = BVOp (r :< Ø)        acosOp+    acos    = liftOp1 acosOp     {-# INLINE acos #-}-    atan  r       = BVOp (r :< Ø)        atanOp+    atan    = liftOp1 atanOp     {-# INLINE atan #-}-    sinh  r       = BVOp (r :< Ø)        sinhOp+    sinh    = liftOp1 sinhOp     {-# INLINE sinh #-}-    cosh  r       = BVOp (r :< Ø)        coshOp+    cosh    = liftOp1 coshOp     {-# INLINE cosh #-}-    tanh  r       = BVOp (r :< Ø)        tanhOp+    tanh    = liftOp1 tanhOp     {-# INLINE tanh #-}-    asinh r       = BVOp (r :< Ø)        asinhOp+    asinh   = liftOp1 asinhOp     {-# INLINE asinh #-}-    acosh r       = BVOp (r :< Ø)        acoshOp+    acosh   = liftOp1 acoshOp     {-# INLINE acosh #-}-    atanh r       = BVOp (r :< Ø)        atanhOp+    atanh   = liftOp1 atanhOp     {-# INLINE atanh #-} +-- Some utility functions to get around a lens dependency+itraverse+    :: forall t a b f. (Traversable t, Monad f)+    => (Int -> a -> f b) -> t a -> f (t b)+itraverse f xs = evalStateT (traverse (StateT . go) xs) 0+  where+    go :: a -> Int -> f (b, Int)+    go x i = (,i+1) <$> f i x+{-# INLINE itraverse #-}++ixi :: Int -> Lens' [a] a+ixi _ _ []     = error "ixi"+ixi 0 f (x:xs) = (:xs) <$> f x+ixi n f (x:xs) = (x:)  <$> ixi (n - 1) f xs+{-# INLINE ixi #-}++ixt :: forall b a. Traversal' b a -> Int -> Lens' b a+ixt t i f xs = stuff <$> ixi i f contents+  where+    contents = xs ^.. t+    stuff    = evalState (traverseOf t (state . const go) xs)+      where+        go :: [a] -> (a,  [a])+        go []     = error "asList"+        go (y:ys) = (y, ys)+{-# INLINE ixt #-}