stripe / rainier

Bayesian inference in Scala.
https://rainier.fit
Apache License 2.0
432 stars 51 forks source link

Decimal implementation with double and fraction reprs #432

Closed avibryant closed 4 years ago

avibryant commented 4 years ago

@sritchie I'd love your eyes on this one

This is a follow up to https://github.com/stripe/rainier/pull/431 . For context, this would not be merged into develop, but into a 0.3-dev branch.

It introduces two representations for constants: one is just a double, and the other is a numerator/denominator which are both long. If you create a Real with an integer or long, it will stay in the fraction form as long as you continue to do operations on it with other ints/longs. As soon as you introduce a double, it will convert into double and stay there (we never check to see if we can convert the other way).

This should catch the most common cases where we otherwise lose precision through floating point errors, without hampering performance of the common case of doubles. From simple benchmarking, this is not appreciably slower than just using doubles directly in Real.

Note: I disabled the CategoricalSuite in cats because it was failing on a -Inf != -Inf, and I wasn't sure how best to address that.

sritchie commented 4 years ago

@avibryant, there are two main typeclasses I use with a categorical distribution that might belong in rainier - Expectation and Decompose. I'm curious if there is a more natural way to accomplish this, so let me describe the two traits.

Two main things make sense to model, sometimes, as categorical distributions:

Expectation

I use this Expectation typeclass to generate expected values, given some Categorical. It's generic so it can work for cats.Id[A] instances, my categorical distribution, rainier's, or some new thing that comes along.

I could back Double out to Decimal or Real now.

Decompose

I added a Decompose typeclass to break some M[A] down into a linear combination of (A, R) pairs... experiment outcome and probability. That feels like an idea that comes up a lot.

trait Decompose[M[_], R] extends Serializable {
  def ring: Ring[R]
  def decompose[A](ma: M[A]): Iterator[(A, R)]
}

What do you think? Is there a better way to get an expected value out of a categorical distribution? My stats vocab is not great (yet!) so it might be staring me in the face.

sritchie commented 4 years ago

@avibryant, is this why those tests are failing? Do we now have a Constant(Double.NEG_INF) or whatever that is in not-pseudocode?

  def eqReal(epsilon: Double): Eq[Real] = {
    val bde = eqDecimal(epsilon)
    Eq.instance { (left, right) =>
      (left, right) match {
        case (Infinity, Infinity) | (NegInfinity, NegInfinity) => true
        case (Constant(a), Constant(b))                        => bde.eqv(a, b)
        // TODO - the Eq instance returned here doesn't test for full
        // DAG equality; this final case needs to be expanded to
        // compare the various NonConstant options.
        case (a, b) => a == b
      }
    }
  }

https://github.com/stripe/rainier/blob/73dee58cd36c0229611f8b12bb90638061ae9abc/rainier-cats/src/main/scala/com/stripe/rainier/cats/package.scala#L125

avibryant commented 4 years ago

@sritchie yes, but I don't think we necessarily want to paper over it by explicitly checking for infinities. Previously we introduced the Infinity and NegInfinity subtypes of Real because BigDecimal cannot represent infinities, and tried to make it impossible to ever end up with Constant(inf). Somehow the change to Decimal is making that possible again here, which we may want to treat as a bug.

sritchie commented 4 years ago

@avibryant yup, looks like the Eq instance giving us trouble here too https://travis-ci.org/stripe/rainier/jobs/626809074#L362

avibryant commented 4 years ago

@sritchie re Expectation and Decompose, the regime Rainier is designed for is one where you cannot compute exact probabilities for outcomes, but instead can only generate samples - so expectations would be estimated by averaging over the sampled outcomes, rather than ever having probabilities explicitly attached. Categorical with constant probabilities (ie where the probabilities are not themselves continuous random variables) is an extreme special case where you have discrete outcomes with known probabilities and so can actually compute probabilities exactly; but that case is the least interesting one to rainier and it wouldn't make sense to have code for it.

sritchie commented 4 years ago

Yeah, that makes sense. This came up in the more boring contrived bits of tabular RL, where you keep cycling through the sequence of all states and updating your estimates based on expectations of what reward and estimated next state value each action might get you.

A more interesting case that I do need to cover next is a case where you have:

And you can only take a single action. You take the action, and eventually you get some return... then you update your estimate if the value of the action by

that is a more interesting use case. Is that only useful as a conceptual model for something that exists in Rainier? We're sampling and updating estimates based on the same, but we still have this problem of aggregating using the probabilities as weights.

I'm not really sure what I'm asking... I can write code that will work, BUT I'd rather shift my mental model and reuse existing concepts. I'd love any pointers!

avibryant commented 4 years ago

@sritchie there's a whole bunch of interesting stuff here to talk about, though we should probably move it to a different medium. The short version, though, is that what you're describing - both because it involves mostly discrete distributions, and because it's about modeling an open-ended process that is probably intricately tied into a simulation - is a really poor fit for Rainier's design choices (Rainier is largely built around the constraints of Hamiltonian Monte Carlo, which is for continuous parameters, fixed model structure, and seeing all the data at once - basically the polar opposite of what you're doing).

In the probabilistic programming world, the techniques that are better suited to this are particle filters, aka Sequential Monte Carlo. There's lot of systems out there geared toward that - in fact, the first version of Rainier was, based on porting this paper to Scala: http://mlg.eng.cam.ac.uk/pub/pdf/SciGhaGor15.pdf Another one that looks really cool is https://birch-lang.org/ (which is its own language), and there's https://probprog.github.io/anglican/ for Clojure.

avibryant commented 4 years ago

@sritchie I tried hard to get a relaxed enough Eq[Decimal] to let the Categorical suite run, but whatever is generating random examples is doing too good a job of exposing floating point numerical instability and it's just not tenable to run the tests.

codecov-io commented 4 years ago

Codecov Report

Merging #432 into 0.3-dev will increase coverage by 1.65%. The diff coverage is 89.33%.

Impacted file tree graph

@@             Coverage Diff             @@
##           0.3-dev     #432      +/-   ##
===========================================
+ Coverage    58.74%   60.39%   +1.65%     
===========================================
  Files           69       70       +1     
  Lines         2642     2886     +244     
  Branches       153      146       -7     
===========================================
+ Hits          1552     1743     +191     
- Misses        1090     1143      +53
Impacted Files Coverage Δ
...main/scala/com/stripe/rainier/compute/ToReal.scala 73.68% <75%> (+0.35%) :arrow_up:
...c/main/scala/com/stripe/rainier/cats/package.scala 64.4% <85.71%> (-20.7%) :arrow_down:
...ain/scala/com/stripe/rainier/compute/Decimal.scala 91.37% <91.22%> (-3.63%) :arrow_down:
...c/main/scala/com/stripe/rainier/core/package.scala 0% <0%> (-100%) :arrow_down:
.../scala/com/stripe/rainier/scalacheck/package.scala 63.15% <0%> (-15.79%) :arrow_down:
...in/scala/com/stripe/rainier/core/Categorical.scala 36.95% <0%> (-7.49%) :arrow_down:
...ain/scala/com/stripe/rainier/ir/DataFunction.scala 100% <0%> (ø) :arrow_up:
...src/main/scala/com/stripe/rainier/compute/Fn.scala 100% <0%> (ø)
...in/scala/com/stripe/rainier/compute/Compiler.scala 92.68% <0%> (+0.09%) :arrow_up:
... and 7 more

Continue to review full report at Codecov.

Legend - Click here to learn more Δ = absolute <relative> (impact), ø = not affected, ? = missing data Powered by Codecov. Last update 9b5ea65...7ff6a6e. Read the comment docs.